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 13618 for NEMO/trunk/src/ICE/icedyn_adv_umx.F90 – NEMO

Ignore:
Timestamp:
2020-10-16T10:07:55+02:00 (4 years ago)
Author:
clem
Message:

trunk: optimization of both Prather and UMx advection schemes. Related to ticket #2552

File:
1 edited

Legend:

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

    r13609 r13618  
    115115      ELSEWHERE                      ; zs_i(:,:,:) = 0._wp 
    116116      END WHERE 
    117       DO jl = 1, jpl 
    118          DO_2D( 0, 0, 0, 0 ) 
    119             zhip_max(ji,jj,jl) = MAX( epsi20, ph_ip(ji,jj,jl), ph_ip(ji+1,jj  ,jl), ph_ip(ji  ,jj+1,jl), & 
    120                &                                               ph_ip(ji-1,jj  ,jl), ph_ip(ji  ,jj-1,jl), & 
    121                &                                               ph_ip(ji+1,jj+1,jl), ph_ip(ji-1,jj-1,jl), & 
    122                &                                               ph_ip(ji+1,jj-1,jl), ph_ip(ji-1,jj+1,jl) ) 
    123             zhi_max (ji,jj,jl) = MAX( epsi20, ph_i (ji,jj,jl), ph_i (ji+1,jj  ,jl), ph_i (ji  ,jj+1,jl), & 
    124                &                                               ph_i (ji-1,jj  ,jl), ph_i (ji  ,jj-1,jl), & 
    125                &                                               ph_i (ji+1,jj+1,jl), ph_i (ji-1,jj-1,jl), & 
    126                &                                               ph_i (ji+1,jj-1,jl), ph_i (ji-1,jj+1,jl) ) 
    127             zhs_max (ji,jj,jl) = MAX( epsi20, ph_s (ji,jj,jl), ph_s (ji+1,jj  ,jl), ph_s (ji  ,jj+1,jl), & 
    128                &                                               ph_s (ji-1,jj  ,jl), ph_s (ji  ,jj-1,jl), & 
    129                &                                               ph_s (ji+1,jj+1,jl), ph_s (ji-1,jj-1,jl), & 
    130                &                                               ph_s (ji+1,jj-1,jl), ph_s (ji-1,jj+1,jl) ) 
    131             zsi_max (ji,jj,jl) = MAX( epsi20, zs_i (ji,jj,jl), zs_i (ji+1,jj  ,jl), zs_i (ji  ,jj+1,jl), & 
    132                &                                               zs_i (ji-1,jj  ,jl), zs_i (ji  ,jj-1,jl), & 
    133                &                                               zs_i (ji+1,jj+1,jl), zs_i (ji-1,jj-1,jl), & 
    134                &                                               zs_i (ji+1,jj-1,jl), zs_i (ji-1,jj+1,jl) ) 
    135          END_2D 
    136       END DO 
     117      CALL icemax3D( ph_i , zhi_max ) 
     118      CALL icemax3D( ph_s , zhs_max ) 
     119      CALL icemax3D( ph_ip, zhip_max) 
     120      CALL icemax3D( zs_i , zsi_max ) 
    137121      CALL lbc_lnk_multi( 'icedyn_adv_umx', zhi_max, 'T', 1._wp, zhs_max, 'T', 1._wp, zhip_max, 'T', 1._wp, zsi_max, 'T', 1._wp ) 
    138122      ! 
     
    147131         ELSEWHERE                      ; ze_s(:,:,jk,:) = 0._wp 
    148132         END WHERE 
    149       END DO 
    150       DO jl = 1, jpl 
    151          DO_3D( 0, 0, 0, 0, 1, nlay_i ) 
    152             zei_max(ji,jj,jk,jl) = MAX( epsi20, ze_i(ji,jj,jk,jl), ze_i(ji+1,jj  ,jk,jl), ze_i(ji  ,jj+1,jk,jl), & 
    153                &                                                   ze_i(ji-1,jj  ,jk,jl), ze_i(ji  ,jj-1,jk,jl), & 
    154                &                                                   ze_i(ji+1,jj+1,jk,jl), ze_i(ji-1,jj-1,jk,jl), & 
    155                &                                                   ze_i(ji+1,jj-1,jk,jl), ze_i(ji-1,jj+1,jk,jl) ) 
    156          END_3D 
    157       END DO 
    158       DO jl = 1, jpl 
    159          DO_3D( 0, 0, 0, 0, 1, nlay_s ) 
    160             zes_max(ji,jj,jk,jl) = MAX( epsi20, ze_s(ji,jj,jk,jl), ze_s(ji+1,jj  ,jk,jl), ze_s(ji  ,jj+1,jk,jl), & 
    161                &                                                   ze_s(ji-1,jj  ,jk,jl), ze_s(ji  ,jj-1,jk,jl), & 
    162                &                                                   ze_s(ji+1,jj+1,jk,jl), ze_s(ji-1,jj-1,jk,jl), & 
    163                &                                                   ze_s(ji+1,jj-1,jk,jl), ze_s(ji-1,jj+1,jk,jl) ) 
    164          END_3D 
    165       END DO 
    166       CALL lbc_lnk( 'icedyn_adv_pra', zei_max, 'T', 1. ) 
    167       CALL lbc_lnk( 'icedyn_adv_pra', zes_max, 'T', 1. ) 
     133      END DO    
     134      CALL icemax4D( ze_i , zei_max ) 
     135      CALL icemax4D( ze_s , zes_max ) 
     136      CALL lbc_lnk( 'icedyn_adv_umx', zei_max, 'T', 1._wp ) 
     137      CALL lbc_lnk( 'icedyn_adv_umx', zes_max, 'T', 1._wp ) 
    168138      ! 
    169139      ! 
     
    386356            ENDIF 
    387357         ENDIF 
     358 
     359         ! --- Lateral boundary conditions --- ! 
     360         IF    ( ln_pnd_LEV .AND. ln_pnd_lids ) THEN 
     361            CALL lbc_lnk_multi( 'icedyn_adv_umx', pa_i,'T',1._wp, pv_i,'T',1._wp, pv_s,'T',1._wp, psv_i,'T',1._wp, poa_i,'T',1._wp & 
     362               &                                , pa_ip,'T',1._wp, pv_ip,'T',1._wp, pv_il,'T',1._wp ) 
     363         ELSEIF( ln_pnd_LEV .AND. .NOT.ln_pnd_lids ) THEN 
     364            CALL lbc_lnk_multi( 'icedyn_adv_umx', pa_i,'T',1._wp, pv_i,'T',1._wp, pv_s,'T',1._wp, psv_i,'T',1._wp, poa_i,'T',1._wp & 
     365               &                                , pa_ip,'T',1._wp, pv_ip,'T',1._wp ) 
     366         ELSE 
     367            CALL lbc_lnk_multi( 'icedyn_adv_umx', pa_i,'T',1._wp, pv_i,'T',1._wp, pv_s,'T',1._wp, psv_i,'T',1._wp, poa_i,'T',1._wp ) 
     368         ENDIF 
     369         CALL lbc_lnk( 'icedyn_adv_umx', pe_i, 'T', 1._wp ) 
     370         CALL lbc_lnk( 'icedyn_adv_umx', pe_s, 'T', 1._wp ) 
    388371         ! 
    389372         !== Open water area ==! 
     
    393376               &                          - ( zudy(ji,jj) - zudy(ji-1,jj) + zvdx(ji,jj) - zvdx(ji,jj-1) ) * r1_e1e2t(ji,jj) * zdt 
    394377         END_2D 
    395          CALL lbc_lnk( 'icedyn_adv_umx', pato_i, 'T',  1.0_wp ) 
     378         CALL lbc_lnk( 'icedyn_adv_umx', pato_i, 'T',  1._wp ) 
    396379         ! 
    397380         ! --- diagnostics --- ! 
     
    462445      !!             work on H (and not V). It is partly related to the multi-category approach 
    463446      !!             Therefore, after advection we limit the thickness to the largest value of the 9-points around (only if ice 
    464       !!             concentration is small). Since we do not limit S and T, large values can occur at the edge but it does not really matter 
    465       !!             since sv_i and e_i are still good. 
     447      !!             concentration is small). We also limit S and T. 
    466448      !!---------------------------------------------------------------------- 
    467449      REAL(wp)                        , INTENT(in   )           ::   pamsk            ! advection of concentration (1) or other tracers (0) 
     
    507489      IF( pamsk == 0._wp ) THEN 
    508490         DO jl = 1, jpl 
    509             DO_2D( 1, 0, 1, 0 ) 
     491            DO_2D( 0, 0, 1, 0 ) 
    510492               IF( ABS( pu(ji,jj) ) > epsi10 ) THEN 
    511493                  zfu_ho (ji,jj,jl) = zfu_ho (ji,jj,jl) * puc    (ji,jj,jl) / pu(ji,jj) 
     
    516498               ENDIF 
    517499               ! 
     500            END_2D 
     501            DO_2D( 1, 0, 0, 0 ) 
    518502               IF( ABS( pv(ji,jj) ) > epsi10 ) THEN 
    519503                  zfv_ho (ji,jj,jl) = zfv_ho (ji,jj,jl) * pvc    (ji,jj,jl) / pv(ji,jj) 
     
    550534      IF( PRESENT( pua_ho ) ) THEN 
    551535         DO jl = 1, jpl 
    552             DO_2D( 1, 0, 1, 0 ) 
    553                pua_ho (ji,jj,jl) = zfu_ho (ji,jj,jl) ; pva_ho (ji,jj,jl) = zfv_ho (ji,jj,jl) 
    554                pua_ups(ji,jj,jl) = zfu_ups(ji,jj,jl) ; pva_ups(ji,jj,jl) = zfv_ups(ji,jj,jl) 
     536            DO_2D( 0, 0, 1, 0 ) 
     537               pua_ho (ji,jj,jl) = zfu_ho (ji,jj,jl) 
     538               pua_ups(ji,jj,jl) = zfu_ups(ji,jj,jl) 
     539            END_2D 
     540            DO_2D( 1, 0, 0, 0 ) 
     541               pva_ho (ji,jj,jl) = zfv_ho (ji,jj,jl) 
     542               pva_ups(ji,jj,jl) = zfv_ups(ji,jj,jl) 
    555543            END_2D 
    556544         END DO 
     
    566554         END_2D 
    567555      END DO 
    568       CALL lbc_lnk( 'icedyn_adv_umx', ptc, 'T',  1.0_wp ) 
    569556      ! 
    570557   END SUBROUTINE adv_umx 
     
    605592            ! 
    606593            DO jl = 1, jpl              !-- flux in x-direction 
    607                DO_2D( 1, 0, 1, 0 ) 
     594               DO_2D( 1, 1, 1, 0 ) 
    608595                  pfu_ups(ji,jj,jl) = MAX( pu(ji,jj), 0._wp ) * pt(ji,jj,jl) + MIN( pu(ji,jj), 0._wp ) * pt(ji+1,jj,jl) 
    609596               END_2D 
     
    611598            ! 
    612599            DO jl = 1, jpl              !-- first guess of tracer from u-flux 
    613                DO_2D( 0, 0, 0, 0 ) 
     600               DO_2D( 1, 1, 0, 0 ) 
    614601                  ztra = - ( pfu_ups(ji,jj,jl) - pfu_ups(ji-1,jj,jl) )              & 
    615602                     &   + ( pu     (ji,jj   ) - pu     (ji-1,jj   ) ) * pt(ji,jj,jl) * (1.-pamsk) 
     
    618605               END_2D 
    619606            END DO 
    620             CALL lbc_lnk( 'icedyn_adv_umx', zpt, 'T', 1.0_wp ) 
    621607            ! 
    622608            DO jl = 1, jpl              !-- flux in y-direction 
    623                DO_2D( 1, 0, 1, 0 ) 
     609               DO_2D( 1, 0, 0, 0 ) 
    624610                  pfv_ups(ji,jj,jl) = MAX( pv(ji,jj), 0._wp ) * zpt(ji,jj,jl) + MIN( pv(ji,jj), 0._wp ) * zpt(ji,jj+1,jl) 
    625611               END_2D 
     
    629615            ! 
    630616            DO jl = 1, jpl              !-- flux in y-direction 
    631                DO_2D( 1, 0, 1, 0 ) 
     617               DO_2D( 1, 0, 1, 1 ) 
    632618                  pfv_ups(ji,jj,jl) = MAX( pv(ji,jj), 0._wp ) * pt(ji,jj,jl) + MIN( pv(ji,jj), 0._wp ) * pt(ji,jj+1,jl) 
    633619               END_2D 
     
    635621            ! 
    636622            DO jl = 1, jpl              !-- first guess of tracer from v-flux 
    637                DO_2D( 0, 0, 0, 0 ) 
     623               DO_2D( 0, 0, 1, 1 ) 
    638624                  ztra = - ( pfv_ups(ji,jj,jl) - pfv_ups(ji,jj-1,jl) )  & 
    639625                     &   + ( pv     (ji,jj   ) - pv     (ji,jj-1   ) ) * pt(ji,jj,jl) * (1.-pamsk) 
     
    642628               END_2D 
    643629            END DO 
    644             CALL lbc_lnk( 'icedyn_adv_umx', zpt, 'T', 1.0_wp ) 
    645630            ! 
    646631            DO jl = 1, jpl              !-- flux in x-direction 
    647                DO_2D( 1, 0, 1, 0 ) 
     632               DO_2D( 0, 0, 1, 0 ) 
    648633                  pfu_ups(ji,jj,jl) = MAX( pu(ji,jj), 0._wp ) * zpt(ji,jj,jl) + MIN( pu(ji,jj), 0._wp ) * zpt(ji+1,jj,jl) 
    649634               END_2D 
     
    694679         ! 
    695680         DO jl = 1, jpl 
    696             DO_2D( 1, 0, 1, 0 ) 
     681            DO_2D( 1, 1, 1, 0 ) 
    697682               pfu_ho(ji,jj,jl) = 0.5_wp * pu(ji,jj) * ( pt(ji,jj,jl) + pt(ji+1,jj  ,jl) ) 
     683            END_2D 
     684            DO_2D( 1, 0, 1, 1 ) 
    698685               pfv_ho(ji,jj,jl) = 0.5_wp * pv(ji,jj) * ( pt(ji,jj,jl) + pt(ji  ,jj+1,jl) ) 
    699686            END_2D 
     
    712699            ! 
    713700            DO jl = 1, jpl              !-- flux in x-direction 
    714                DO_2D( 1, 0, 1, 0 ) 
     701               DO_2D( 1, 1, 1, 0 ) 
    715702                  pfu_ho(ji,jj,jl) = 0.5_wp * pu(ji,jj) * ( pt(ji,jj,jl) + pt(ji+1,jj,jl) ) 
    716703               END_2D 
     
    719706 
    720707            DO jl = 1, jpl              !-- first guess of tracer from u-flux 
    721                DO_2D( 0, 0, 0, 0 ) 
     708               DO_2D( 1, 1, 0, 0 ) 
    722709                  ztra = - ( pfu_ho(ji,jj,jl) - pfu_ho(ji-1,jj,jl) )              & 
    723710                     &   + ( pu    (ji,jj   ) - pu    (ji-1,jj   ) ) * pt(ji,jj,jl) * (1.-pamsk) 
     
    726713               END_2D 
    727714            END DO 
    728             CALL lbc_lnk( 'icedyn_adv_umx', zpt, 'T', 1.0_wp ) 
    729715 
    730716            DO jl = 1, jpl              !-- flux in y-direction 
    731                DO_2D( 1, 0, 1, 0 ) 
     717               DO_2D( 1, 0, 0, 0 ) 
    732718                  pfv_ho(ji,jj,jl) = 0.5_wp * pv(ji,jj) * ( zpt(ji,jj,jl) + zpt(ji,jj+1,jl) ) 
    733719               END_2D 
     
    738724            ! 
    739725            DO jl = 1, jpl              !-- flux in y-direction 
    740                DO_2D( 1, 0, 1, 0 ) 
     726               DO_2D( 1, 0, 1, 1 ) 
    741727                  pfv_ho(ji,jj,jl) = 0.5_wp * pv(ji,jj) * ( pt(ji,jj,jl) + pt(ji,jj+1,jl) ) 
    742728               END_2D 
     
    745731            ! 
    746732            DO jl = 1, jpl              !-- first guess of tracer from v-flux 
    747                DO_2D( 0, 0, 0, 0 ) 
     733               DO_2D( 0, 0, 1, 1 ) 
    748734                  ztra = - ( pfv_ho(ji,jj,jl) - pfv_ho(ji,jj-1,jl) )  & 
    749735                     &   + ( pv    (ji,jj   ) - pv    (ji,jj-1   ) ) * pt(ji,jj,jl) * (1.-pamsk) 
     
    752738               END_2D 
    753739            END DO 
    754             CALL lbc_lnk( 'icedyn_adv_umx', zpt, 'T', 1.0_wp ) 
    755740            ! 
    756741            DO jl = 1, jpl              !-- flux in x-direction 
    757                DO_2D( 1, 0, 1, 0 ) 
     742               DO_2D( 0, 0, 1, 0 ) 
    758743                  pfu_ho(ji,jj,jl) = 0.5_wp * pu(ji,jj) * ( zpt(ji,jj,jl) + zpt(ji+1,jj,jl) ) 
    759744               END_2D 
     
    912897         !         
    913898         DO jl = 1, jpl 
    914             DO_2D( 1, 0, 1, 0 ) 
     899            DO_2D( 0, 0, 1, 0 ) 
    915900               pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * (                                pt(ji+1,jj,jl) + pt(ji,jj,jl)   & 
    916901                  &                                         - SIGN( 1._wp, pu(ji,jj) ) * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) ) 
     
    921906         ! 
    922907         DO jl = 1, jpl 
    923             DO_2D( 1, 0, 1, 0 ) 
     908            DO_2D( 0, 0, 1, 0 ) 
    924909               zcu  = pu(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 
    925910               pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * (                                pt(ji+1,jj,jl) + pt(ji,jj,jl)   & 
     
    931916         ! 
    932917         DO jl = 1, jpl 
    933             DO_2D( 1, 0, 1, 0 ) 
     918            DO_2D( 0, 0, 1, 0 ) 
    934919               zcu  = pu(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 
    935920               zdx2 = e1u(ji,jj) * e1u(ji,jj) 
     
    945930         ! 
    946931         DO jl = 1, jpl 
    947             DO_2D( 1, 0, 1, 0 ) 
     932            DO_2D( 0, 0, 1, 0 ) 
    948933               zcu  = pu(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 
    949934               zdx2 = e1u(ji,jj) * e1u(ji,jj) 
     
    959944         ! 
    960945         DO jl = 1, jpl 
    961             DO_2D( 1, 0, 1, 0 ) 
     946            DO_2D( 0, 0, 1, 0 ) 
    962947               zcu  = pu(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 
    963948               zdx2 = e1u(ji,jj) * e1u(ji,jj) 
     
    980965      IF( ll_neg ) THEN 
    981966         DO jl = 1, jpl 
    982             DO_2D( 1, 0, 1, 0 ) 
     967            DO_2D( 0, 0, 1, 0 ) 
    983968               IF( pt_u(ji,jj,jl) < 0._wp .OR. ( imsk_small(ji,jj,jl) == 0 .AND. pamsk == 0. ) ) THEN 
    984969                  pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * (                                pt(ji+1,jj,jl) + pt(ji,jj,jl)   & 
     
    10481033      CASE( 1 )                                                !==  1st order central TIM  ==! (Eq. 21) 
    10491034         DO jl = 1, jpl 
    1050             DO_2D( 1, 0, 1, 0 ) 
     1035            DO_2D( 1, 0, 0, 0 ) 
    10511036               pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * (                                pt(ji,jj+1,jl) + pt(ji,jj,jl)   & 
    10521037                  &                                         - SIGN( 1._wp, pv(ji,jj) ) * ( pt(ji,jj+1,jl) - pt(ji,jj,jl) ) ) 
     
    10561041      CASE( 2 )                                                !==  2nd order central TIM  ==! (Eq. 23) 
    10571042         DO jl = 1, jpl 
    1058             DO_2D( 1, 0, 1, 0 ) 
     1043            DO_2D( 1, 0, 0, 0 ) 
    10591044               zcv  = pv(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 
    10601045               pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * (                                pt(ji,jj+1,jl) + pt(ji,jj,jl)   & 
     
    10651050      CASE( 3 )                                                !==  3rd order central TIM  ==! (Eq. 24) 
    10661051         DO jl = 1, jpl 
    1067             DO_2D( 1, 0, 1, 0 ) 
     1052            DO_2D( 1, 0, 0, 0 ) 
    10681053               zcv  = pv(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 
    10691054               zdy2 = e2v(ji,jj) * e2v(ji,jj) 
     
    10781063      CASE( 4 )                                                !==  4th order central TIM  ==! (Eq. 27) 
    10791064         DO jl = 1, jpl 
    1080             DO_2D( 1, 0, 1, 0 ) 
     1065            DO_2D( 1, 0, 0, 0 ) 
    10811066               zcv  = pv(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 
    10821067               zdy2 = e2v(ji,jj) * e2v(ji,jj) 
     
    10911076      CASE( 5 )                                                !==  5th order central TIM  ==! (Eq. 29) 
    10921077         DO jl = 1, jpl 
    1093             DO_2D( 1, 0, 1, 0 ) 
     1078            DO_2D( 1, 0, 0, 0 ) 
    10941079               zcv  = pv(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 
    10951080               zdy2 = e2v(ji,jj) * e2v(ji,jj) 
     
    11121097      IF( ll_neg ) THEN 
    11131098         DO jl = 1, jpl 
    1114             DO_2D( 1, 0, 1, 0 ) 
     1099            DO_2D( 1, 0, 0, 0 ) 
    11151100               IF( pt_v(ji,jj,jl) < 0._wp .OR. ( jmsk_small(ji,jj,jl) == 0 .AND. pamsk == 0. ) ) THEN 
    11161101                  pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * (                              ( pt(ji,jj+1,jl) + pt(ji,jj,jl) )  & 
     
    11221107      !                                                     !-- High order flux in j-direction  --! 
    11231108      DO jl = 1, jpl 
    1124          DO_2D( 1, 0, 1, 0 ) 
     1109         DO_2D( 1, 0, 0, 0 ) 
    11251110            pfv_ho(ji,jj,jl) = pv(ji,jj) * pt_v(ji,jj,jl) 
    11261111         END_2D 
     
    11581143      ! -------------------------------------------------- 
    11591144      DO jl = 1, jpl 
    1160          DO_2D( 1, 0, 1, 0 ) 
     1145         DO_2D( 0, 0, 1, 0 ) 
    11611146            pfu_ho(ji,jj,jl) = pfu_ho(ji,jj,jl) - pfu_ups(ji,jj,jl) 
     1147         END_2D 
     1148         DO_2D( 1, 0, 0, 0 ) 
    11621149            pfv_ho(ji,jj,jl) = pfv_ho(ji,jj,jl) - pfv_ups(ji,jj,jl) 
    11631150         END_2D 
     
    12651252      ! --------------------------------- 
    12661253      DO jl = 1, jpl 
    1267          DO_2D( 1, 0, 1, 0 ) 
     1254         DO_2D( 0, 0, 1, 0 ) 
    12681255            zau = MIN( 1._wp , zbetdo(ji,jj,jl) , zbetup(ji+1,jj,jl) ) 
    12691256            zbu = MIN( 1._wp , zbetup(ji,jj,jl) , zbetdo(ji+1,jj,jl) ) 
     
    12761263         END_2D 
    12771264 
    1278          DO_2D( 1, 0, 1, 0 ) 
     1265         DO_2D( 1, 0, 0, 0 ) 
    12791266            zav = MIN( 1._wp , zbetdo(ji,jj,jl) , zbetup(ji,jj+1,jl) ) 
    12801267            zbv = MIN( 1._wp , zbetup(ji,jj,jl) , zbetdo(ji,jj+1,jl) ) 
     
    16331620   END SUBROUTINE Hsnow 
    16341621 
     1622   SUBROUTINE icemax3D( pice , pmax ) 
     1623      !!--------------------------------------------------------------------- 
     1624      !!                   ***  ROUTINE icemax3D ***                      
     1625      !! ** Purpose :  compute the max of the 9 points around 
     1626      !!---------------------------------------------------------------------- 
     1627      REAL(wp), DIMENSION(:,:,:)      , INTENT(in ) ::   pice   ! input 
     1628      REAL(wp), DIMENSION(:,:,:)      , INTENT(out) ::   pmax   ! output 
     1629      REAL(wp), DIMENSION(2:jpim1,jpj)              ::   zmax   ! temporary array 
     1630      INTEGER  ::   ji, jj, jl   ! dummy loop indices 
     1631      !!---------------------------------------------------------------------- 
     1632      DO jl = 1, jpl 
     1633         DO jj = Njs0-1, Nje0+1     
     1634            DO ji = Nis0, Nie0 
     1635               zmax(ji,jj) = MAX( epsi20, pice(ji,jj,jl), pice(ji-1,jj,jl), pice(ji+1,jj,jl) ) 
     1636            END DO 
     1637         END DO 
     1638         DO jj = Njs0, Nje0     
     1639            DO ji = Nis0, Nie0 
     1640               pmax(ji,jj,jl) = MAX( epsi20, zmax(ji,jj), zmax(ji,jj-1), zmax(ji,jj+1) ) 
     1641            END DO 
     1642         END DO 
     1643      END DO 
     1644   END SUBROUTINE icemax3D 
     1645 
     1646   SUBROUTINE icemax4D( pice , pmax ) 
     1647      !!--------------------------------------------------------------------- 
     1648      !!                   ***  ROUTINE icemax4D ***                      
     1649      !! ** Purpose :  compute the max of the 9 points around 
     1650      !!---------------------------------------------------------------------- 
     1651      REAL(wp), DIMENSION(:,:,:,:)    , INTENT(in ) ::   pice   ! input 
     1652      REAL(wp), DIMENSION(:,:,:,:)    , INTENT(out) ::   pmax   ! output 
     1653      REAL(wp), DIMENSION(2:jpim1,jpj)              ::   zmax   ! temporary array 
     1654      INTEGER  ::   jlay, ji, jj, jk, jl   ! dummy loop indices 
     1655      !!---------------------------------------------------------------------- 
     1656      jlay = SIZE( pice , 3 )   ! size of input arrays 
     1657      DO jl = 1, jpl 
     1658         DO jk = 1, jlay 
     1659            DO jj = Njs0-1, Nje0+1     
     1660               DO ji = Nis0, Nie0 
     1661                  zmax(ji,jj) = MAX( epsi20, pice(ji,jj,jk,jl), pice(ji-1,jj,jk,jl), pice(ji+1,jj,jk,jl) ) 
     1662               END DO 
     1663            END DO 
     1664            DO jj = Njs0, Nje0     
     1665               DO ji = Nis0, Nie0 
     1666                  pmax(ji,jj,jk,jl) = MAX( epsi20, zmax(ji,jj), zmax(ji,jj-1), zmax(ji,jj+1) ) 
     1667               END DO 
     1668            END DO 
     1669         END DO 
     1670      END DO 
     1671   END SUBROUTINE icemax4D 
    16351672 
    16361673#else 
Note: See TracChangeset for help on using the changeset viewer.