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 10331 for NEMO/branches – NEMO

Changeset 10331 for NEMO/branches


Ignore:
Timestamp:
2018-11-19T10:50:43+01:00 (5 years ago)
Author:
clem
Message:

this is the best shot we have for now concerning this very tricky issue of crazy ice thicknesses after advection

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_adv_umx.F90

    r10315 r10331  
    3737   REAL(wp), DIMENSION(:,:), ALLOCATABLE :: z1_a, z1_a_ups, zua_ups, zva_ups 
    3838    
     39   ! rachid trick (only for nonosc) 
     40   LOGICAL :: ll_rachid = .TRUE. 
     41 
    3942   ! alternate directions for upstream 
    4043   ! clem: it gives results for Lipscomb test that are the same as "ll_upsxy=false" 
     
    4447   ! prelimiter 
    4548   ! clem: use it to avoid overshoot in H 
    46    LOGICAL :: ll_prelimiter = .TRUE. 
    4749   LOGICAL :: ll_prelimiter_zalesak = .TRUE.  ! from: Zalesak(1979) eq. 14 => better than Devore especially in 2D 
    48    LOGICAL :: ll_prelimiter_devore  = .FALSE. ! from: Devore eq. 11 
     50   LOGICAL :: ll_prelimiter_devore  = .FALSE.  ! from: Devore eq. 11 
    4951 
    5052   ! iterate on the limiter (only nonosc) 
     
    154156 
    155157            ! 1/A_ups 
    156             WHERE( za_ups(:,:) > epsi20 )   ;   z1_a_ups(:,:) = 1._wp / za_ups(:,:) 
    157             ELSEWHERE                       ;   z1_a_ups(:,:) = 0. 
    158             END WHERE 
    159  
    160             !!zua_ho = zudy ; zva_ho = zvdx 
    161             !!CALL adv_umx( kn_umx, kt, zdt, zudy, zvdx, zua_ho , zva_ho , zcu_box, zcv_box, pv_i(:,:,jl), pv_i(:,:,jl) ) 
     158!!            IF( .NOT. ll_rachid ) THEN 
     159               WHERE( za_ups(:,:) > epsi20 )   ;   z1_a_ups(:,:) = 1._wp / za_ups(:,:) 
     160               ELSEWHERE                       ;   z1_a_ups(:,:) = 0. 
     161               END WHERE 
     162!!            ELSE 
     163!!               WHERE( pa_i(:,:,jl) > epsi20 )   ;   z1_a_ups(:,:) = 1._wp / pa_i(:,:,jl) 
     164!!               ELSEWHERE                        ;   z1_a_ups(:,:) = 0. 
     165!!               END WHERE               
     166!!            ENDIF 
     167!! 
     168!!            IF( ll_rachid )   zua_ups = zua_ho ; zva_ups = zva_ho 
    162169            zamsk = 0._wp ; zh_i(:,:) = pv_i (:,:,jl) * z1_a(:,:) 
    163170            CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy, zvdx, zua_ho , zva_ho , zcu_box, zcv_box, zh_i(:,:), pv_i (:,:,jl), zua_ups, zva_ups )              ! Ice volume 
     
    193200               END WHERE 
    194201 
     202!!            IF( ll_rachid )   zua_ups = zua_ho ; zva_ups = zva_ho 
    195203               zamsk = 0._wp ; zh_ip(:,:) = pv_ip(:,:,jl) * z1_a(:,:) 
    196204               CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy, zvdx, zua_ho , zva_ho , zcu_box, zcv_box, zh_ip(:,:), pv_ip(:,:,jl), zua_ups, zva_ups )              ! mp volume 
     
    256264      !  upstream (_ups) advection with initial mass fluxes 
    257265      ! --------------------------------------------------- 
     266      IF( ll_rachid )    zfu_ups=0.; zfv_ups=0. 
    258267      IF( .NOT. ll_upsxy ) THEN 
    259268 
     
    320329          
    321330      ENDIF 
     331 
     332      ! Rachid trick 
     333      ! ------------ 
     334      IF( ll_rachid .AND. kn_limiter /= 1 )  & 
     335         & CALL ctl_stop( 'STOP', 'icedyn_adv_umx: ll_rachid incompatible with limiters other than nonosc' ) 
     336      IF( ll_rachid ) THEN 
     337         call lbc_lnk (zfu_ups,'U',-1.) 
     338         call lbc_lnk (zfv_ups,'V',-1.) 
     339      IF( pamsk == 0. ) THEN 
     340         WHERE( ABS( puc(:,:) ) > 0._wp )   ;   zfu_ups(:,:) = zfu_ups(:,:) / puc(:,:) * pu(:,:) 
     341         ELSEWHERE                          ;   zfu_ups(:,:) = 0._wp 
     342         END WHERE 
     343 
     344         WHERE( ABS( pvc(:,:) ) > 0._wp )   ;   zfv_ups(:,:) = zfv_ups(:,:) / pvc(:,:) * pv(:,:) 
     345         ELSEWHERE                          ;   zfv_ups(:,:) = 0._wp 
     346         END WHERE          
     347      ENDIF 
     348      ENDIF 
     349 
    322350      ! guess after content field with upstream scheme 
    323351      DO jj = 2, jpjm1 
     
    325353            ztra          = - (   zfu_ups(ji,jj) - zfu_ups(ji-1,jj  )   & 
    326354               &                + zfv_ups(ji,jj) - zfv_ups(ji  ,jj-1) ) * r1_e1e2t(ji,jj) 
    327             zt_ups(ji,jj) = ( ptc(ji,jj) + pdt * ztra ) * tmask(ji,jj,1) 
     355            IF( ll_rachid ) THEN 
     356               zt_ups(ji,jj) = ( pt (ji,jj) + pdt * ztra ) * tmask(ji,jj,1) 
     357            ELSE 
     358               zt_ups(ji,jj) = ( ptc(ji,jj) + pdt * ztra ) * tmask(ji,jj,1) 
     359            ENDIF 
    328360         END DO 
    329361      END DO 
     
    336368      CASE ( 20 )                          !== centered second order ==! 
    337369         ! 
    338          CALL cen2( kn_limiter, jt, kt, pdt, pt, pu, pv, puc, pvc, ptc, zfu_ho, zfv_ho,  & 
     370         CALL cen2( pamsk, kn_limiter, jt, kt, pdt, pt, pu, pv, puc, pvc, ptc, zfu_ho, zfv_ho,  & 
    339371            &       zt_ups, zfu_ups, zfv_ups ) 
    340372         ! 
     
    377409   END SUBROUTINE adv_umx 
    378410 
    379    SUBROUTINE cen2( kn_limiter, jt, kt, pdt, pt, pu, pv, puc, pvc, ptc, pfu_ho, pfv_ho, & 
     411   SUBROUTINE cen2( pamsk, kn_limiter, jt, kt, pdt, pt, pu, pv, puc, pvc, ptc, pfu_ho, pfv_ho, & 
    380412      &             pt_ups, pfu_ups, pfv_ups ) 
    381413      !!--------------------------------------------------------------------- 
     
    389421      !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74.  
    390422      !!---------------------------------------------------------------------- 
     423      REAL(wp)                    , INTENT(in   ) ::   pamsk            ! advection of concentration (1) or other tracers (0) 
    391424      INTEGER                     , INTENT(in   ) ::   kn_limiter       ! limiter 
    392425      INTEGER                     , INTENT(in   ) ::   jt               ! number of sub-iteration 
     
    399432      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) ::   pfu_ho, pfv_ho   ! high order fluxes  
    400433      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) ::   pt_ups           ! upstream guess of tracer content  
    401       REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) ::   pfu_ups, pfv_ups ! upstream fluxes  
     434      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pfu_ups, pfv_ups ! upstream fluxes  
    402435      ! 
    403436      INTEGER  ::   ji, jj    ! dummy loop indices 
     
    415448         END DO 
    416449         IF    ( kn_limiter == 1 ) THEN 
    417             CALL nonosc_2d( pdt, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 
     450            CALL nonosc_2d( pamsk, pdt, pu, puc, pv, pvc, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 
    418451         ELSEIF( kn_limiter == 2 ) THEN 
    419452            CALL limiter_x( pdt, pu, puc, pt, pfu_ho ) 
     
    483516 
    484517         ENDIF 
    485          IF( kn_limiter == 1 )   CALL nonosc_2d( pdt, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 
     518         IF( kn_limiter == 1 )   CALL nonosc_2d( pamsk, pdt, pu, puc, pv, pvc, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 
    486519          
    487520      ENDIF 
     
    516549      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) ::   pfu_ho, pfv_ho   ! high order fluxes  
    517550      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) ::   pt_ups           ! upstream guess of tracer content  
    518       REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) ::   pfu_ups, pfv_ups ! upstream fluxes  
     551      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pfu_ups, pfv_ups ! upstream fluxes  
    519552      ! 
    520553      INTEGER  ::   ji, jj    ! dummy loop indices 
     
    570603      IF( kn_limiter == 1 ) THEN 
    571604         IF( .NOT. ll_limiter_it2 ) THEN 
    572             CALL nonosc_2d ( pdt, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 
     605            CALL nonosc_2d ( pamsk, pdt, pu, puc, pv, pvc, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 
    573606         ELSE 
    574607            zzfu_ho(:,:) = pfu_ho(:,:) 
    575608            zzfv_ho(:,:) = pfv_ho(:,:) 
    576609            ! 1st iteration of nonosc (limit the flux with the upstream solution) 
    577             CALL nonosc_2d ( pdt, ptc, pt_ups, pfu_ups, pfv_ups, zzfu_ho, zzfv_ho ) 
     610            CALL nonosc_2d ( pamsk, pdt, pu, puc, pv, pvc, ptc, pt_ups, pfu_ups, pfv_ups, zzfu_ho, zzfv_ho ) 
    578611            ! guess after content field with high order 
    579612            DO jj = 2, jpjm1 
     
    585618            CALL lbc_lnk( zzt, 'T', 1. ) 
    586619            ! 2nd iteration of nonosc (limit the flux with the limited high order solution) 
    587             CALL nonosc_2d ( pdt, ptc, zzt, zzfu_ho, zzfv_ho, pfu_ho, pfv_ho ) 
     620            CALL nonosc_2d ( pamsk, pdt, pu, puc, pv, pvc, ptc, zzt, zzfu_ho, zzfv_ho, pfu_ho, pfv_ho ) 
    588621         ENDIF 
    589622      ENDIF 
     
    842875      
    843876 
    844    SUBROUTINE nonosc_2d( pdt, ptc, pt_low, pfu_low, pfv_low, pfu_ho, pfv_ho ) 
     877   SUBROUTINE nonosc_2d( pamsk, pdt, pu, puc, pv, pvc, ptc, pt_low, pfu_low, pfv_low, pfu_ho, pfv_ho ) 
    845878      !!--------------------------------------------------------------------- 
    846879      !!                    ***  ROUTINE nonosc  *** 
     
    855888      !!       in-space based differencing for fluid 
    856889      !!---------------------------------------------------------------------- 
     890      REAL(wp)                     , INTENT(in   ) ::   pamsk            ! advection of concentration (1) or other tracers (0) 
    857891      REAL(wp)                     , INTENT(in   ) ::   pdt              ! tracer time-step 
     892      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pu               ! ice i-velocity => u*e2 
     893      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   puc              ! ice i-velocity *A => u*e2*a 
     894      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pv               ! ice j-velocity => v*e1 
     895      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pvc              ! ice j-velocity *A => v*e1*a 
    858896      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   ptc, pt_low      ! before field & upstream guess of after field 
    859       REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pfv_low, pfu_low ! upstream flux 
     897      REAL(wp), DIMENSION (jpi,jpj), INTENT(inout) ::   pfv_low, pfu_low ! upstream flux 
    860898      REAL(wp), DIMENSION (jpi,jpj), INTENT(inout) ::   pfv_ho, pfu_ho   ! monotonic flux 
    861899      ! 
     
    868906      zsml = epsi20 
    869907 
     908      ! Rachid trick (upstream already done in macho) 
     909      ! ------------ 
     910      IF( ll_rachid ) THEN 
     911      IF( pamsk == 0. ) THEN 
     912         WHERE( ABS( puc(:,:) ) > 0._wp )   ;   pfu_ho(:,:) = pfu_ho(:,:) / puc(:,:) * pu(:,:) 
     913         ELSEWHERE                          ;   pfu_ho(:,:) = 0._wp 
     914         END WHERE 
     915 
     916         WHERE( ABS( pvc(:,:) ) > 0._wp )   ;   pfv_ho(:,:) = pfv_ho(:,:) / pvc(:,:) * pv(:,:) 
     917         ELSEWHERE                          ;   pfv_ho(:,:) = 0._wp 
     918         END WHERE          
     919      ENDIF 
     920      ENDIF 
     921       
    870922      ! antidiffusive flux : high order minus low order 
    871923      ! -------------------------------------------------- 
     
    885937      !                        |      |      |        |    * 
    886938      !            t_low :       i-1     i       i+1       i+2    
    887       IF( ll_prelimiter ) THEN 
    888  
     939      IF( ll_prelimiter_zalesak ) THEN 
     940          
    889941         DO jj = 2, jpjm1 
    890942            DO ji = fs_2, fs_jpim1  
     
    895947         CALL lbc_lnk_multi( zti_low, 'T', 1., ztj_low, 'T', 1. ) 
    896948 
    897          IF( ll_prelimiter_zalesak ) THEN 
    898  
    899             !! this does not work 
    900 !!            DO jj = 2, jpjm1 
    901 !!               DO ji = fs_2, fs_jpim1 
    902 !!                  IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., pt_low (ji+1,jj  ) - pt_low (ji  ,jj) ) .AND.     & 
    903 !!               &      SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., pt_low (ji  ,jj+1) - pt_low (ji  ,jj) )           & 
    904 !!               &    ) THEN 
    905 !!                     IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., zti_low(ji+1,jj ) - zti_low(ji  ,jj) ) .AND.   & 
    906 !!               &         SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., ztj_low(ji,jj+1 ) - ztj_low(ji  ,jj) )         & 
    907 !!               &       ) THEN 
    908 !!                        pfu_ho(ji,jj) = 0.  ;   pfv_ho(ji,jj) = 0. 
    909 !!                     ENDIF 
    910 !!                     IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., pt_low (ji  ,jj) - pt_low (ji-1,jj  ) ) .AND.  & 
    911 !!               &         SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., pt_low (ji  ,jj) - pt_low (ji  ,jj-1) )        & 
    912 !!               &       )  THEN 
    913 !!                        pfu_ho(ji,jj) = 0.  ;   pfv_ho(ji,jj) = 0. 
    914 !!                     ENDIF 
    915 !!                  ENDIF 
    916 !!                END DO 
    917 !!            END DO             
    918  
    919             DO jj = 2, jpjm1 
    920                DO ji = fs_2, fs_jpim1 
    921                   IF ( pfu_ho(ji,jj) * ( pt_low(ji+1,jj) - pt_low(ji,jj) ) <= 0. .AND.  & 
    922                      & pfv_ho(ji,jj) * ( pt_low(ji,jj+1) - pt_low(ji,jj) ) <= 0. ) THEN 
    923                      ! 
    924                      IF(  pfu_ho(ji,jj) * ( zti_low(ji+1,jj) - zti_low(ji,jj) ) <= 0 .AND.  & 
    925                         & pfv_ho(ji,jj) * ( ztj_low(ji,jj+1) - ztj_low(ji,jj) ) <= 0)  pfu_ho(ji,jj)=0. ; pfv_ho(ji,jj)=0. 
    926                      ! 
    927                      IF(  pfu_ho(ji,jj) * ( pt_low(ji  ,jj) - pt_low(ji-1,jj) ) <= 0 .AND.  & 
    928                         & pfv_ho(ji,jj) * ( pt_low(ji  ,jj) - pt_low(ji,jj-1) ) <= 0)  pfu_ho(ji,jj)=0. ; pfv_ho(ji,jj)=0. 
    929                      ! 
    930                   ENDIF 
    931                END DO 
    932             END DO 
    933             CALL lbc_lnk_multi( pfu_ho, 'U', -1., pfv_ho, 'V', -1. )   ! lateral boundary cond. 
    934  
    935          ELSEIF( ll_prelimiter_devore ) THEN 
    936             z1_dt = 1._wp / pdt 
    937             DO jj = 2, jpjm1 
    938                DO ji = fs_2, fs_jpim1 
    939                   zsign = SIGN( 1., pt_low(ji+1,jj) - pt_low(ji,jj) ) 
    940                   pfu_ho(ji,jj) =  zsign * MAX( 0. , MIN( ABS(pfu_ho(ji,jj)) , & 
    941                      &                          zsign * ( pt_low (ji  ,jj) - pt_low (ji-1,jj) ) * e1e2t(ji  ,jj) * z1_dt , & 
    942                      &                          zsign * ( zti_low(ji+1,jj) - zti_low(ji  ,jj) ) * e1e2t(ji+1,jj) * z1_dt ) ) 
    943                    
    944                   zsign = SIGN( 1., pt_low(ji,jj+1) - pt_low(ji,jj) ) 
    945                   pfv_ho(ji,jj) =  zsign * MAX( 0. , MIN( ABS(pfv_ho(ji,jj)) , & 
    946                      &                          zsign * ( pt_low (ji,jj  ) - pt_low (ji,jj-1) ) * e1e2t(ji,jj  ) * z1_dt , & 
    947                      &                          zsign * ( ztj_low(ji,jj+1) - ztj_low(ji,jj  ) ) * e1e2t(ji,jj+1) * z1_dt ) ) 
    948                END DO 
    949             END DO 
    950             CALL lbc_lnk_multi( pfu_ho, 'U', -1., pfv_ho, 'V', -1. )   ! lateral boundary cond. 
     949 
     950         !! this does not work 
     951         !!            DO jj = 2, jpjm1 
     952         !!               DO ji = fs_2, fs_jpim1 
     953         !!                  IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., pt_low (ji+1,jj  ) - pt_low (ji  ,jj) ) .AND.     & 
     954         !!               &      SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., pt_low (ji  ,jj+1) - pt_low (ji  ,jj) )           & 
     955         !!               &    ) THEN 
     956         !!                     IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., zti_low(ji+1,jj ) - zti_low(ji  ,jj) ) .AND.   & 
     957         !!               &         SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., ztj_low(ji,jj+1 ) - ztj_low(ji  ,jj) )         & 
     958         !!               &       ) THEN 
     959         !!                        pfu_ho(ji,jj) = 0.  ;   pfv_ho(ji,jj) = 0. 
     960         !!                     ENDIF 
     961         !!                     IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., pt_low (ji  ,jj) - pt_low (ji-1,jj  ) ) .AND.  & 
     962         !!               &         SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., pt_low (ji  ,jj) - pt_low (ji  ,jj-1) )        & 
     963         !!               &       )  THEN 
     964         !!                        pfu_ho(ji,jj) = 0.  ;   pfv_ho(ji,jj) = 0. 
     965         !!                     ENDIF 
     966         !!                  ENDIF 
     967         !!                END DO 
     968         !!            END DO             
     969 
     970         DO jj = 2, jpjm1 
     971            DO ji = fs_2, fs_jpim1 
     972               IF ( pfu_ho(ji,jj) * ( pt_low(ji+1,jj) - pt_low(ji,jj) ) <= 0. .AND.  & 
     973                  & pfv_ho(ji,jj) * ( pt_low(ji,jj+1) - pt_low(ji,jj) ) <= 0. ) THEN 
     974                  ! 
     975                  IF(  pfu_ho(ji,jj) * ( zti_low(ji+1,jj) - zti_low(ji,jj) ) <= 0 .AND.  & 
     976                     & pfv_ho(ji,jj) * ( ztj_low(ji,jj+1) - ztj_low(ji,jj) ) <= 0)  pfu_ho(ji,jj)=0. ; pfv_ho(ji,jj)=0. 
     977                  ! 
     978                  IF(  pfu_ho(ji,jj) * ( pt_low(ji  ,jj) - pt_low(ji-1,jj) ) <= 0 .AND.  & 
     979                     & pfv_ho(ji,jj) * ( pt_low(ji  ,jj) - pt_low(ji,jj-1) ) <= 0)  pfu_ho(ji,jj)=0. ; pfv_ho(ji,jj)=0. 
     980                  ! 
     981               ENDIF 
     982            END DO 
     983         END DO 
     984         CALL lbc_lnk_multi( pfu_ho, 'U', -1., pfv_ho, 'V', -1. )   ! lateral boundary cond. 
     985 
     986      ELSEIF( ll_prelimiter_devore ) THEN 
     987         DO jj = 2, jpjm1 
     988            DO ji = fs_2, fs_jpim1  
     989               zti_low(ji,jj)= pt_low(ji+1,jj  ) 
     990               ztj_low(ji,jj)= pt_low(ji  ,jj+1) 
     991            END DO 
     992         END DO 
     993         CALL lbc_lnk_multi( zti_low, 'T', 1., ztj_low, 'T', 1. ) 
     994 
     995         z1_dt = 1._wp / pdt 
     996         DO jj = 2, jpjm1 
     997            DO ji = fs_2, fs_jpim1 
     998               zsign = SIGN( 1., pt_low(ji+1,jj) - pt_low(ji,jj) ) 
     999               pfu_ho(ji,jj) =  zsign * MAX( 0. , MIN( ABS(pfu_ho(ji,jj)) , & 
     1000                  &                          zsign * ( pt_low (ji  ,jj) - pt_low (ji-1,jj) ) * e1e2t(ji  ,jj) * z1_dt , & 
     1001                  &                          zsign * ( zti_low(ji+1,jj) - zti_low(ji  ,jj) ) * e1e2t(ji+1,jj) * z1_dt ) ) 
     1002 
     1003               zsign = SIGN( 1., pt_low(ji,jj+1) - pt_low(ji,jj) ) 
     1004               pfv_ho(ji,jj) =  zsign * MAX( 0. , MIN( ABS(pfv_ho(ji,jj)) , & 
     1005                  &                          zsign * ( pt_low (ji,jj  ) - pt_low (ji,jj-1) ) * e1e2t(ji,jj  ) * z1_dt , & 
     1006                  &                          zsign * ( ztj_low(ji,jj+1) - ztj_low(ji,jj  ) ) * e1e2t(ji,jj+1) * z1_dt ) ) 
     1007            END DO 
     1008         END DO 
     1009         CALL lbc_lnk_multi( pfu_ho, 'U', -1., pfv_ho, 'V', -1. )   ! lateral boundary cond. 
    9511010             
    952          ENDIF 
     1011      ENDIF 
    9531012          
    954       ENDIF 
    9551013       
    9561014      ! Search local extrema 
     
    9821040         DO ji = fs_2, fs_jpim1   ! vector opt. 
    9831041            ! 
    984             !zup  = MAX( zbup(ji,jj), zbup(ji-1,jj  ), zbup(ji+1,jj  ), zbup(ji  ,jj-1), zbup(ji  ,jj+1) )  ! search max/min in neighbourhood 
    985             !zdo  = MIN( zbdo(ji,jj), zbdo(ji-1,jj  ), zbdo(ji+1,jj  ), zbdo(ji  ,jj-1), zbdo(ji  ,jj+1) ) 
     1042!            zup  = MAX( zbup(ji,jj), zbup(ji-1,jj  ), zbup(ji+1,jj  ), zbup(ji  ,jj-1), zbup(ji  ,jj+1) )  ! search max/min in neighbourhood 
     1043!            zdo  = MIN( zbdo(ji,jj), zbdo(ji-1,jj  ), zbdo(ji+1,jj  ), zbdo(ji  ,jj-1), zbdo(ji  ,jj+1) ) 
    9861044            ! 
    9871045            zup  = MAX( zbup(ji,jj), zbup(ji-1,jj  ), zbup(ji+1,jj  ), zbup(ji  ,jj-1), zbup(ji  ,jj+1), &  ! search max/min in neighbourhood 
     
    10141072      CALL lbc_lnk_multi( zbetup, 'T', 1., zbetdo, 'T', 1. )   ! lateral boundary cond. (unchanged sign) 
    10151073 
     1074      ! Rachid trick 
     1075      ! ------------ 
     1076      IF( ll_rachid ) THEN 
     1077      IF( pamsk == 0. ) THEN 
     1078         WHERE( ABS( pu(:,:) ) > 0._wp ) 
     1079            pfu_ho (:,:) = pfu_ho (:,:) * puc(:,:) / pu(:,:) 
     1080            pfu_low(:,:) = pfu_low(:,:) * puc(:,:) / pu(:,:) 
     1081         ELSEWHERE 
     1082            pfu_ho (:,:) = 0._wp 
     1083            pfu_low(:,:) = 0._wp 
     1084         END WHERE 
     1085 
     1086         WHERE( ABS( pv(:,:) ) > 0._wp ) 
     1087            pfv_ho (:,:) = pfv_ho (:,:) * pvc(:,:) / pv(:,:) 
     1088            pfv_low(:,:) = pfv_low(:,:) * pvc(:,:) / pv(:,:) 
     1089         ELSEWHERE 
     1090            pfv_ho (:,:) = 0._wp 
     1091            pfv_low(:,:) = 0._wp 
     1092         END WHERE          
     1093      ENDIF 
     1094      ENDIF 
     1095       
    10161096      ! monotonic flux in the y direction 
    10171097      ! --------------------------------- 
Note: See TracChangeset for help on using the changeset viewer.