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

Changeset 14013


Ignore:
Timestamp:
2020-12-02T16:23:27+01:00 (3 years ago)
Author:
jchanut
Message:

Manually merge with Christian's branch agrif top sponge/interp #2129

Location:
NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST/agrif_oce_interp.F90

    r14000 r14013  
    940940         ELSE 
    941941          
     942            IF ( Agrif_Parent(ln_zps) ) THEN ! Account for partial cells  
     943                                             ! linear vertical interpolation 
     944               DO jj=j1,j2 
     945                  DO ji=i1,i2 
     946                     ! 
     947                     N_in  = mbkt(ji,jj) 
     948                     N_out = mbkt(ji,jj) 
     949                     z_in(1) = ptab(ji,jj,1,n2) 
     950                     tabin(1,1:jpts) = ptab(ji,jj,1,1:jpts) 
     951                     DO jk=2, N_in 
     952                        z_in(jk) = ptab(ji,jj,jk,n2) 
     953                        tabin(jk,1:jpts) = ptab(ji,jj,jk,1:jpts) 
     954                     END DO 
     955                     IF (.NOT.ln_linssh) z_in(1:N_in) = z_in(1:N_in) - ptab(ji,jj,k2,n2) 
     956                     z_out(1) = 0.5_wp * e3t(ji,jj,1,Krhs_a) 
     957                     DO jk=2, N_out 
     958                        z_out(jk) = z_out(jk-1) + 0.5_wp * (e3t(ji,jj,jk-1,Krhs_a) + e3t(ji,jj,jk,Krhs_a)) 
     959                     END DO 
     960                     IF (.NOT.ln_linssh) z_out(1:N_out) = z_out(1:N_out) - ssh(ji,jj,Krhs_a) 
     961                     CALL remap_linear(tabin(1:N_in,1:jpts),z_in(1:N_in),ptab(ji,jj,1:N_out,1:jpts), & 
     962                                   &   z_out(1:N_out),N_in,N_out,jpts)   
     963                  END DO 
     964               END DO 
     965            ENDIF 
     966 
    942967            DO jn =1, jpts 
    943  
    944                IF ( Agrif_Parent(ln_zps) ) THEN ! Account for partial cells  
    945                                                 ! linear vertical interpolation 
    946                   DO jj=j1,j2 
    947                      DO ji=i1,i2 
    948                         ! 
    949                         N_in  = mbkt(ji,jj) 
    950                         N_out = mbkt(ji,jj) 
    951                         z_in(1) = ptab(ji,jj,1,n2) 
    952                         tabin(1,1:jpts) = ptab(ji,jj,1,1:jpts) 
    953                         DO jk=2, N_in 
    954                            z_in(jk) = ptab(ji,jj,jk,n2) 
    955                            tabin(jk,1:jpts) = ptab(ji,jj,jk,1:jpts) 
    956                         END DO 
    957                         IF (.NOT.ln_linssh) z_in(1:N_in) = z_in(1:N_in) - ptab(ji,jj,k2,n2) 
    958                         z_out(1) = 0.5_wp * e3t(ji,jj,1,Krhs_a) 
    959                         DO jk=2, N_out 
    960                            z_out(jk) = z_out(jk-1) + 0.5_wp * (e3t(ji,jj,jk-1,Krhs_a) + e3t(ji,jj,jk,Krhs_a)) 
    961                         END DO 
    962                         IF (.NOT.ln_linssh) z_out(1:N_out) = z_out(1:N_out) - ssh(ji,jj,Krhs_a) 
    963                         CALL remap_linear(tabin(1:N_in,1:jpts),z_in(1:N_in),ptab(ji,jj,1:N_out,1:jpts), & 
    964                                       &   z_out(1:N_out),N_in,N_out,jpts)   
    965                      END DO 
    966                   END DO 
    967  
    968                ENDIF 
    969  
    970968               ts(i1:i2,j1:j2,1:jpk,jn,Krhs_a) = ptab(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk) 
    971  
    972969            END DO 
    973970         ENDIF 
  • NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST/agrif_top_interp.F90

    r13352 r14013  
    4343      Agrif_SpecialValue    = 0._wp 
    4444      Agrif_UseSpecialValue = .TRUE. 
    45       l_vremap = ln_vert_remap 
     45      l_vremap              = ln_vert_remap 
    4646      ! 
    4747      CALL Agrif_Bc_variable( trn_id, procname=interptrn ) 
     
    6464      INTEGER  :: item 
    6565      ! vertical interpolation: 
    66       REAL(wp) :: zhtot 
    67       REAL(wp), DIMENSION(k1:k2,1:jptra) :: tabin 
    68       REAL(wp), DIMENSION(k1:k2) :: h_in, z_in 
     66      REAL(wp) :: zhtot, zwgt 
     67      REAL(wp), DIMENSION(k1:k2,1:jptra) :: tabin, tabin_i 
     68      REAL(wp), DIMENSION(k1:k2) :: z_in, h_in_i, z_in_i 
    6969      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out 
    7070      !!---------------------------------------------------------------------- 
     
    8585         END DO 
    8686 
    87          IF( l_vremap .OR. l_ini_child) THEN 
    88             ! Interpolate thicknesses 
     87         IF( l_vremap .OR. l_ini_child .OR. ln_zps ) THEN 
     88            ! Fill cell depths (i.e. gdept) to be interpolated 
    8989            ! Warning: these are masked, hence extrapolated prior interpolation. 
    90             DO jk=k1,k2 
    91                DO jj=j1,j2 
    92                   DO ji=i1,i2 
    93                       ptab(ji,jj,jk,jptra+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) 
    94  
     90            DO jj=j1,j2 
     91               DO ji=i1,i2 
     92                  ptab(ji,jj,k1,jptra+1) = 0.5_wp * tmask(ji,jj,k1) * e3t(ji,jj,k1,Kmm_a) 
     93                  DO jk=k1+1,k2 
     94                     ptab(ji,jj,jk,jptra+1) = tmask(ji,jj,jk) * & 
     95                        & ( ptab(ji,jj,jk-1,jptra+1) + 0.5_wp * (e3t(ji,jj,jk-1,Kmm_a)+e3t(ji,jj,jk,Kmm_a)) ) 
    9596                  END DO 
    9697               END DO 
    9798            END DO 
    9899 
    99             ! Extrapolate thicknesses in partial bottom cells: 
    100             ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 
    101             IF (ln_zps) THEN 
    102                DO jj=j1,j2 
    103                   DO ji=i1,i2 
    104                       jk = mbkt(ji,jj) 
    105                       ptab(ji,jj,jk,jptra+1) = 0._wp 
    106                   END DO 
    107                END DO            
    108             END IF 
    109          
    110100            ! Save ssh at last level: 
    111101            IF (.NOT.ln_linssh) THEN 
    112102               ptab(i1:i2,j1:j2,k2,jptra+1) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1)  
    113             ELSE 
    114                ptab(i1:i2,j1:j2,k2,jptra+1) = 0._wp 
    115103            END IF       
    116104         ENDIF 
     
    126114            DO jj=j1,j2 
    127115               DO ji=i1,i2 
    128                   tr(ji,jj,:,:,Krhs_a) = 0.                   
     116                  tr(ji,jj,:,:,Krhs_a) = 0.   
     117                  ! 
     118                  ! Build vertical grids: 
    129119                  N_in = mbkt_parent(ji,jj) 
    130                   zhtot = 0._wp 
    131                   DO jk=1,N_in !k2 = jpk of parent grid 
    132                      IF (jk==N_in) THEN 
    133                         h_in(jk) = ht0_parent(ji,jj) + ptab(ji,jj,k2,n2) - zhtot 
    134                      ELSE 
    135                         h_in(jk) = ptab(ji,jj,jk,n2) 
    136                      ENDIF 
    137                      zhtot = zhtot + h_in(jk) 
    138                      tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1) 
    139                   END DO 
    140                   z_in(1) = 0.5_wp * h_in(1) - zhtot + ht0_parent(ji,jj) 
     120                  ! Input grid (account for partial cells if any): 
     121                  DO jk=1,N_in 
     122                     z_in(jk) = ptab(ji,jj,jk,n2) - ptab(ji,jj,k2,n2) 
     123                     tabin(jk,1:jptra) = ptab(ji,jj,jk,1:jptra) 
     124                  END DO 
     125                   
     126                  ! Intermediate grid: 
     127                  DO jk = 1, N_in 
     128                     h_in_i(jk) = e3t0_parent(ji,jj,jk) * &  
     129                       &       (1._wp + ptab(ji,jj,k2,n2)/(ht0_parent(ji,jj)*ssmask(ji,jj) + 1._wp - ssmask(ji,jj))) 
     130                  END DO 
     131                  z_in_i(1) = 0.5_wp * h_in_i(1) 
    141132                  DO jk=2,N_in 
    142                      z_in(jk) = z_in(jk-1) + 0.5_wp * (h_in(jk-1)+h_in(jk)) 
    143                   END DO 
    144  
    145                   N_out = 0 
    146                   DO jk=1,jpk ! jpk of child grid 
    147                      IF (tmask(ji,jj,jk) == 0._wp) EXIT  
    148                      N_out = N_out + 1 
     133                     z_in_i(jk) = z_in_i(jk-1) + 0.5_wp * ( h_in_i(jk) + h_in_i(jk-1) ) 
     134                  END DO 
     135                  z_in_i(1:N_in) = z_in_i(1:N_in)  - ptab(ji,jj,k2,n2) 
     136 
     137                  ! Output (Child) grid: 
     138                  N_out = mbkt(ji,jj) 
     139                  DO jk=1,N_out 
    149140                     h_out(jk) = e3t(ji,jj,jk,Krhs_a) 
    150141                  END DO 
    151  
    152                   z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + ht_0(ji,jj) 
     142                  z_out(1) = 0.5_wp * h_out(1) 
    153143                  DO jk=2,N_out 
    154                      z_out(jk) = z_out(jk-1) + 0.5_wp * (h_out(jk-1)+h_out(jk)) 
    155                   END DO 
     144                     z_out(jk) = z_out(jk-1) + 0.5_wp * ( h_out(jk)+h_out(jk-1) ) 
     145                  END DO 
     146                  IF (.NOT.ln_linssh) z_out(1:N_out) = z_out(1:N_out)  - ssh(ji,jj,Krhs_a)                
    156147 
    157148                  IF (N_in*N_out > 0) THEN 
     
    159150                        CALL remap_linear(tabin(1:N_in,1:jptra),z_in(1:N_in),tr(ji,jj,1:N_out,1:jptra,Krhs_a),          & 
    160151                                      &   z_out(1:N_out),N_in,N_out,jptra)   
    161                      ELSE  
    162                         CALL reconstructandremap(tabin(1:N_in,1:jptra),h_in(1:N_in),tr(ji,jj,1:N_out,1:jptra,Krhs_a),   & 
    163                                       &   h_out(1:N_out),N_in,N_out,jptra)   
     152                     ELSE   
     153                        CALL remap_linear(tabin(1:N_in,1:jptra),z_in(1:N_in),tabin_i(1:N_in,1:jptra),                     & 
     154                                     &   z_in_i(1:N_in),N_in,N_in,jptra) 
     155                        CALL reconstructandremap(tabin_i(1:N_in,1:jptra),h_in_i(1:N_in),tr(ji,jj,1:N_out,1:jptra,Krhs_a), & 
     156                                      &   h_out(1:N_out),N_in,N_out,jptra)    
    164157                     ENDIF 
    165158                  ENDIF 
     
    169162  
    170163         ELSE 
    171           
     164 
     165            IF ( Agrif_Parent(ln_zps) ) THEN ! Account for partial cells  
     166                                             ! linear vertical interpolation 
     167               DO jj=j1,j2 
     168                  DO ji=i1,i2 
     169                     ! 
     170                     N_in  = mbkt(ji,jj) 
     171                     N_out = mbkt(ji,jj) 
     172                     z_in(1) = ptab(ji,jj,1,n2) 
     173                     tabin(1,1:jptra) = ptab(ji,jj,1,1:jptra) 
     174                     DO jk=2, N_in 
     175                        z_in(jk) = ptab(ji,jj,jk,n2) 
     176                        tabin(jk,1:jptra) = ptab(ji,jj,jk,1:jptra) 
     177                     END DO 
     178                     IF (.NOT.ln_linssh) z_in(1:N_in) = z_in(1:N_in) - ptab(ji,jj,k2,n2) 
     179                     z_out(1) = 0.5_wp * e3t(ji,jj,1,Krhs_a) 
     180                     DO jk=2, N_out 
     181                        z_out(jk) = z_out(jk-1) + 0.5_wp * (e3t(ji,jj,jk-1,Krhs_a) + e3t(ji,jj,jk,Krhs_a)) 
     182                     END DO 
     183                     IF (.NOT.ln_linssh) z_out(1:N_out) = z_out(1:N_out) - ssh(ji,jj,Krhs_a) 
     184                     CALL remap_linear(tabin(1:N_in,1:jptra),z_in(1:N_in),ptab(ji,jj,1:N_out,1:jptra), & 
     185                                   &   z_out(1:N_out),N_in,N_out,jptra)   
     186                  END DO 
     187               END DO 
     188 
     189            ENDIF 
     190 
    172191            DO jn=1, jptra 
    173192                tr(i1:i2,j1:j2,1:jpk,jn,Krhs_a)=ptab(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk)  
  • NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST/agrif_top_sponge.F90

    r13987 r14013  
    4545      ! 
    4646#if defined SPONGE_TOP 
    47 !! Assume persistence  
     47      !! Assume persistence: 
    4848      zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot()) 
    49       CALL Agrif_sponge 
     49 
    5050      Agrif_SpecialValue    = 0._wp 
    5151      Agrif_UseSpecialValue = .TRUE. 
     52      l_vremap              = ln_vert_remap 
    5253      tabspongedone_trn     = .FALSE. 
     54      ! 
    5355      CALL Agrif_Bc_Variable( trn_sponge_id, calledweight=zcoef, procname=interptrn_sponge ) 
     56      ! 
    5457      Agrif_UseSpecialValue = .FALSE. 
     58      l_vremap              = .FALSE. 
    5559#endif 
    5660      ! 
     
    5862 
    5963 
    60    SUBROUTINE interptrn_sponge( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) 
    61       !!---------------------------------------------------------------------- 
    62       !!                   *** ROUTINE interptrn_sponge *** 
     64   SUBROUTINE interptrn_sponge( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before)  
     65      !!---------------------------------------------------------------------- 
     66      !!                 *** ROUTINE interptrn_sponge *** 
    6367      !!---------------------------------------------------------------------- 
    6468      INTEGER                                     , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2 
     
    6872      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
    6973      INTEGER  ::   iku, ikv 
    70       REAL(wp) ::   ztra, zabe1, zabe2, zbtr 
    71       REAL(wp), DIMENSION(i1:i2,j1:j2,jpk)  ::  ztu, ztv 
    72       REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2) ::   trbdiff 
     74      REAL(wp) :: ztra, zabe1, zabe2, zbtr, zhtot 
     75      REAL(wp), DIMENSION(i1-1:i2,j1-1:j2,jpk) :: ztu, ztv 
     76      REAL(wp), DIMENSION(i1:i2,j1:j2,jpk,n1:n2) ::trbdiff 
    7377      ! vertical interpolation: 
    7478      REAL(wp), DIMENSION(i1:i2,j1:j2,jpk,n1:n2) ::tabres_child 
    75       REAL(wp), DIMENSION(k1:k2,1:jptra) :: tabin 
    76       REAL(wp), DIMENSION(k1:k2) :: h_in 
    77       REAL(wp), DIMENSION(1:jpk) :: h_out 
     79      REAL(wp), DIMENSION(k1:k2,n1:n2-1) :: tabin, tabin_i 
     80      REAL(wp), DIMENSION(k1:k2) :: z_in, z_in_i, h_in_i 
     81      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out 
    7882      INTEGER :: N_in, N_out 
    7983      !!---------------------------------------------------------------------- 
     
    9094         END DO 
    9195 
    92 # if defined key_vertical 
    93          DO jk=k1,k2 
    94             DO jj=j1,j2 
    95                DO ji=i1,i2 
    96                   tabres(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kbb_a) 
    97                END DO 
    98             END DO 
    99          END DO 
    100 # endif 
    101       ELSE       
    102 # if defined key_vertical 
    103          tabres_child(:,:,:,:) = 0. 
    104          DO jj=j1,j2 
    105             DO ji=i1,i2 
    106                N_in = 0 
    107                DO jk=k1,k2 !k2 = jpk of parent grid 
    108                   IF (tabres(ji,jj,jk,n2) == 0) EXIT 
    109                   N_in = N_in + 1 
    110                   tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1) 
    111                   h_in(N_in) = tabres(ji,jj,jk,n2) 
    112                END DO 
    113                N_out = 0 
    114                DO jk=1,jpk ! jpk of child grid 
    115                   IF (tmask(ji,jj,jk) == 0) EXIT  
    116                   N_out = N_out + 1 
    117                   h_out(jk) = e3t(ji,jj,jk,Kbb_a) !Child grid scale factors. Could multiply by e1e2t here instead of division above 
    118                ENDDO 
    119                IF (N_in > 0) THEN 
    120                   CALL reconstructandremap(tabin(1:N_in,1:jptra),h_in,tabres_child(ji,jj,1:N_out,1:jptra),h_out,N_in,N_out,jptra) 
    121                ENDIF 
    122             ENDDO 
    123          ENDDO 
    124 # endif 
    125  
    126          DO jj=j1,j2 
    127             DO ji=i1,i2 
    128                DO jk=1,jpkm1 
    129 # if defined key_vertical 
    130                   trbdiff(ji,jj,jk,1:jptra) = ( tr(ji,jj,jk,1:jptra,Kbb_a) - tabres_child(ji,jj,jk,1:jptra) ) * tmask(ji,jj,jk) 
    131 # else 
    132                   trbdiff(ji,jj,jk,1:jptra) = ( tr(ji,jj,jk,1:jptra,Kbb_a) - tabres(ji,jj,jk,1:jptra)) * tmask(ji,jj,jk) 
    133 # endif 
    134                ENDDO 
    135             ENDDO 
    136          ENDDO 
    137  
    138          DO jn = 1, jptra 
     96         IF ( l_vremap.OR.ln_zps ) THEN 
     97 
     98            ! Fill cell depths (i.e. gdept) to be interpolated 
     99            ! Warning: these are masked, hence extrapolated prior interpolation. 
     100            DO jj=j1,j2 
     101               DO ji=i1,i2 
     102                  tabres(ji,jj,k1,jptra+1) = 0.5_wp * tmask(ji,jj,k1) * e3t(ji,jj,k1,Kbb_a) 
     103                  DO jk=k1+1,k2 
     104                     tabres(ji,jj,jk,jptra+1) = tmask(ji,jj,jk) * & 
     105                        & ( tabres(ji,jj,jk-1,jptra+1) + 0.5_wp * (e3t(ji,jj,jk-1,Kbb_a)+e3t(ji,jj,jk,Kbb_a)) ) 
     106                  END DO 
     107               END DO 
     108            END DO 
     109 
     110            ! Save ssh at last level: 
     111            IF ( .NOT.ln_linssh ) THEN 
     112               tabres(i1:i2,j1:j2,k2,jptra+1) = ssh(i1:i2,j1:j2,Kbb_a)*tmask(i1:i2,j1:j2,1)  
     113            END IF   
     114     
     115         END IF 
     116 
     117      ELSE  
     118         ! 
     119         IF ( l_vremap ) THEN 
     120 
     121            IF (ln_linssh) tabres(i1:i2,j1:j2,k2,n2) = 0._wp 
     122 
     123            DO jj=j1,j2 
     124               DO ji=i1,i2 
     125 
     126                  tabres_child(ji,jj,:,:) = 0._wp  
     127                  ! Build vertical grids: 
     128                  N_in = mbkt_parent(ji,jj) 
     129                  ! Input grid (account for partial cells if any): 
     130                  DO jk=1,N_in 
     131                     z_in(jk) = tabres(ji,jj,jk,n2) - tabres(ji,jj,k2,n2) 
     132                     tabin(jk,1:jptra) = tabres(ji,jj,jk,1:jptra) 
     133                  END DO 
     134                   
     135                  ! Intermediate grid: 
     136                  DO jk = 1, N_in 
     137                     h_in_i(jk) = e3t0_parent(ji,jj,jk) * &  
     138                       &       (1._wp + tabres(ji,jj,k2,n2)/(ht0_parent(ji,jj)*ssmask(ji,jj) + 1._wp - ssmask(ji,jj))) 
     139                  END DO 
     140                  z_in_i(1) = 0.5_wp * h_in_i(1) 
     141                  DO jk=2,N_in 
     142                     z_in_i(jk) = z_in_i(jk-1) + 0.5_wp * ( h_in_i(jk) + h_in_i(jk-1) ) 
     143                  END DO 
     144                  z_in_i(1:N_in) = z_in_i(1:N_in)  - tabres(ji,jj,k2,n2) 
     145 
     146                  ! Output (Child) grid: 
     147                  N_out = mbkt(ji,jj) 
     148                  DO jk=1,N_out 
     149                     h_out(jk) = e3t(ji,jj,jk,Kbb_a) 
     150                  END DO 
     151                  z_out(1) = 0.5_wp * h_out(1) 
     152                  DO jk=2,N_out 
     153                     z_out(jk) = z_out(jk-1) + 0.5_wp * ( h_out(jk)+h_out(jk-1) ) 
     154                  END DO 
     155                  IF (.NOT.ln_linssh) z_out(1:N_out) = z_out(1:N_out)  - ssh(ji,jj,Kbb_a) 
     156 
     157                  ! Account for small differences in the free-surface 
     158                  IF ( sum(h_out(1:N_out)) > sum(h_in_i(1:N_in) )) THEN 
     159                     h_out(1) = h_out(1)  - ( sum(h_out(1:N_out))-sum(h_in_i(1:N_in)) ) 
     160                  ELSE 
     161                     h_in_i(1)= h_in_i(1) - ( sum(h_in_i(1:N_in))-sum(h_out(1:N_out)) ) 
     162                  END IF 
     163                  IF (N_in*N_out > 0) THEN 
     164                     CALL remap_linear(tabin(1:N_in,1:jptra),z_in(1:N_in),tabin_i(1:N_in,1:jptra),z_in_i(1:N_in),N_in,N_in,jptra) 
     165                     CALL reconstructandremap(tabin_i(1:N_in,1:jptra),h_in_i(1:N_in),tabres_child(ji,jj,1:N_out,1:jptra),h_out(1:N_out),N_in,N_out,jptra) 
     166!                     CALL remap_linear(tabin(1:N_in,1:jptra),z_in(1:N_in),tabres_child(ji,jj,1:N_out,1:jptra),z_out(1:N_in),N_in,N_out,jptra)   
     167                  ENDIF 
     168               END DO 
     169            END DO 
     170 
     171            DO jj=j1,j2 
     172               DO ji=i1,i2 
     173                  DO jk=1,jpkm1 
     174                     trbdiff(ji,jj,jk,1:jptra) = (tr(ji,jj,jk,1:jptra,Kbb_a) - tabres_child(ji,jj,jk,1:jptra)) * tmask(ji,jj,jk) 
     175                  END DO 
     176               END DO 
     177            END DO 
     178 
     179         ELSE 
     180 
     181            IF ( Agrif_Parent(ln_zps) ) THEN ! Account for partial cells 
     182 
     183               DO jj=j1,j2 
     184                  DO ji=i1,i2 
     185                     ! 
     186                     N_in  = mbkt(ji,jj)  
     187                     N_out = mbkt(ji,jj)  
     188                     z_in(1) = tabres(ji,jj,1,n2) 
     189                     tabin(1,1:jptra) = tabres(ji,jj,1,1:jptra) 
     190                     DO jk=2, N_in 
     191                        z_in(jk) = tabres(ji,jj,jk,n2) 
     192                        tabin(jk,1:jptra) = tabres(ji,jj,jk,1:jptra) 
     193                     END DO  
     194                     IF (.NOT.ln_linssh) z_in(1:N_in) = z_in(1:N_in) - tabres(ji,jj,k2,n2) 
     195 
     196                     z_out(1) = 0.5_wp * e3t(ji,jj,1,Kbb_a) 
     197                     DO jk=2, N_out 
     198                        z_out(jk) = z_out(jk-1) + 0.5_wp * (e3t(ji,jj,jk-1,Kbb_a) + e3t(ji,jj,jk,Kbb_a))  
     199                     END DO  
     200                     IF (.NOT.ln_linssh) z_out(1:N_out) = z_out(1:N_out) - ssh(ji,jj,Kbb_a) 
     201 
     202                     CALL remap_linear(tabin(1:N_in,1:jptra), z_in(1:N_in), tabres(ji,jj,1:N_out,1:jptra), & 
     203                                         &   z_out(1:N_out), N_in, N_out, jptra) 
     204                  END DO 
     205               END DO 
     206            ENDIF 
     207 
     208            DO jj=j1,j2 
     209               DO ji=i1,i2 
     210                  DO jk=1,jpkm1 
     211                     trbdiff(ji,jj,jk,1:jptra) = (tr(ji,jj,jk,1:jptra,Kbb_a) - tabres(ji,jj,jk,1:jptra))*tmask(ji,jj,jk) 
     212                  END DO 
     213               END DO 
     214            END DO 
     215 
     216         END IF 
     217 
     218         DO jn = 1, jptra             
    139219            DO jk = 1, jpkm1 
    140                ztu(i1:i2,j1:j2,jk) = 0._wp 
     220               ztu(i1-1:i2,j1-1:j2,jk) = 0._wp 
    141221               DO jj = j1,j2 
    142222                  DO ji = i1,i2-1 
    143                      zabe1 = rn_sponge_tra * r1_Dt * fspu(ji,jj) * umask(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm_a) 
    144                      ztu(ji,jj,jk) = zabe1 * ( trbdiff(ji+1,jj,jk,jn) - trbdiff(ji,jj,jk,jn) ) 
    145                   END DO 
    146                END DO 
    147                ztv(i1:i2,j1:j2,jk) = 0._wp 
     223                     zabe1 = rn_sponge_tra * r1_Dt * umask(ji,jj,jk) * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) 
     224                     ztu(ji,jj,jk) = zabe1 * fspu(ji,jj) * ( trbdiff(ji+1,jj  ,jk,jn) - trbdiff(ji,jj,jk,jn) )  
     225                  END DO 
     226               END DO 
     227               ztv(i1-1:i2,j1-1:j2,jk) = 0._wp 
    148228               DO ji = i1,i2 
    149229                  DO jj = j1,j2-1 
    150                      zabe2 = rn_sponge_tra * r1_Dt * fspv(ji,jj) * vmask(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm_a) 
    151                      ztv(ji,jj,jk) = zabe2 * ( trbdiff(ji,jj+1,jk,jn) - trbdiff(ji,jj,jk,jn) ) 
     230                     zabe2 = rn_sponge_tra * r1_Dt * vmask(ji,jj,jk) * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm_a) 
     231                     ztv(ji,jj,jk) = zabe2 * fspv(ji,jj) * ( trbdiff(ji  ,jj+1,jk,jn) - trbdiff(ji,jj,jk,jn) ) 
    152232                  END DO 
    153233               END DO 
     
    166246            END DO 
    167247            ! 
     248! JC: there is something wrong with the Laplacian in corners 
    168249            DO jk = 1, jpkm1 
    169                DO jj = j1+1,j2-1 
    170                   DO ji = i1+1,i2-1 
    171                      IF( .NOT. tabspongedone_trn(ji,jj) ) THEN 
     250               DO jj = j1,j2 
     251                  DO ji = i1,i2 
     252                     IF (.NOT. tabspongedone_trn(ji,jj)) THEN  
    172253                        zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm_a) 
    173                                                 ! horizontal diffusive trends 
    174                         ztra = zbtr * (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji,jj-1,jk)  ) & 
    175                              &  - rn_trelax_tra * r1_Dt * fspt(ji,jj) * trbdiff(ji,jj,jk,jn) 
     254                        ! horizontal diffusive trends 
     255                        ztra = zbtr * ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk)   &  
     256                          &           + ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) & 
     257                          &   - rn_trelax_tra * r1_Dt * fspt(ji,jj) * trbdiff(ji,jj,jk,jn) 
    176258                        ! add it to the general tracer trends 
    177259                        tr(ji,jj,jk,jn,Krhs_a) = tr(ji,jj,jk,jn,Krhs_a) + ztra 
     
    179261                  END DO 
    180262               END DO 
    181             END DO 
     263 
     264            END DO 
     265            ! 
    182266         END DO 
    183267         ! 
    184          tabspongedone_trn(i1+1:i2-1,j1+1:j2-1) = .TRUE. 
     268         tabspongedone_trn(i1:i2,j1:j2) = .TRUE. 
     269         ! 
    185270      ENDIF 
    186       !                  
     271      ! 
    187272   END SUBROUTINE interptrn_sponge 
    188273 
Note: See TracChangeset for help on using the changeset viewer.