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 13613 for NEMO/trunk – NEMO

Changeset 13613 for NEMO/trunk


Ignore:
Timestamp:
2020-10-15T15:42:25+02:00 (4 years ago)
Author:
smasson
Message:

trunk: icedyn_adv_pra.F90 optimisation see #2552

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/ICE/icedyn_adv_pra.F90

    r13609 r13613  
    288288                  CALL adv_x( zdt , zudy , 0._wp , zarea , z0vl , sxvl , sxxvl , syvl , syyvl , sxyvl )  
    289289               ENDIF 
    290            ENDIF 
    291             ! 
     290            ENDIF 
     291            ! 
     292         ENDIF 
     293          
     294         ! --- Lateral boundary conditions --- ! 
     295         !     caution: for gradients (sx and sy) the sign changes 
     296         CALL lbc_lnk_multi( 'icedyn_adv_pra', z0ice , 'T', 1._wp, sxice , 'T', -1._wp, syice , 'T', -1._wp  & ! ice volume 
     297            &                                , sxxice, 'T', 1._wp, syyice, 'T',  1._wp, sxyice, 'T',  1._wp  & 
     298            &                                , z0snw , 'T', 1._wp, sxsn  , 'T', -1._wp, sysn  , 'T', -1._wp  & ! snw volume 
     299            &                                , sxxsn , 'T', 1._wp, syysn , 'T',  1._wp, sxysn , 'T',  1._wp  ) 
     300         CALL lbc_lnk_multi( 'icedyn_adv_pra', z0smi , 'T', 1._wp, sxsal , 'T', -1._wp, sysal , 'T', -1._wp  & ! ice salinity 
     301            &                                , sxxsal, 'T', 1._wp, syysal, 'T',  1._wp, sxysal, 'T',  1._wp  & 
     302            &                                , z0ai  , 'T', 1._wp, sxa   , 'T', -1._wp, sya   , 'T', -1._wp  & ! ice concentration 
     303            &                                , sxxa  , 'T', 1._wp, syya  , 'T',  1._wp, sxya  , 'T',  1._wp  ) 
     304         CALL lbc_lnk_multi( 'icedyn_adv_pra', z0oi  , 'T', 1._wp, sxage , 'T', -1._wp, syage , 'T', -1._wp  & ! ice age 
     305            &                                , sxxage, 'T', 1._wp, syyage, 'T',  1._wp, sxyage, 'T',  1._wp  ) 
     306         CALL lbc_lnk_multi( 'icedyn_adv_pra', z0es  , 'T', 1._wp, sxc0  , 'T', -1._wp, syc0  , 'T', -1._wp  & ! snw enthalpy 
     307            &                                , sxxc0 , 'T', 1._wp, syyc0 , 'T',  1._wp, sxyc0 , 'T',  1._wp  )  
     308         CALL lbc_lnk_multi( 'icedyn_adv_pra', z0ei  , 'T', 1._wp, sxe   , 'T', -1._wp, sye   , 'T', -1._wp  & ! ice enthalpy 
     309            &                                , sxxe  , 'T', 1._wp, syye  , 'T',  1._wp, sxye  , 'T',  1._wp  ) 
     310         IF ( ln_pnd_LEV ) THEN 
     311            CALL lbc_lnk_multi( 'icedyn_adv_pra', z0ap , 'T', 1._wp, sxap , 'T', -1._wp, syap , 'T', -1._wp  & ! melt pond fraction 
     312               &                                , sxxap, 'T', 1._wp, syyap, 'T',  1._wp, sxyap, 'T',  1._wp  & 
     313               &                                , z0vp , 'T', 1._wp, sxvp , 'T', -1._wp, syvp , 'T', -1._wp  & ! melt pond volume 
     314               &                                , sxxvp, 'T', 1._wp, syyvp, 'T',  1._wp, sxyvp, 'T',  1._wp  )  
     315            IF ( ln_pnd_lids ) THEN 
     316               CALL lbc_lnk_multi( 'icedyn_adv_pra', z0vl ,'T', 1._wp, sxvl ,'T', -1._wp, syvl ,'T', -1._wp  & ! melt pond lid volume 
     317                  &                                , sxxvl,'T', 1._wp, syyvl,'T',  1._wp, sxyvl,'T',  1._wp  )  
     318            ENDIF 
    292319         ENDIF 
    293320 
     
    368395      !!  
    369396      INTEGER  ::   ji, jj, jl, jcat                     ! dummy loop indices 
     397      INTEGER  ::   jj0                                  ! dummy loop indices 
    370398      REAL(wp) ::   zs1max, zslpmax, ztemp               ! local scalars 
    371399      REAL(wp) ::   zs1new, zalf , zalfq , zbt           !   -      - 
     
    375403      REAL(wp), DIMENSION(jpi,jpj) ::   zalg, zalg1, zalg1q         !  -      - 
    376404      !----------------------------------------------------------------------- 
     405      ! in order to avoid lbc_lnk (communications): 
     406      !    jj loop must be 1:jpj   if adv_x is called first 
     407      !                and 2:jpj-1 if adv_x is called second 
     408      jj0 = NINT(pcrh) 
    377409      ! 
    378410      jcat = SIZE( ps0 , 3 )   ! size of input arrays 
     
    381413         ! 
    382414         ! Limitation of moments.                                            
    383          DO_2D( 0, 0, 1, 1 ) 
    384             !  Initialize volumes of boxes  (=area if adv_x first called, =psm otherwise)                                      
    385             psm (ji,jj,jl) = MAX( pcrh * e1e2t(ji,jj) + ( 1.0 - pcrh ) * psm(ji,jj,jl) , epsi20 ) 
    386             ! 
    387             zslpmax = MAX( 0._wp, ps0(ji,jj,jl) ) 
    388             zs1max  = 1.5 * zslpmax 
    389             zs1new  = MIN( zs1max, MAX( -zs1max, psx(ji,jj,jl) ) ) 
    390             zs2new  = MIN(  2.0 * zslpmax - 0.3334 * ABS( zs1new ),      & 
    391                &            MAX( ABS( zs1new ) - zslpmax, psxx(ji,jj,jl) )  ) 
    392             rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask 
    393  
    394             ps0 (ji,jj,jl) = zslpmax   
    395             psx (ji,jj,jl) = zs1new         * rswitch 
    396             psxx(ji,jj,jl) = zs2new         * rswitch 
    397             psy (ji,jj,jl) = psy (ji,jj,jl) * rswitch 
    398             psyy(ji,jj,jl) = psyy(ji,jj,jl) * rswitch 
    399             psxy(ji,jj,jl) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj,jl) ) ) * rswitch 
    400          END_2D 
    401  
    402          !  Calculate fluxes and moments between boxes i<-->i+1               
    403          DO_2D( 0, 0, 1, 1 )                   !  Flux from i to i+1 WHEN u GT 0 
    404             zbet(ji,jj)  =  MAX( 0._wp, SIGN( 1._wp, put(ji,jj) ) ) 
    405             zalf         =  MAX( 0._wp, put(ji,jj) ) * pdt / psm(ji,jj,jl) 
    406             zalfq        =  zalf * zalf 
    407             zalf1        =  1.0 - zalf 
    408             zalf1q       =  zalf1 * zalf1 
    409             ! 
    410             zfm (ji,jj)  =  zalf  *   psm (ji,jj,jl) 
    411             zf0 (ji,jj)  =  zalf  * ( ps0 (ji,jj,jl) + zalf1 * ( psx(ji,jj,jl) + (zalf1 - zalf) * psxx(ji,jj,jl) ) ) 
    412             zfx (ji,jj)  =  zalfq * ( psx (ji,jj,jl) + 3.0 * zalf1 * psxx(ji,jj,jl) ) 
    413             zfxx(ji,jj)  =  zalf  *   psxx(ji,jj,jl) * zalfq 
    414             zfy (ji,jj)  =  zalf  * ( psy (ji,jj,jl) + zalf1 * psxy(ji,jj,jl) ) 
    415             zfxy(ji,jj)  =  zalfq *   psxy(ji,jj,jl) 
    416             zfyy(ji,jj)  =  zalf  *   psyy(ji,jj,jl) 
    417  
    418             !  Readjust moments remaining in the box. 
    419             psm (ji,jj,jl)  =  psm (ji,jj,jl) - zfm(ji,jj) 
    420             ps0 (ji,jj,jl)  =  ps0 (ji,jj,jl) - zf0(ji,jj) 
    421             psx (ji,jj,jl)  =  zalf1q * ( psx(ji,jj,jl) - 3.0 * zalf * psxx(ji,jj,jl) ) 
    422             psxx(ji,jj,jl)  =  zalf1  * zalf1q * psxx(ji,jj,jl) 
    423             psy (ji,jj,jl)  =  psy (ji,jj,jl) - zfy(ji,jj) 
    424             psyy(ji,jj,jl)  =  psyy(ji,jj,jl) - zfyy(ji,jj) 
    425             psxy(ji,jj,jl)  =  zalf1q * psxy(ji,jj,jl) 
    426          END_2D 
    427  
    428          DO_2D( 0, 0, 1, 0 )                   !  Flux from i+1 to i when u LT 0. 
    429             zalf          = MAX( 0._wp, -put(ji,jj) ) * pdt / psm(ji+1,jj,jl)  
    430             zalg  (ji,jj) = zalf 
    431             zalfq         = zalf * zalf 
    432             zalf1         = 1.0 - zalf 
    433             zalg1 (ji,jj) = zalf1 
    434             zalf1q        = zalf1 * zalf1 
    435             zalg1q(ji,jj) = zalf1q 
    436             ! 
    437             zfm   (ji,jj) = zfm (ji,jj) + zalf  *    psm (ji+1,jj,jl) 
    438             zf0   (ji,jj) = zf0 (ji,jj) + zalf  * (  ps0 (ji+1,jj,jl) & 
    439                &                                   - zalf1 * ( psx(ji+1,jj,jl) - (zalf1 - zalf ) * psxx(ji+1,jj,jl) ) ) 
    440             zfx   (ji,jj) = zfx (ji,jj) + zalfq * (  psx (ji+1,jj,jl) - 3.0 * zalf1 * psxx(ji+1,jj,jl) ) 
    441             zfxx  (ji,jj) = zfxx(ji,jj) + zalf  *    psxx(ji+1,jj,jl) * zalfq 
    442             zfy   (ji,jj) = zfy (ji,jj) + zalf  * (  psy (ji+1,jj,jl) - zalf1 * psxy(ji+1,jj,jl) ) 
    443             zfxy  (ji,jj) = zfxy(ji,jj) + zalfq *    psxy(ji+1,jj,jl) 
    444             zfyy  (ji,jj) = zfyy(ji,jj) + zalf  *    psyy(ji+1,jj,jl) 
    445          END_2D 
    446  
    447          DO_2D( 0, 0, 0, 0 )                   !  Readjust moments remaining in the box. 
    448             zbt  =       zbet(ji-1,jj) 
    449             zbt1 = 1.0 - zbet(ji-1,jj) 
    450             ! 
    451             psm (ji,jj,jl) = zbt * psm(ji,jj,jl) + zbt1 * ( psm(ji,jj,jl) - zfm(ji-1,jj) ) 
    452             ps0 (ji,jj,jl) = zbt * ps0(ji,jj,jl) + zbt1 * ( ps0(ji,jj,jl) - zf0(ji-1,jj) ) 
    453             psx (ji,jj,jl) = zalg1q(ji-1,jj) * ( psx(ji,jj,jl) + 3.0 * zalg(ji-1,jj) * psxx(ji,jj,jl) ) 
    454             psxx(ji,jj,jl) = zalg1 (ji-1,jj) * zalg1q(ji-1,jj) * psxx(ji,jj,jl) 
    455             psy (ji,jj,jl) = zbt * psy (ji,jj,jl) + zbt1 * ( psy (ji,jj,jl) - zfy (ji-1,jj) ) 
    456             psyy(ji,jj,jl) = zbt * psyy(ji,jj,jl) + zbt1 * ( psyy(ji,jj,jl) - zfyy(ji-1,jj) ) 
    457             psxy(ji,jj,jl) = zalg1q(ji-1,jj) * psxy(ji,jj,jl) 
    458          END_2D 
    459  
    460          !   Put the temporary moments into appropriate neighboring boxes.     
    461          DO_2D( 0, 0, 0, 0 )                   !   Flux from i to i+1 IF u GT 0. 
    462             zbt  =       zbet(ji-1,jj) 
    463             zbt1 = 1.0 - zbet(ji-1,jj) 
    464             psm(ji,jj,jl) = zbt * ( psm(ji,jj,jl) + zfm(ji-1,jj) ) + zbt1 * psm(ji,jj,jl) 
    465             zalf          = zbt * zfm(ji-1,jj) / psm(ji,jj,jl) 
    466             zalf1         = 1.0 - zalf 
    467             ztemp         = zalf * ps0(ji,jj,jl) - zalf1 * zf0(ji-1,jj) 
    468             ! 
    469             ps0 (ji,jj,jl) =  zbt  * ( ps0(ji,jj,jl) + zf0(ji-1,jj) ) + zbt1 * ps0(ji,jj,jl) 
    470             psx (ji,jj,jl) =  zbt  * ( zalf * zfx(ji-1,jj) + zalf1 * psx(ji,jj,jl) + 3.0 * ztemp ) + zbt1 * psx(ji,jj,jl) 
    471             psxx(ji,jj,jl) =  zbt  * ( zalf * zalf * zfxx(ji-1,jj) + zalf1 * zalf1 * psxx(ji,jj,jl)                             & 
    472                &                     + 5.0 * ( zalf * zalf1 * ( psx (ji,jj,jl) - zfx(ji-1,jj) ) - ( zalf1 - zalf ) * ztemp )  ) & 
    473                &            + zbt1 * psxx(ji,jj,jl) 
    474             psxy(ji,jj,jl) =  zbt  * ( zalf * zfxy(ji-1,jj) + zalf1 * psxy(ji,jj,jl)             & 
    475                &                     + 3.0 * (- zalf1*zfy(ji-1,jj)  + zalf * psy(ji,jj,jl) ) )   & 
    476                &            + zbt1 * psxy(ji,jj,jl) 
    477             psy (ji,jj,jl) =  zbt  * ( psy (ji,jj,jl) + zfy (ji-1,jj) ) + zbt1 * psy (ji,jj,jl) 
    478             psyy(ji,jj,jl) =  zbt  * ( psyy(ji,jj,jl) + zfyy(ji-1,jj) ) + zbt1 * psyy(ji,jj,jl) 
    479          END_2D 
    480  
    481          DO_2D( 0, 0, 0, 0 )                   !  Flux from i+1 to i IF u LT 0. 
    482             zbt  =       zbet(ji,jj) 
    483             zbt1 = 1.0 - zbet(ji,jj) 
    484             psm(ji,jj,jl) = zbt * psm(ji,jj,jl) + zbt1 * ( psm(ji,jj,jl) + zfm(ji,jj) ) 
    485             zalf          = zbt1 * zfm(ji,jj) / psm(ji,jj,jl) 
    486             zalf1         = 1.0 - zalf 
    487             ztemp         = - zalf * ps0(ji,jj,jl) + zalf1 * zf0(ji,jj) 
    488             ! 
    489             ps0 (ji,jj,jl) = zbt * ps0 (ji,jj,jl) + zbt1 * ( ps0(ji,jj,jl) + zf0(ji,jj) ) 
    490             psx (ji,jj,jl) = zbt * psx (ji,jj,jl) + zbt1 * ( zalf * zfx(ji,jj) + zalf1 * psx(ji,jj,jl) + 3.0 * ztemp ) 
    491             psxx(ji,jj,jl) = zbt * psxx(ji,jj,jl) + zbt1 * ( zalf * zalf * zfxx(ji,jj) + zalf1 * zalf1 * psxx(ji,jj,jl) & 
    492                &                                           + 5.0 * ( zalf * zalf1 * ( - psx(ji,jj,jl) + zfx(ji,jj) )    & 
    493                &                                           + ( zalf1 - zalf ) * ztemp ) ) 
    494             psxy(ji,jj,jl) = zbt * psxy(ji,jj,jl) + zbt1 * ( zalf * zfxy(ji,jj) + zalf1 * psxy(ji,jj,jl)  & 
    495                &                                           + 3.0 * ( zalf1 * zfy(ji,jj) - zalf * psy(ji,jj,jl) ) ) 
    496             psy (ji,jj,jl) = zbt * psy (ji,jj,jl) + zbt1 * ( psy (ji,jj,jl) + zfy (ji,jj) ) 
    497             psyy(ji,jj,jl) = zbt * psyy(ji,jj,jl) + zbt1 * ( psyy(ji,jj,jl) + zfyy(ji,jj) ) 
    498          END_2D 
    499  
     415         DO jj = Njs0 - jj0, Nje0 + jj0 
     416             
     417            DO ji = Nis0 - 1, Nie0 + 1 
     418 
     419               !  Initialize volumes of boxes  (=area if adv_x first called, =psm otherwise)                                      
     420               psm (ji,jj,jl) = MAX( pcrh * e1e2t(ji,jj) + ( 1.0 - pcrh ) * psm(ji,jj,jl) , epsi20 ) 
     421               ! 
     422               zslpmax = MAX( 0._wp, ps0(ji,jj,jl) ) 
     423               zs1max  = 1.5 * zslpmax 
     424               zs1new  = MIN( zs1max, MAX( -zs1max, psx(ji,jj,jl) ) ) 
     425               zs2new  = MIN( 2.0 * zslpmax - 0.3334 * ABS( zs1new ), MAX( ABS( zs1new ) - zslpmax, psxx(ji,jj,jl) ) ) 
     426               rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask 
     427                
     428               ps0 (ji,jj,jl) = zslpmax   
     429               psx (ji,jj,jl) = zs1new         * rswitch 
     430               psxx(ji,jj,jl) = zs2new         * rswitch 
     431               psy (ji,jj,jl) = psy (ji,jj,jl) * rswitch 
     432               psyy(ji,jj,jl) = psyy(ji,jj,jl) * rswitch 
     433               psxy(ji,jj,jl) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj,jl) ) ) * rswitch 
     434                
     435               !  Calculate fluxes and moments between boxes i<-->i+1               
     436               !                                !  Flux from i to i+1 WHEN u GT 0  
     437               zbet(ji,jj)  =  MAX( 0._wp, SIGN( 1._wp, put(ji,jj) ) ) 
     438               zalf         =  MAX( 0._wp, put(ji,jj) ) * pdt / psm(ji,jj,jl) 
     439               zalfq        =  zalf * zalf 
     440               zalf1        =  1.0 - zalf 
     441               zalf1q       =  zalf1 * zalf1 
     442               ! 
     443               zfm (ji,jj)  =  zalf  *   psm (ji,jj,jl) 
     444               zf0 (ji,jj)  =  zalf  * ( ps0 (ji,jj,jl) + zalf1 * ( psx(ji,jj,jl) + (zalf1 - zalf) * psxx(ji,jj,jl) ) ) 
     445               zfx (ji,jj)  =  zalfq * ( psx (ji,jj,jl) + 3.0 * zalf1 * psxx(ji,jj,jl) ) 
     446               zfxx(ji,jj)  =  zalf  *   psxx(ji,jj,jl) * zalfq 
     447               zfy (ji,jj)  =  zalf  * ( psy (ji,jj,jl) + zalf1 * psxy(ji,jj,jl) ) 
     448               zfxy(ji,jj)  =  zalfq *   psxy(ji,jj,jl) 
     449               zfyy(ji,jj)  =  zalf  *   psyy(ji,jj,jl) 
     450                
     451               !                                !  Readjust moments remaining in the box. 
     452               psm (ji,jj,jl)  =  psm (ji,jj,jl) - zfm(ji,jj) 
     453               ps0 (ji,jj,jl)  =  ps0 (ji,jj,jl) - zf0(ji,jj) 
     454               psx (ji,jj,jl)  =  zalf1q * ( psx(ji,jj,jl) - 3.0 * zalf * psxx(ji,jj,jl) ) 
     455               psxx(ji,jj,jl)  =  zalf1  * zalf1q * psxx(ji,jj,jl) 
     456               psy (ji,jj,jl)  =  psy (ji,jj,jl) - zfy(ji,jj) 
     457               psyy(ji,jj,jl)  =  psyy(ji,jj,jl) - zfyy(ji,jj) 
     458               psxy(ji,jj,jl)  =  zalf1q * psxy(ji,jj,jl) 
     459            END DO 
     460             
     461            DO ji = Nis0 - 1, Nie0 
     462               !                                !  Flux from i+1 to i when u LT 0. 
     463               zalf          = MAX( 0._wp, -put(ji,jj) ) * pdt / psm(ji+1,jj,jl)  
     464               zalg  (ji,jj) = zalf 
     465               zalfq         = zalf * zalf 
     466               zalf1         = 1.0 - zalf 
     467               zalg1 (ji,jj) = zalf1 
     468               zalf1q        = zalf1 * zalf1 
     469               zalg1q(ji,jj) = zalf1q 
     470               ! 
     471               zfm   (ji,jj) = zfm (ji,jj) + zalf  *    psm (ji+1,jj,jl) 
     472               zf0   (ji,jj) = zf0 (ji,jj) + zalf  * (  ps0 (ji+1,jj,jl) & 
     473                  &                                   - zalf1 * ( psx(ji+1,jj,jl) - (zalf1 - zalf ) * psxx(ji+1,jj,jl) ) ) 
     474               zfx   (ji,jj) = zfx (ji,jj) + zalfq * (  psx (ji+1,jj,jl) - 3.0 * zalf1 * psxx(ji+1,jj,jl) ) 
     475               zfxx  (ji,jj) = zfxx(ji,jj) + zalf  *    psxx(ji+1,jj,jl) * zalfq 
     476               zfy   (ji,jj) = zfy (ji,jj) + zalf  * (  psy (ji+1,jj,jl) - zalf1 * psxy(ji+1,jj,jl) ) 
     477               zfxy  (ji,jj) = zfxy(ji,jj) + zalfq *    psxy(ji+1,jj,jl) 
     478               zfyy  (ji,jj) = zfyy(ji,jj) + zalf  *    psyy(ji+1,jj,jl) 
     479            END DO 
     480 
     481            DO ji = Nis0, Nie0 
     482               !                                !  Readjust moments remaining in the box. 
     483               zbt  =       zbet(ji-1,jj) 
     484               zbt1 = 1.0 - zbet(ji-1,jj) 
     485               ! 
     486               psm (ji,jj,jl) = zbt * psm(ji,jj,jl) + zbt1 * ( psm(ji,jj,jl) - zfm(ji-1,jj) ) 
     487               ps0 (ji,jj,jl) = zbt * ps0(ji,jj,jl) + zbt1 * ( ps0(ji,jj,jl) - zf0(ji-1,jj) ) 
     488               psx (ji,jj,jl) = zalg1q(ji-1,jj) * ( psx(ji,jj,jl) + 3.0 * zalg(ji-1,jj) * psxx(ji,jj,jl) ) 
     489               psxx(ji,jj,jl) = zalg1 (ji-1,jj) * zalg1q(ji-1,jj) * psxx(ji,jj,jl) 
     490               psy (ji,jj,jl) = zbt * psy (ji,jj,jl) + zbt1 * ( psy (ji,jj,jl) - zfy (ji-1,jj) ) 
     491               psyy(ji,jj,jl) = zbt * psyy(ji,jj,jl) + zbt1 * ( psyy(ji,jj,jl) - zfyy(ji-1,jj) ) 
     492               psxy(ji,jj,jl) = zalg1q(ji-1,jj) * psxy(ji,jj,jl) 
     493                
     494               !   Put the temporary moments into appropriate neighboring boxes.     
     495               !                                !   Flux from i to i+1 IF u GT 0. 
     496               zbt  =       zbet(ji-1,jj) 
     497               zbt1 = 1.0 - zbet(ji-1,jj) 
     498               psm(ji,jj,jl) = zbt * ( psm(ji,jj,jl) + zfm(ji-1,jj) ) + zbt1 * psm(ji,jj,jl) 
     499               zalf          = zbt * zfm(ji-1,jj) / psm(ji,jj,jl) 
     500               zalf1         = 1.0 - zalf 
     501               ztemp         = zalf * ps0(ji,jj,jl) - zalf1 * zf0(ji-1,jj) 
     502               ! 
     503               ps0 (ji,jj,jl) =  zbt  * ( ps0(ji,jj,jl) + zf0(ji-1,jj) ) + zbt1 * ps0(ji,jj,jl) 
     504               psx (ji,jj,jl) =  zbt  * ( zalf * zfx(ji-1,jj) + zalf1 * psx(ji,jj,jl) + 3.0 * ztemp ) + zbt1 * psx(ji,jj,jl) 
     505               psxx(ji,jj,jl) =  zbt  * ( zalf * zalf * zfxx(ji-1,jj) + zalf1 * zalf1 * psxx(ji,jj,jl)                             & 
     506                  &                     + 5.0 * ( zalf * zalf1 * ( psx (ji,jj,jl) - zfx(ji-1,jj) ) - ( zalf1 - zalf ) * ztemp )  ) & 
     507                  &            + zbt1 * psxx(ji,jj,jl) 
     508               psxy(ji,jj,jl) =  zbt  * ( zalf * zfxy(ji-1,jj) + zalf1 * psxy(ji,jj,jl)             & 
     509                  &                     + 3.0 * (- zalf1*zfy(ji-1,jj)  + zalf * psy(ji,jj,jl) ) )   & 
     510                  &            + zbt1 * psxy(ji,jj,jl) 
     511               psy (ji,jj,jl) =  zbt  * ( psy (ji,jj,jl) + zfy (ji-1,jj) ) + zbt1 * psy (ji,jj,jl) 
     512               psyy(ji,jj,jl) =  zbt  * ( psyy(ji,jj,jl) + zfyy(ji-1,jj) ) + zbt1 * psyy(ji,jj,jl) 
     513                
     514               !                                !  Flux from i+1 to i IF u LT 0. 
     515               zbt  =       zbet(ji,jj) 
     516               zbt1 = 1.0 - zbet(ji,jj) 
     517               psm(ji,jj,jl) = zbt * psm(ji,jj,jl) + zbt1 * ( psm(ji,jj,jl) + zfm(ji,jj) ) 
     518               zalf          = zbt1 * zfm(ji,jj) / psm(ji,jj,jl) 
     519               zalf1         = 1.0 - zalf 
     520               ztemp         = - zalf * ps0(ji,jj,jl) + zalf1 * zf0(ji,jj) 
     521               ! 
     522               ps0 (ji,jj,jl) = zbt * ps0 (ji,jj,jl) + zbt1 * ( ps0(ji,jj,jl) + zf0(ji,jj) ) 
     523               psx (ji,jj,jl) = zbt * psx (ji,jj,jl) + zbt1 * ( zalf * zfx(ji,jj) + zalf1 * psx(ji,jj,jl) + 3.0 * ztemp ) 
     524               psxx(ji,jj,jl) = zbt * psxx(ji,jj,jl) + zbt1 * ( zalf * zalf * zfxx(ji,jj) + zalf1 * zalf1 * psxx(ji,jj,jl) & 
     525                  &                                           + 5.0 * ( zalf * zalf1 * ( - psx(ji,jj,jl) + zfx(ji,jj) )    & 
     526                  &                                           + ( zalf1 - zalf ) * ztemp ) ) 
     527               psxy(ji,jj,jl) = zbt * psxy(ji,jj,jl) + zbt1 * ( zalf * zfxy(ji,jj) + zalf1 * psxy(ji,jj,jl)  & 
     528                  &                                           + 3.0 * ( zalf1 * zfy(ji,jj) - zalf * psy(ji,jj,jl) ) ) 
     529               psy (ji,jj,jl) = zbt * psy (ji,jj,jl) + zbt1 * ( psy (ji,jj,jl) + zfy (ji,jj) ) 
     530               psyy(ji,jj,jl) = zbt * psyy(ji,jj,jl) + zbt1 * ( psyy(ji,jj,jl) + zfyy(ji,jj) ) 
     531            END DO 
     532            ! 
     533         END DO 
     534         ! 
    500535      END DO 
    501  
    502       !-- Lateral boundary conditions 
    503       CALL lbc_lnk_multi( 'icedyn_adv_pra', psm(:,:,1:jcat) , 'T',  1.0_wp, ps0 , 'T',  1.0_wp   & 
    504          &                                , psx             , 'T', -1.0_wp, psy , 'T', -1.0_wp   &   ! caution gradient ==> the sign changes 
    505          &                                , psxx            , 'T',  1.0_wp, psyy, 'T',  1.0_wp , psxy, 'T',  1.0_wp ) 
    506       ! 
     536      !       
    507537   END SUBROUTINE adv_x 
    508538 
     
    525555      !! 
    526556      INTEGER  ::   ji, jj, jl, jcat                     ! dummy loop indices 
     557      INTEGER  ::   ji0                                  ! dummy loop indices 
    527558      REAL(wp) ::   zs1max, zslpmax, ztemp               ! temporary scalars 
    528559      REAL(wp) ::   zs1new, zalf , zalfq , zbt           !    -         - 
     
    532563      REAL(wp), DIMENSION(jpi,jpj) ::   zalg, zalg1, zalg1q     !  -      - 
    533564      !--------------------------------------------------------------------- 
     565      ! in order to avoid lbc_lnk (communications): 
     566      !    ji loop must be 1:jpi   if adv_y is called first 
     567      !                and 2:jpi-1 if adv_y is called second 
     568      ji0 = NINT(pcrh) 
    534569      ! 
    535570      jcat = SIZE( ps0 , 3 )   ! size of input arrays 
     
    538573         ! 
    539574         ! Limitation of moments. 
    540          DO_2D( 1, 1, 0, 0 ) 
     575         DO_2D( 1, 1, ji0, ji0 ) 
    541576            !  Initialize volumes of boxes (=area if adv_x first called, =psm otherwise) 
    542577            psm(ji,jj,jl) = MAX(  pcrh * e1e2t(ji,jj) + ( 1.0 - pcrh ) * psm(ji,jj,jl) , epsi20  ) 
     
    545580            zs1max  = 1.5 * zslpmax 
    546581            zs1new  = MIN( zs1max, MAX( -zs1max, psy(ji,jj,jl) ) ) 
    547             zs2new  = MIN(  ( 2.0 * zslpmax - 0.3334 * ABS( zs1new ) ),   & 
    548                &             MAX( ABS( zs1new )-zslpmax, psyy(ji,jj,jl) )  ) 
     582            zs2new  = MIN( ( 2.0 * zslpmax - 0.3334 * ABS( zs1new ) ), MAX( ABS( zs1new )-zslpmax, psyy(ji,jj,jl) ) ) 
    549583            rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask 
    550584            ! 
     
    555589            psyy(ji,jj,jl) = zs2new         * rswitch 
    556590            psxy(ji,jj,jl) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj,jl) ) ) * rswitch 
    557          END_2D 
    558591  
    559          !  Calculate fluxes and moments between boxes j<-->j+1               
    560          DO_2D( 1, 1, 0, 0 )                   !  Flux from j to j+1 WHEN v GT 0 
     592            !  Calculate fluxes and moments between boxes j<-->j+1               
     593            !                                !  Flux from j to j+1 WHEN v GT 0    
    561594            zbet(ji,jj)  =  MAX( 0._wp, SIGN( 1._wp, pvt(ji,jj) ) ) 
    562595            zalf         =  MAX( 0._wp, pvt(ji,jj) ) * pdt / psm(ji,jj,jl) 
     
    573606            zfxx(ji,jj)  =  zalf  * psxx(ji,jj,jl) 
    574607            ! 
    575             !  Readjust moments remaining in the box. 
     608            !                                !  Readjust moments remaining in the box. 
    576609            psm (ji,jj,jl)  =  psm (ji,jj,jl) - zfm(ji,jj) 
    577610            ps0 (ji,jj,jl)  =  ps0 (ji,jj,jl) - zf0(ji,jj) 
     
    583616         END_2D 
    584617         ! 
    585          DO_2D( 1, 0, 0, 0 )                   !  Flux from j+1 to j when v LT 0. 
     618         DO_2D( 1, 0, ji0, ji0 ) 
     619            !                                !  Flux from j+1 to j when v LT 0. 
    586620            zalf          = MAX( 0._wp, -pvt(ji,jj) ) * pdt / psm(ji,jj+1,jl)  
    587621            zalg  (ji,jj) = zalf 
     
    602636         END_2D 
    603637 
    604          !  Readjust moments remaining in the box.  
    605          DO_2D( 0, 0, 0, 0 ) 
     638         DO_2D( 0, 0, ji0, ji0 ) 
     639            !                                !  Readjust moments remaining in the box. 
    606640            zbt  =         zbet(ji,jj-1) 
    607641            zbt1 = ( 1.0 - zbet(ji,jj-1) ) 
     
    614648            psxx(ji,jj,jl) = zbt * psxx(ji,jj,jl) + zbt1 * ( psxx(ji,jj,jl) - zfxx(ji,jj-1) ) 
    615649            psxy(ji,jj,jl) = zalg1q(ji,jj-1) * psxy(ji,jj,jl) 
    616          END_2D 
    617  
    618          !   Put the temporary moments into appropriate neighboring boxes.     
    619          DO_2D( 0, 0, 0, 0 )                   !  Flux from j to j+1 IF v GT 0. 
     650 
     651            !   Put the temporary moments into appropriate neighboring boxes.     
     652            !                                !   Flux from j to j+1 IF v GT 0. 
    620653            zbt  =       zbet(ji,jj-1) 
    621654            zbt1 = 1.0 - zbet(ji,jj-1) 
     
    636669            psx (ji,jj,jl) =   zbt * ( psx (ji,jj,jl) + zfx (ji,jj-1) ) + zbt1 * psx (ji,jj,jl) 
    637670            psxx(ji,jj,jl) =   zbt * ( psxx(ji,jj,jl) + zfxx(ji,jj-1) ) + zbt1 * psxx(ji,jj,jl) 
    638          END_2D 
    639  
    640          DO_2D( 0, 0, 0, 0 )                   !  Flux from j+1 to j IF v LT 0. 
     671 
     672            !                                !  Flux from j+1 to j IF v LT 0. 
    641673            zbt  =       zbet(ji,jj) 
    642674            zbt1 = 1.0 - zbet(ji,jj) 
     
    656688            psxx(ji,jj,jl) = zbt * psxx(ji,jj,jl) + zbt1 * ( psxx(ji,jj,jl) + zfxx(ji,jj) ) 
    657689         END_2D 
    658  
     690         ! 
    659691      END DO 
    660  
    661       !-- Lateral boundary conditions 
    662       CALL lbc_lnk_multi( 'icedyn_adv_pra', psm(:,:,1:jcat) , 'T',  1.0_wp, ps0 , 'T',  1.0_wp   & 
    663          &                                , psx             , 'T', -1.0_wp, psy , 'T', -1.0_wp   &   ! caution gradient ==> the sign changes 
    664          &                                , psxx            , 'T',  1.0_wp, psyy, 'T',  1.0_wp , psxy, 'T',  1.0_wp ) 
    665692      ! 
    666693   END SUBROUTINE adv_y 
     
    890917            ! 
    891918            !                                                        ! ice thickness 
    892             CALL iom_get( numrir, jpdom_auto, 'sxice' , sxice ) 
    893             CALL iom_get( numrir, jpdom_auto, 'syice' , syice ) 
     919            CALL iom_get( numrir, jpdom_auto, 'sxice' , sxice , psgn = -1._wp ) 
     920            CALL iom_get( numrir, jpdom_auto, 'syice' , syice , psgn = -1._wp ) 
    894921            CALL iom_get( numrir, jpdom_auto, 'sxxice', sxxice ) 
    895922            CALL iom_get( numrir, jpdom_auto, 'syyice', syyice ) 
    896923            CALL iom_get( numrir, jpdom_auto, 'sxyice', sxyice ) 
    897924            !                                                        ! snow thickness 
    898             CALL iom_get( numrir, jpdom_auto, 'sxsn'  , sxsn  ) 
    899             CALL iom_get( numrir, jpdom_auto, 'sysn'  , sysn  ) 
     925            CALL iom_get( numrir, jpdom_auto, 'sxsn'  , sxsn  , psgn = -1._wp ) 
     926            CALL iom_get( numrir, jpdom_auto, 'sysn'  , sysn  , psgn = -1._wp ) 
    900927            CALL iom_get( numrir, jpdom_auto, 'sxxsn' , sxxsn  ) 
    901928            CALL iom_get( numrir, jpdom_auto, 'syysn' , syysn  ) 
    902929            CALL iom_get( numrir, jpdom_auto, 'sxysn' , sxysn  ) 
    903930            !                                                        ! ice concentration 
    904             CALL iom_get( numrir, jpdom_auto, 'sxa'   , sxa    ) 
    905             CALL iom_get( numrir, jpdom_auto, 'sya'   , sya    ) 
     931            CALL iom_get( numrir, jpdom_auto, 'sxa'   , sxa   , psgn = -1._wp ) 
     932            CALL iom_get( numrir, jpdom_auto, 'sya'   , sya   , psgn = -1._wp ) 
    906933            CALL iom_get( numrir, jpdom_auto, 'sxxa'  , sxxa   ) 
    907934            CALL iom_get( numrir, jpdom_auto, 'syya'  , syya   ) 
    908935            CALL iom_get( numrir, jpdom_auto, 'sxya'  , sxya   ) 
    909936            !                                                        ! ice salinity 
    910             CALL iom_get( numrir, jpdom_auto, 'sxsal' , sxsal ) 
    911             CALL iom_get( numrir, jpdom_auto, 'sysal' , sysal ) 
     937            CALL iom_get( numrir, jpdom_auto, 'sxsal' , sxsal , psgn = -1._wp ) 
     938            CALL iom_get( numrir, jpdom_auto, 'sysal' , sysal , psgn = -1._wp ) 
    912939            CALL iom_get( numrir, jpdom_auto, 'sxxsal', sxxsal ) 
    913940            CALL iom_get( numrir, jpdom_auto, 'syysal', syysal ) 
    914941            CALL iom_get( numrir, jpdom_auto, 'sxysal', sxysal ) 
    915942            !                                                        ! ice age 
    916             CALL iom_get( numrir, jpdom_auto, 'sxage' , sxage ) 
    917             CALL iom_get( numrir, jpdom_auto, 'syage' , syage ) 
     943            CALL iom_get( numrir, jpdom_auto, 'sxage' , sxage , psgn = -1._wp ) 
     944            CALL iom_get( numrir, jpdom_auto, 'syage' , syage , psgn = -1._wp ) 
    918945            CALL iom_get( numrir, jpdom_auto, 'sxxage', sxxage ) 
    919946            CALL iom_get( numrir, jpdom_auto, 'syyage', syyage ) 
     
    922949            DO jk = 1, nlay_s 
    923950               WRITE(zchar1,'(I2.2)') jk 
    924                znam = 'sxc0'//'_l'//zchar1  ;   CALL iom_get( numrir, jpdom_auto, znam , z3d )   ;   sxc0 (:,:,jk,:) = z3d(:,:,:) 
    925                znam = 'syc0'//'_l'//zchar1  ;   CALL iom_get( numrir, jpdom_auto, znam , z3d )   ;   syc0 (:,:,jk,:) = z3d(:,:,:) 
     951               znam = 'sxc0'//'_l'//zchar1  ;   CALL iom_get( numrir, jpdom_auto, znam , z3d, psgn = -1._wp )   ;   sxc0 (:,:,jk,:) = z3d(:,:,:) 
     952               znam = 'syc0'//'_l'//zchar1  ;   CALL iom_get( numrir, jpdom_auto, znam , z3d, psgn = -1._wp )   ;   syc0 (:,:,jk,:) = z3d(:,:,:) 
    926953               znam = 'sxxc0'//'_l'//zchar1 ;   CALL iom_get( numrir, jpdom_auto, znam , z3d )   ;   sxxc0(:,:,jk,:) = z3d(:,:,:) 
    927954               znam = 'syyc0'//'_l'//zchar1 ;   CALL iom_get( numrir, jpdom_auto, znam , z3d )   ;   syyc0(:,:,jk,:) = z3d(:,:,:) 
     
    931958            DO jk = 1, nlay_i 
    932959               WRITE(zchar1,'(I2.2)') jk 
    933                znam = 'sxe'//'_l'//zchar1   ;   CALL iom_get( numrir, jpdom_auto, znam , z3d )   ;   sxe (:,:,jk,:) = z3d(:,:,:) 
    934                znam = 'sye'//'_l'//zchar1   ;   CALL iom_get( numrir, jpdom_auto, znam , z3d )   ;   sye (:,:,jk,:) = z3d(:,:,:) 
     960               znam = 'sxe'//'_l'//zchar1   ;   CALL iom_get( numrir, jpdom_auto, znam , z3d, psgn = -1._wp )   ;   sxe (:,:,jk,:) = z3d(:,:,:) 
     961               znam = 'sye'//'_l'//zchar1   ;   CALL iom_get( numrir, jpdom_auto, znam , z3d, psgn = -1._wp )   ;   sye (:,:,jk,:) = z3d(:,:,:) 
    935962               znam = 'sxxe'//'_l'//zchar1  ;   CALL iom_get( numrir, jpdom_auto, znam , z3d )   ;   sxxe(:,:,jk,:) = z3d(:,:,:) 
    936963               znam = 'syye'//'_l'//zchar1  ;   CALL iom_get( numrir, jpdom_auto, znam , z3d )   ;   syye(:,:,jk,:) = z3d(:,:,:) 
     
    940967            IF( ln_pnd_LEV ) THEN                                    ! melt pond fraction 
    941968               IF( iom_varid( numror, 'sxap', ldstop = .FALSE. ) > 0 ) THEN 
    942                   CALL iom_get( numrir, jpdom_auto, 'sxap' , sxap ) 
    943                   CALL iom_get( numrir, jpdom_auto, 'syap' , syap ) 
     969                  CALL iom_get( numrir, jpdom_auto, 'sxap' , sxap , psgn = -1._wp ) 
     970                  CALL iom_get( numrir, jpdom_auto, 'syap' , syap , psgn = -1._wp ) 
    944971                  CALL iom_get( numrir, jpdom_auto, 'sxxap', sxxap ) 
    945972                  CALL iom_get( numrir, jpdom_auto, 'syyap', syyap ) 
    946973                  CALL iom_get( numrir, jpdom_auto, 'sxyap', sxyap ) 
    947974                  !                                                     ! melt pond volume 
    948                   CALL iom_get( numrir, jpdom_auto, 'sxvp' , sxvp ) 
    949                   CALL iom_get( numrir, jpdom_auto, 'syvp' , syvp ) 
     975                  CALL iom_get( numrir, jpdom_auto, 'sxvp' , sxvp , psgn = -1._wp ) 
     976                  CALL iom_get( numrir, jpdom_auto, 'syvp' , syvp , psgn = -1._wp ) 
    950977                  CALL iom_get( numrir, jpdom_auto, 'sxxvp', sxxvp ) 
    951978                  CALL iom_get( numrir, jpdom_auto, 'syyvp', syyvp ) 
     
    958985               IF ( ln_pnd_lids ) THEN                               ! melt pond lid volume 
    959986                  IF( iom_varid( numror, 'sxvl', ldstop = .FALSE. ) > 0 ) THEN 
    960                      CALL iom_get( numrir, jpdom_auto, 'sxvl' , sxvl ) 
    961                      CALL iom_get( numrir, jpdom_auto, 'syvl' , syvl ) 
     987                     CALL iom_get( numrir, jpdom_auto, 'sxvl' , sxvl , psgn = -1._wp ) 
     988                     CALL iom_get( numrir, jpdom_auto, 'syvl' , syvl , psgn = -1._wp ) 
    962989                     CALL iom_get( numrir, jpdom_auto, 'sxxvl', sxxvl ) 
    963990                     CALL iom_get( numrir, jpdom_auto, 'syyvl', syyvl ) 
Note: See TracChangeset for help on using the changeset viewer.