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

Changeset 8565


Ignore:
Timestamp:
2017-09-27T12:09:10+02:00 (7 years ago)
Author:
clem
Message:

trying to respect naming convention

Location:
branches/2017/dev_r8183_ICEMODEL/NEMOGCM
Files:
14 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/CONFIG/SHARED/namelist_ice_lim3_ref

    r8531 r8565  
    6868      rn_pstar      =   2.0e+04       !     ice strength thickness parameter [N/m2] 
    6969      rn_crhg       =   20.0          !     ice strength conc. parameter (-) 
    70    ln_str_R75       = .false.         !  ice strength param.: Rothrock_75 => P = Cf*coeff*integral(wr.h^2)     
    71       rn_perdg      =   17.0          !     ridging work divided by pot. energy change in ridging 
    7270                   ! -- ice_rdgrft -- ! 
    7371   rn_csrdg         =   0.5           !  fraction of shearing energy contributing to ridging 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice1D.F90

    r8564 r8565  
    2727   !: are the variables corresponding to 2d vectors 
    2828 
    29 !! please, DOCTOR naming convention   starting by i means LOCAL integer 
    30 !!         ===>>>  rename idxice  as nidice or nidthd, or nidx_thd or midx ... ? 
    31    INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   idxice !: selected points for ice thermo 
    32    INTEGER , PUBLIC                                  ::   nidx   !  number of selected points 
     29   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   nptidx !: selected points for ice thermo 
     30   INTEGER , PUBLIC                                  ::   npti   !  number of selected points 
    3331 
    3432   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   qlead_1d      
     
    175173 
    176174      ii = 1 
    177       ALLOCATE( idxice   (jpij) ,   & 
     175      ALLOCATE( nptidx   (jpij) ,   & 
    178176         &      qlead_1d (jpij) , ftr_ice_1d(jpij) , qsr_ice_1d (jpij) ,   & 
    179177         &      fr1_i0_1d(jpij) , fr2_i0_1d (jpij) , qns_ice_1d(jpij)  ,   & 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icealb.F90

    r8534 r8565  
    2929   REAL(wp), PUBLIC, PARAMETER ::   rn_alb_oce = 0.066   !: ocean or lead albedo (Pegau and Paulson, Ann. Glac. 2001) 
    3030 
    31    REAL(wp) , PARAMETER ::   rc1    = 0.05    ! snow thickness (only for nn_ice_alb=0) 
    32    REAL(wp) , PARAMETER ::   rc2    = 0.10    !  "        " 
    33    REAL(wp) , PARAMETER ::   rcloud = 0.06    ! cloud effect on albedo (only-for nn_ice_alb=0) 
    34    REAL(wp) , PARAMETER ::   r1_c1 = 1. / rc1 
    35    REAL(wp) , PARAMETER ::   r1_c2 = 1. / rc2 
     31   REAL(wp) , PARAMETER ::   ppc1    = 0.05    ! snow thickness (only for nn_ice_alb=0) 
     32   REAL(wp) , PARAMETER ::   ppc2    = 0.10    !  "        " 
     33   REAL(wp) , PARAMETER ::   ppcloud = 0.06    ! cloud effect on albedo (only-for nn_ice_alb=0) 
     34   REAL(wp) , PARAMETER ::   pp1_c1 = 1. / ppc1 
     35   REAL(wp) , PARAMETER ::   pp1_c2 = 1. / ppc2 
    3636   ! 
    3737   ! ** albedo namelist (namalb) 
     
    154154               DO ji = 1, jpi 
    155155                  !                    ! Freezing snow 
    156                   ! no effect of underlying ice layer IF snow thickness > c1. Albedo does not depend on snow thick if > rc2 
    157                   zswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ( ph_snw(ji,jj,jl) - rc1 ) ) ) 
     156                  ! no effect of underlying ice layer IF snow thickness > c1. Albedo does not depend on snow thick if > ppc2 
     157                  zswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ( ph_snw(ji,jj,jl) - ppc1 ) ) ) 
    158158                  zalb_sf = ( 1._wp - zswitch ) * (  zalb_it(ji,jj,jl)  & 
    159                      &                                   + ph_snw(ji,jj,jl) * ( rn_alb_sdry - zalb_it(ji,jj,jl) ) * r1_c1  )   & 
     159                     &                                   + ph_snw(ji,jj,jl) * ( rn_alb_sdry - zalb_it(ji,jj,jl) ) * pp1_c1  )   & 
    160160                     &    +           zswitch   * rn_alb_sdry   
    161161                     ! 
    162162                  !                    ! Melting snow 
    163                   ! no effect of underlying ice layer. Albedo does not depend on snow thick IF > rc2 
    164                   zswitch = MAX( 0._wp , SIGN( 1._wp , ph_snw(ji,jj,jl) - rc2 ) ) 
    165                   zalb_sm = ( 1._wp - zswitch ) * ( rn_alb_imlt + ph_snw(ji,jj,jl) * ( rn_alb_smlt - rn_alb_imlt ) * r1_c2 )   & 
     163                  ! no effect of underlying ice layer. Albedo does not depend on snow thick IF > ppc2 
     164                  zswitch = MAX( 0._wp , SIGN( 1._wp , ph_snw(ji,jj,jl) - ppc2 ) ) 
     165                  zalb_sm = ( 1._wp - zswitch ) * ( rn_alb_imlt + ph_snw(ji,jj,jl) * ( rn_alb_smlt - rn_alb_imlt ) * pp1_c2 )   & 
    166166                     &      +         zswitch   *   rn_alb_smlt  
    167167                     ! 
     
    178178         END DO 
    179179         ! 
    180          pa_ice_os(:,:,:) = pa_ice_cs(:,:,:) + rcloud       ! Oberhuber correction for overcast sky 
     180         pa_ice_os(:,:,:) = pa_ice_cs(:,:,:) + ppcloud       ! Oberhuber correction for overcast sky 
    181181         ! 
    182182         ! 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icectl.F90

    r8564 r8565  
    7171      REAL(wp) ::   zvtrp, zetrp 
    7272      REAL(wp) ::   zarea, zv_sill, zs_sill, zt_sill 
    73       REAL(wp), PARAMETER ::   zconv = 1.e-9 ! convert W to GW and kg to Mt 
     73      REAL(wp), PARAMETER ::   ppconv = 1.e-9 ! convert W to GW and kg to Mt 
    7474      !!------------------------------------------------------------------- 
    7575      ! 
     
    7979            &                    wfx_opw(:,:) + wfx_res(:,:)     + wfx_dyn(:,:)     + wfx_lam(:,:)     + wfx_ice_sub(:,:) +  & 
    8080            &                    wfx_snw_sni(:,:) + wfx_snw_sum(:,:) + wfx_snw_dyn(:,:) + wfx_snw_sub(:,:) + wfx_spr(:,:)    & 
    81             &                  ) * e1e2t(:,:) ) * zconv 
     81            &                  ) * e1e2t(:,:) ) * ppconv 
    8282         ! 
    8383         !                          ! salt flux 
    8484         pdiag_fs = glob_sum(  ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) +  & 
    8585            &                    sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) + sfx_sub(:,:) + sfx_lam(:,:)    & 
    86             &                  ) *  e1e2t(:,:) ) * zconv  
     86            &                  ) *  e1e2t(:,:) ) * ppconv  
    8787         ! 
    8888         !                          ! heat flux 
    8989         pdiag_ft = glob_sum(  ( hfx_sum(:,:) + hfx_bom(:,:) + hfx_bog(:,:) + hfx_dif(:,:) + hfx_opw(:,:) + hfx_snw(:,:)  &  
    9090            &                  - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) - hfx_sub(:,:) - hfx_spr(:,:)   & 
    91             &                  ) *  e1e2t(:,:) ) * zconv 
    92  
    93          pdiag_v = glob_sum( SUM( v_i * rhoic + v_s * rhosn, dim=3 ) * e1e2t * zconv ) 
    94  
    95          pdiag_s = glob_sum( SUM( sv_i * rhoic            , dim=3 ) * e1e2t * zconv ) 
     91            &                  ) *  e1e2t(:,:) ) * ppconv 
     92 
     93         pdiag_v = glob_sum( SUM( v_i * rhoic + v_s * rhosn, dim=3 ) * e1e2t * ppconv ) 
     94 
     95         pdiag_s = glob_sum( SUM( sv_i * rhoic            , dim=3 ) * e1e2t * ppconv ) 
    9696 
    9797         pdiag_t = glob_sum( (  SUM( SUM( e_i(:,:,1:nlay_i,:), dim=4 ), dim=3 )     & 
    98             &                 + SUM( SUM( e_s(:,:,1:nlay_s,:), dim=4 ), dim=3 )   ) * e1e2t ) * zconv 
     98            &                 + SUM( SUM( e_s(:,:,1:nlay_s,:), dim=4 ), dim=3 )   ) * e1e2t ) * ppconv 
    9999 
    100100      ELSEIF( icount == 1 ) THEN 
     
    104104            &                wfx_opw(:,:) + wfx_res(:,:)     + wfx_dyn(:,:)     + wfx_lam(:,:)     + wfx_ice_sub(:,:) +  & 
    105105            &                wfx_snw_sni(:,:) + wfx_snw_sum(:,:) + wfx_snw_dyn(:,:) + wfx_snw_sub(:,:) + wfx_spr(:,:)        & 
    106             &              ) * e1e2t(:,:) ) * zconv - pdiag_fv 
     106            &              ) * e1e2t(:,:) ) * ppconv - pdiag_fv 
    107107 
    108108         ! salt flux 
    109109         zfs  = glob_sum(  ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) +  & 
    110110            &                sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) + sfx_sub(:,:) + sfx_lam(:,:)    &  
    111             &              ) * e1e2t(:,:) ) * zconv - pdiag_fs 
     111            &              ) * e1e2t(:,:) ) * ppconv - pdiag_fs 
    112112 
    113113         ! heat flux 
    114114         zft  = glob_sum(  ( hfx_sum(:,:) + hfx_bom(:,:) + hfx_bog(:,:) + hfx_dif(:,:) + hfx_opw(:,:) + hfx_snw(:,:)  &  
    115115            &              - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) - hfx_sub(:,:) - hfx_spr(:,:)   & 
    116             &              ) * e1e2t(:,:) ) * zconv - pdiag_ft 
     116            &              ) * e1e2t(:,:) ) * ppconv - pdiag_ft 
    117117  
    118118         ! outputs 
    119          zv = ( ( glob_sum( SUM( v_i * rhoic + v_s * rhosn, dim=3 ) * e1e2t  ) * zconv  & 
     119         zv = ( ( glob_sum( SUM( v_i * rhoic + v_s * rhosn, dim=3 ) * e1e2t  ) * ppconv  & 
    120120            &     - pdiag_v ) * r1_rdtice - zfv ) * rday 
    121121 
    122          zs = ( ( glob_sum( SUM( sv_i * rhoic            , dim=3 ) * e1e2t ) * zconv  & 
     122         zs = ( ( glob_sum( SUM( sv_i * rhoic            , dim=3 ) * e1e2t ) * ppconv  & 
    123123            &     - pdiag_s ) * r1_rdtice + zfs ) * rday 
    124124 
    125125         zt = ( glob_sum( (  SUM( SUM( e_i(:,:,1:nlay_i,:), dim=4 ), dim=3 )   & 
    126             &                + SUM( SUM( e_s(:,:,1:nlay_s,:), dim=4 ), dim=3 ) ) * e1e2t ) * zconv   & 
     126            &                + SUM( SUM( e_s(:,:,1:nlay_s,:), dim=4 ), dim=3 ) ) * e1e2t ) * ppconv   & 
    127127            &   - pdiag_t ) * r1_rdtice + zft 
    128128 
    129129         ! zvtrp and zetrp must be close to 0 if the advection scheme is conservative 
    130          zvtrp = glob_sum( ( diag_trp_vi * rhoic + diag_trp_vs * rhosn ) * e1e2t  ) * zconv * rday  
    131          zetrp = glob_sum( ( diag_trp_ei         + diag_trp_es         ) * e1e2t  ) * zconv 
     130         zvtrp = glob_sum( ( diag_trp_vi * rhoic + diag_trp_vs * rhosn ) * e1e2t  ) * ppconv * rday  
     131         zetrp = glob_sum( ( diag_trp_ei         + diag_trp_es         ) * e1e2t  ) * ppconv 
    132132 
    133133         zvmin = glob_min( v_i ) 
     
    136136 
    137137         ! set threshold values and calculate the ice area (+epsi10 to set a threshold > 0 when there is no ice)  
    138          zarea   = glob_sum( SUM( a_i + epsi10, dim=3 ) * e1e2t ) * zconv ! in 1.e9 m2 
     138         zarea   = glob_sum( SUM( a_i + epsi10, dim=3 ) * e1e2t ) * ppconv ! in 1.e9 m2 
    139139         zv_sill = zarea * 2.5e-5 
    140140         zs_sill = zarea * 25.e-5 
     
    176176      REAL(wp)                        :: zhfx, zsfx, zvfx 
    177177      REAL(wp)                        :: zarea, zv_sill, zs_sill, zt_sill 
    178       REAL(wp), PARAMETER             :: zconv = 1.e-9 ! convert W to GW and kg to Mt 
     178      REAL(wp), PARAMETER             :: ppconv = 1.e-9 ! convert W to GW and kg to Mt 
    179179      !!------------------------------------------------------------------- 
    180180 
    181181      ! water flux 
    182       zvfx  = glob_sum( ( wfx_ice + wfx_snw + wfx_spr + wfx_sub + diag_vice + diag_vsnw ) * e1e2t ) * zconv * rday 
     182      zvfx  = glob_sum( ( wfx_ice + wfx_snw + wfx_spr + wfx_sub + diag_vice + diag_vsnw ) * e1e2t ) * ppconv * rday 
    183183 
    184184      ! salt flux 
    185       zsfx  = glob_sum( ( sfx + diag_sice ) * e1e2t ) * zconv * rday 
     185      zsfx  = glob_sum( ( sfx + diag_sice ) * e1e2t ) * ppconv * rday 
    186186 
    187187      ! heat flux 
    188188      zhfx  = glob_sum( ( hfx_in - hfx_out - diag_heat - diag_trp_ei - diag_trp_es   & 
    189189      !  &              - SUM( qevap_ice * a_i_b, dim=3 )                           & !!clem: I think this line must be commented (but need check) 
    190          &              ) * e1e2t ) * zconv 
     190         &              ) * e1e2t ) * ppconv 
    191191 
    192192      ! set threshold values and calculate the ice area (+epsi10 to set a threshold > 0 when there is no ice)  
    193       zarea   = glob_sum( SUM( a_i + epsi10, dim=3 ) * e1e2t ) * zconv ! in 1.e9 m2 
     193      zarea   = glob_sum( SUM( a_i + epsi10, dim=3 ) * e1e2t ) * ppconv ! in 1.e9 m2 
    194194      zv_sill = zarea * 2.5e-5 
    195195      zs_sill = zarea * 25.e-5 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icedyn_rdgrft.F90

    r8564 r8565  
    6161   REAL(wp) ::   rn_pstar         ! determines ice strength, Hibler JPO79 
    6262   REAL(wp) ::   rn_crhg          ! determines changes in ice strength 
    63    LOGICAL  ::   ln_str_R75       ! ice strength parameterization (Rothrock75) 
    64    REAL(wp) ::   rn_perdg         ! ridging work divided by pot. energy change in ridging 
    6563   REAL(wp) ::   rn_csrdg         ! fraction of shearing energy contributing to ridging             
    6664   LOGICAL  ::   ln_partf_lin     ! participation function linear (Thorndike et al. (1975)) 
     
    130128      !! 
    131129      INTEGER  ::   ji, jj, jk, jl             ! dummy loop index 
    132       INTEGER  ::   niter, iterate_ridging     ! local integer  
    133       INTEGER  ::   nidx2                      ! local integer 
     130      INTEGER  ::   iter, iterate_ridging      ! local integer  
     131      INTEGER  ::   ipti                       ! local integer 
    134132      REAL(wp) ::   zfac                       ! local scalar 
    135       INTEGER , DIMENSION(jpij) ::   idxice2       ! compute ridge/raft or not 
     133      INTEGER , DIMENSION(jpij) ::   iptidx        ! compute ridge/raft or not 
    136134      REAL(wp), DIMENSION(jpij) ::   zdivu_adv     ! divu as implied by transport scheme  (1/s) 
    137135      REAL(wp), DIMENSION(jpij) ::   zdivu, zdelt  ! 1D divu_i & delta_i 
    138136      ! 
    139       INTEGER, PARAMETER ::   nitermax = 20     
     137      INTEGER, PARAMETER ::   jp_itermax = 20     
    140138      !!------------------------------------------------------------------- 
    141139      ! controls 
     
    154152      ! 0) Identify grid cells with ice 
    155153      !-------------------------------- 
    156       nidx = 0   ;   idxice(:) = 0 
     154      npti = 0   ;   nptidx(:) = 0 
    157155      DO jj = 1, jpj 
    158156         DO ji = 1, jpi 
    159157            IF ( SUM(a_i(ji,jj,:)) > 0._wp ) THEN 
    160                nidx           = nidx + 1 
    161                idxice( nidx ) = (jj - 1) * jpi + ji 
     158               npti           = npti + 1 
     159               nptidx( npti ) = (jj - 1) * jpi + ji 
    162160            ENDIF 
    163161         END DO 
    164162      END DO 
    165163 
    166       IF( nidx > 0 ) THEN 
     164      IF( npti > 0 ) THEN 
    167165         
    168166         ! just needed here 
    169          CALL tab_2d_1d( nidx, idxice(1:nidx), zdivu(1:nidx), divu_i(:,:) ) 
    170          CALL tab_2d_1d( nidx, idxice(1:nidx), zdelt(1:nidx), delta_i(:,:) ) 
     167         CALL tab_2d_1d( npti, nptidx(1:npti), zdivu(1:npti), divu_i(:,:) ) 
     168         CALL tab_2d_1d( npti, nptidx(1:npti), zdelt(1:npti), delta_i(:,:) ) 
    171169         ! needed here and in the iteration loop 
    172          CALL tab_3d_2d( nidx, idxice(1:nidx), a_i_2d  (1:nidx,1:jpl), a_i(:,:,:) ) 
    173          CALL tab_3d_2d( nidx, idxice(1:nidx), v_i_2d  (1:nidx,1:jpl), v_i(:,:,:) ) 
    174          CALL tab_2d_1d( nidx, idxice(1:nidx), ato_i_1d(1:nidx)      , ato_i(:,:) ) 
    175  
    176          DO ji = 1, nidx 
     170         CALL tab_3d_2d( npti, nptidx(1:npti), a_i_2d  (1:npti,1:jpl), a_i(:,:,:) ) 
     171         CALL tab_3d_2d( npti, nptidx(1:npti), v_i_2d  (1:npti,1:jpl), v_i(:,:,:) ) 
     172         CALL tab_2d_1d( npti, nptidx(1:npti), ato_i_1d(1:npti)      , ato_i(:,:) ) 
     173 
     174         DO ji = 1, npti 
    177175            !-----------------------------------------------------------------------------! 
    178176            ! 2) Dynamical inputs (closing rate, divu_adv, opning) 
     
    219217         !------------------------------------------------ 
    220218         CALL rdgrft_prep( a_i_2d, v_i_2d, ato_i_1d ) 
    221          nidx2 = 0 ; idxice2(:) = 0 
    222          DO ji = 1, nidx 
     219         ipti = 0 ; iptidx(:) = 0 
     220         DO ji = 1, npti 
    223221            IF( SUM( ABS(apartf(ji,1:jpl)) ) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN 
    224                nidx2 = nidx2 + 1 
    225                idxice2    (nidx2)   = idxice     (ji) 
    226                a_i_2d     (nidx2,:) = a_i_2d     (ji,:)  ! adjust to new indices 
    227                v_i_2d     (nidx2,:) = v_i_2d     (ji,:)  ! adjust to new indices 
    228                ato_i_1d   (nidx2)   = ato_i_1d   (ji)    ! adjust to new indices 
    229                closing_net(nidx2)   = closing_net(ji)    ! adjust to new indices 
    230                zdivu_adv  (nidx2)   = zdivu_adv  (ji)    ! adjust to new indices 
    231                opning     (nidx2)   = opning     (ji)    ! adjust to new indices 
     222               ipti = ipti + 1 
     223               iptidx     (ipti)   = nptidx     (ji) 
     224               a_i_2d     (ipti,:) = a_i_2d     (ji,:)  ! adjust to new indices 
     225               v_i_2d     (ipti,:) = v_i_2d     (ji,:)  ! adjust to new indices 
     226               ato_i_1d   (ipti)   = ato_i_1d   (ji)    ! adjust to new indices 
     227               closing_net(ipti)   = closing_net(ji)    ! adjust to new indices 
     228               zdivu_adv  (ipti)   = zdivu_adv  (ji)    ! adjust to new indices 
     229               opning     (ipti)   = opning     (ji)    ! adjust to new indices 
    232230            ENDIF 
    233231         END DO 
    234          idxice(:) = idxice2(:) 
    235          nidx      = nidx2 
     232         nptidx(:) = iptidx(:) 
     233         npti      = ipti 
    236234 
    237235      ENDIF 
     
    240238      ! Start ridging 
    241239      !-------------- 
    242       IF( nidx > 0 ) THEN 
     240      IF( npti > 0 ) THEN 
    243241          
    244242         !----------- 
     
    246244         !----------- 
    247245         ! fields used but not modified 
    248          CALL tab_2d_1d( nidx, idxice(1:nidx), sss_1d(1:nidx), sss_m(:,:) ) 
    249          CALL tab_2d_1d( nidx, idxice(1:nidx), sst_1d(1:nidx), sst_m(:,:) ) 
     246         CALL tab_2d_1d( npti, nptidx(1:npti), sss_1d(1:npti), sss_m(:,:) ) 
     247         CALL tab_2d_1d( npti, nptidx(1:npti), sst_1d(1:npti), sst_m(:,:) ) 
    250248         ! the following fields are modified in this routine 
    251          !!CALL tab_2d_1d( nidx, idxice(1:nidx), ato_i_1d(1:nidx), ato_i(:,:) ) 
    252          !!CALL tab_3d_2d( nidx, idxice(1:nidx), a_i_2d(1:nidx,1:jpl), a_i(:,:,:) ) 
    253          !!CALL tab_3d_2d( nidx, idxice(1:nidx), v_i_2d  (1:nidx,1:jpl), v_i  (:,:,:) ) 
    254          CALL tab_3d_2d( nidx, idxice(1:nidx), v_s_2d (1:nidx,1:jpl), v_s (:,:,:) ) 
    255          CALL tab_3d_2d( nidx, idxice(1:nidx), sv_i_2d(1:nidx,1:jpl), sv_i(:,:,:) ) 
    256          CALL tab_3d_2d( nidx, idxice(1:nidx), oa_i_2d(1:nidx,1:jpl), oa_i(:,:,:) ) 
     249         !!CALL tab_2d_1d( npti, nptidx(1:npti), ato_i_1d(1:npti), ato_i(:,:) ) 
     250         !!CALL tab_3d_2d( npti, nptidx(1:npti), a_i_2d(1:npti,1:jpl), a_i(:,:,:) ) 
     251         !!CALL tab_3d_2d( npti, nptidx(1:npti), v_i_2d  (1:npti,1:jpl), v_i  (:,:,:) ) 
     252         CALL tab_3d_2d( npti, nptidx(1:npti), v_s_2d (1:npti,1:jpl), v_s (:,:,:) ) 
     253         CALL tab_3d_2d( npti, nptidx(1:npti), sv_i_2d(1:npti,1:jpl), sv_i(:,:,:) ) 
     254         CALL tab_3d_2d( npti, nptidx(1:npti), oa_i_2d(1:npti,1:jpl), oa_i(:,:,:) ) 
    257255         IF ( nn_pnd_scheme > 0 ) THEN 
    258             CALL tab_3d_2d( nidx, idxice(1:nidx), a_ip_2d(1:nidx,1:jpl), a_ip(:,:,:) ) 
    259             CALL tab_3d_2d( nidx, idxice(1:nidx), v_ip_2d(1:nidx,1:jpl), v_ip(:,:,:) ) 
     256            CALL tab_3d_2d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip(:,:,:) ) 
     257            CALL tab_3d_2d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip(:,:,:) ) 
    260258         ENDIF 
    261259         DO jl = 1, jpl 
    262260            DO jk = 1, nlay_s 
    263                CALL tab_2d_1d( nidx, idxice(1:nidx), ze_s_2d(1:nidx,jk,jl), e_s(:,:,jk,jl) ) 
     261               CALL tab_2d_1d( npti, nptidx(1:npti), ze_s_2d(1:npti,jk,jl), e_s(:,:,jk,jl) ) 
    264262            END DO 
    265263            DO jk = 1, nlay_i 
    266                CALL tab_2d_1d( nidx, idxice(1:nidx), ze_i_2d(1:nidx,jk,jl), e_i(:,:,jk,jl) ) 
    267             END DO 
    268          END DO 
    269          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_dyn_1d    (1:nidx), sfx_dyn    (:,:) ) 
    270          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_bri_1d    (1:nidx), sfx_bri    (:,:) ) 
    271          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_dyn_1d    (1:nidx), wfx_dyn    (:,:) ) 
    272          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_dyn_1d    (1:nidx), hfx_dyn    (:,:) ) 
    273          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_snw_dyn_1d(1:nidx), wfx_snw_dyn(:,:) ) 
     264               CALL tab_2d_1d( npti, nptidx(1:npti), ze_i_2d(1:npti,jk,jl), e_i(:,:,jk,jl) ) 
     265            END DO 
     266         END DO 
     267         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_dyn_1d    (1:npti), sfx_dyn    (:,:) ) 
     268         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_bri_1d    (1:npti), sfx_bri    (:,:) ) 
     269         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_dyn_1d    (1:npti), wfx_dyn    (:,:) ) 
     270         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_dyn_1d    (1:npti), hfx_dyn    (:,:) ) 
     271         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_snw_dyn_1d(1:npti), wfx_snw_dyn(:,:) ) 
    274272         IF ( nn_pnd_scheme > 0 ) THEN 
    275             CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_pnd_1d(1:nidx), wfx_pnd(:,:) ) 
     273            CALL tab_2d_1d( npti, nptidx(1:npti), wfx_pnd_1d(1:npti), wfx_pnd(:,:) ) 
    276274         ENDIF 
    277275          
     
    279277         ! 3) Ridging iteration 
    280278         !-----------------------------------------------------------------------------! 
    281          niter           = 1 
     279         iter            = 1 
    282280         iterate_ridging = 1       
    283          DO WHILE( iterate_ridging > 0 .AND. niter < nitermax ) 
     281         DO WHILE( iterate_ridging > 0 .AND. iter < jp_itermax ) 
    284282 
    285283            CALL rdgrft_prep( a_i_2d, v_i_2d, ato_i_1d ) 
     
    294292            ! rates were reduced above), ridge again with new rates. 
    295293            iterate_ridging = 0 
    296             DO ji = 1, nidx 
     294            DO ji = 1, npti 
    297295               zfac = 1._wp - ( ato_i_1d(ji) + SUM( a_i_2d(ji,:) ) ) 
    298296               IF( ABS( zfac ) < epsi10 ) THEN 
     
    308306            END DO 
    309307            ! 
    310             niter = niter + 1 
    311             IF( niter  >  nitermax )    CALL ctl_warn( 'icedyn_rdgrft: non-converging ridging scheme' ) 
     308            iter = iter + 1 
     309            IF( iter  >  jp_itermax )    CALL ctl_warn( 'icedyn_rdgrft: non-converging ridging scheme' ) 
    312310            ! 
    313311         END DO 
     
    316314         ! 1D <==> 2D 
    317315         !----------- 
    318          CALL tab_1d_2d( nidx, idxice(1:nidx), ato_i_1d(1:nidx), ato_i(:,:) ) 
    319          CALL tab_2d_3d( nidx, idxice(1:nidx), a_i_2d (1:nidx,1:jpl), a_i (:,:,:) ) 
    320          CALL tab_2d_3d( nidx, idxice(1:nidx), v_i_2d (1:nidx,1:jpl), v_i (:,:,:) ) 
    321          CALL tab_2d_3d( nidx, idxice(1:nidx), v_s_2d (1:nidx,1:jpl), v_s (:,:,:) ) 
    322          CALL tab_2d_3d( nidx, idxice(1:nidx), sv_i_2d(1:nidx,1:jpl), sv_i(:,:,:) ) 
    323          CALL tab_2d_3d( nidx, idxice(1:nidx), oa_i_2d(1:nidx,1:jpl), oa_i(:,:,:) ) 
     316         CALL tab_1d_2d( npti, nptidx(1:npti), ato_i_1d(1:npti), ato_i(:,:) ) 
     317         CALL tab_2d_3d( npti, nptidx(1:npti), a_i_2d (1:npti,1:jpl), a_i (:,:,:) ) 
     318         CALL tab_2d_3d( npti, nptidx(1:npti), v_i_2d (1:npti,1:jpl), v_i (:,:,:) ) 
     319         CALL tab_2d_3d( npti, nptidx(1:npti), v_s_2d (1:npti,1:jpl), v_s (:,:,:) ) 
     320         CALL tab_2d_3d( npti, nptidx(1:npti), sv_i_2d(1:npti,1:jpl), sv_i(:,:,:) ) 
     321         CALL tab_2d_3d( npti, nptidx(1:npti), oa_i_2d(1:npti,1:jpl), oa_i(:,:,:) ) 
    324322         IF ( nn_pnd_scheme > 0 ) THEN 
    325             CALL tab_2d_3d( nidx, idxice(1:nidx), a_ip_2d(1:nidx,1:jpl), a_ip(:,:,:) ) 
    326             CALL tab_2d_3d( nidx, idxice(1:nidx), v_ip_2d(1:nidx,1:jpl), v_ip(:,:,:) ) 
     323            CALL tab_2d_3d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip(:,:,:) ) 
     324            CALL tab_2d_3d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip(:,:,:) ) 
    327325         ENDIF 
    328326         DO jl = 1, jpl 
    329327            DO jk = 1, nlay_s 
    330                CALL tab_1d_2d( nidx, idxice(1:nidx), ze_s_2d(1:nidx,jk,jl), e_s(:,:,jk,jl) ) 
     328               CALL tab_1d_2d( npti, nptidx(1:npti), ze_s_2d(1:npti,jk,jl), e_s(:,:,jk,jl) ) 
    331329            END DO 
    332330            DO jk = 1, nlay_i 
    333                CALL tab_1d_2d( nidx, idxice(1:nidx), ze_i_2d(1:nidx,jk,jl), e_i(:,:,jk,jl) ) 
    334             END DO 
    335          END DO 
    336          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_dyn_1d    (1:nidx), sfx_dyn    (:,:) ) 
    337          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_bri_1d    (1:nidx), sfx_bri    (:,:) ) 
    338          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_dyn_1d    (1:nidx), wfx_dyn    (:,:) ) 
    339          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_dyn_1d    (1:nidx), hfx_dyn    (:,:) ) 
    340          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_snw_dyn_1d(1:nidx), wfx_snw_dyn(:,:) ) 
     331               CALL tab_1d_2d( npti, nptidx(1:npti), ze_i_2d(1:npti,jk,jl), e_i(:,:,jk,jl) ) 
     332            END DO 
     333         END DO 
     334         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_dyn_1d    (1:npti), sfx_dyn    (:,:) ) 
     335         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_bri_1d    (1:npti), sfx_bri    (:,:) ) 
     336         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_dyn_1d    (1:npti), wfx_dyn    (:,:) ) 
     337         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_dyn_1d    (1:npti), hfx_dyn    (:,:) ) 
     338         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_dyn_1d(1:npti), wfx_snw_dyn(:,:) ) 
    341339         IF ( nn_pnd_scheme > 0 ) THEN 
    342             CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_pnd_1d(1:nidx), wfx_pnd(:,:) ) 
     340            CALL tab_1d_2d( npti, nptidx(1:npti), wfx_pnd_1d(1:npti), wfx_pnd(:,:) ) 
    343341         ENDIF 
    344342 
    345       ENDIF ! nidx > 0 
     343      ENDIF ! npti > 0 
    346344    
    347345      CALL ice_var_agg( 1 )  
     
    378376 
    379377      !                       ! Ice thickness needed for rafting 
    380       WHERE( pa_i(1:nidx,:) > 0._wp )   ;   zhi(1:nidx,:) = pv_i(1:nidx,:) / pa_i(1:nidx,:) 
    381       ELSEWHERE                         ;   zhi(1:nidx,:) = 0._wp 
     378      WHERE( pa_i(1:npti,:) > 0._wp )   ;   zhi(1:npti,:) = pv_i(1:npti,:) / pa_i(1:npti,:) 
     379      ELSEWHERE                         ;   zhi(1:npti,:) = 0._wp 
    382380      END WHERE 
    383381 
     
    398396      ! Compute total area of ice plus open water. 
    399397      ! This is in general not equal to one because of divergence during transport 
    400       zasum(1:nidx) = pato_i(1:nidx) + SUM( pa_i(1:nidx,:), dim=2 ) 
    401       ! 
    402       WHERE( zasum(1:nidx) > 0._wp )   ;   z1_asum(1:nidx) = 1._wp / zasum(1:nidx) 
    403       ELSEWHERE                        ;   z1_asum(1:nidx) = 0._wp 
     398      zasum(1:npti) = pato_i(1:npti) + SUM( pa_i(1:npti,:), dim=2 ) 
     399      ! 
     400      WHERE( zasum(1:npti) > 0._wp )   ;   z1_asum(1:npti) = 1._wp / zasum(1:npti) 
     401      ELSEWHERE                        ;   z1_asum(1:npti) = 0._wp 
    404402      END WHERE 
    405403      ! Compute cumulative thickness distribution function 
     
    407405      ! where zGsum(n) is the fractional area in categories 0 to n. 
    408406      ! initial value (in h = 0) equals open water area 
    409       zGsum(1:nidx,-1) = 0._wp 
    410       zGsum(1:nidx,0 ) = pato_i(1:nidx) * z1_asum(1:nidx) 
     407      zGsum(1:npti,-1) = 0._wp 
     408      zGsum(1:npti,0 ) = pato_i(1:npti) * z1_asum(1:npti) 
    411409      DO jl = 1, jpl 
    412          zGsum(1:nidx,jl) = ( pato_i(1:nidx) + SUM( pa_i(1:nidx,1:jl), dim=2 ) ) * z1_asum(1:nidx)  ! sum(1:jl) is ok (and not jpl) 
     410         zGsum(1:npti,jl) = ( pato_i(1:npti) + SUM( pa_i(1:npti,1:jl), dim=2 ) ) * z1_asum(1:npti)  ! sum(1:jl) is ok (and not jpl) 
    413411      END DO 
    414412      ! 
    415413      IF( ln_partf_lin ) THEN          !--- Linear formulation (Thorndike et al., 1975) 
    416414         DO jl = 0, jpl     
    417             DO ji = 1, nidx 
     415            DO ji = 1, npti 
    418416               IF    ( zGsum(ji,jl)   < rn_gstar ) THEN 
    419417                  apartf(ji,jl) = z1_gstar * ( zGsum(ji,jl) - zGsum(ji,jl-1) ) * & 
     
    432430         zfac = 1._wp / ( 1._wp - EXP(-z1_astar) ) 
    433431         DO jl = -1, jpl 
    434             DO ji = 1, nidx 
     432            DO ji = 1, npti 
    435433               zGsum(ji,jl) = EXP( -zGsum(ji,jl) * z1_astar ) * zfac 
    436434            END DO 
    437435         END DO 
    438436         DO jl = 0, jpl 
    439             DO ji = 1, nidx 
     437            DO ji = 1, npti 
    440438               apartf(ji,jl) = zGsum(ji,jl-1) - zGsum(ji,jl) 
    441439            END DO 
     
    447445      IF( ln_rafting .AND. ln_ridging ) THEN             !- ridging & rafting 
    448446         DO jl = 1, jpl 
    449             DO ji = 1, nidx 
     447            DO ji = 1, npti 
    450448               aridge(ji,jl) = ( 1._wp + TANH ( rn_craft * ( zhi(ji,jl) - rn_hraft ) ) ) * 0.5_wp * apartf(ji,jl) 
    451449               araft (ji,jl) = apartf(ji,jl) - aridge(ji,jl) 
     
    454452      ELSEIF( ln_ridging .AND. .NOT. ln_rafting ) THEN   !- ridging alone 
    455453         DO jl = 1, jpl 
    456             DO ji = 1, nidx 
     454            DO ji = 1, npti 
    457455               aridge(ji,jl) = apartf(ji,jl) 
    458456               araft (ji,jl) = 0._wp 
     
    461459      ELSEIF( ln_rafting .AND. .NOT. ln_ridging ) THEN   !- rafting alone    
    462460         DO jl = 1, jpl 
    463             DO ji = 1, nidx 
     461            DO ji = 1, npti 
    464462               aridge(ji,jl) = 0._wp 
    465463               araft (ji,jl) = apartf(ji,jl) 
     
    468466      ELSE                                               !- no ridging & no rafting 
    469467         DO jl = 1, jpl 
    470             DO ji = 1, nidx 
     468            DO ji = 1, npti 
    471469               aridge(ji,jl) = 0._wp 
    472470               araft (ji,jl) = 0._wp          
     
    502500 
    503501      zfac = 1._wp / hi_hrft 
    504       zaksum(1:nidx) = apartf(1:nidx,0) 
     502      zaksum(1:npti) = apartf(1:npti,0) 
    505503      ! Transfer function 
    506504      DO jl = 1, jpl !all categories have a specific transfer function 
    507          DO ji = 1, nidx 
     505         DO ji = 1, npti 
    508506            IF ( apartf(ji,jl) > 0._wp ) THEN 
    509507               zhmean         = MAX( SQRT( rn_hstar * zhi(ji,jl) ), zhi(ji,jl) * hrdg_hi_min ) 
     
    530528      !  closing rate to a gross closing rate.   
    531529      ! NOTE: 0 < aksum <= 1 
    532       WHERE( zaksum(1:nidx) > 0._wp )   ;   closing_gross(1:nidx) = closing_net(1:nidx) / zaksum(1:nidx) 
    533       ELSEWHERE                         ;   closing_gross(1:nidx) = 0._wp 
     530      WHERE( zaksum(1:npti) > 0._wp )   ;   closing_gross(1:npti) = closing_net(1:npti) / zaksum(1:npti) 
     531      ELSEWHERE                         ;   closing_gross(1:npti) = 0._wp 
    534532      END WHERE 
    535533       
    536       DO ji = 1, nidx   
     534      DO ji = 1, npti   
    537535         ! correction to closing rate and opening if closing rate is excessive 
    538536         !--------------------------------------------------------------------- 
     
    552550      ! would be removed.  Reduce the opening rate proportionately. 
    553551      DO jl = 1, jpl 
    554          DO ji = 1, nidx 
     552         DO ji = 1, npti 
    555553            zfac = apartf(ji,jl) * closing_gross(ji) * rdt_ice 
    556554!!            IF( zfac > pa_i(ji,jl) .AND. zfac > 0._wp ) THEN 
     
    594592      ! 1) Compute change in open water area due to closing and opening. 
    595593      !------------------------------------------------------------------------------- 
    596       DO ji = 1, nidx 
     594      DO ji = 1, npti 
    597595         ato_i_1d(ji) = MAX( 0._wp, ato_i_1d(ji) + ( opning(ji) - apartf(ji,0) * closing_gross(ji) ) * rdt_ice ) 
    598596      END DO 
     
    603601      DO jl1 = 1, jpl 
    604602 
    605          CALL tab_2d_1d( nidx, idxice(1:nidx), s_i_1d(1:nidx), s_i(:,:,jl1) ) 
    606  
    607          DO ji = 1, nidx 
     603         CALL tab_2d_1d( npti, nptidx(1:npti), s_i_1d(1:npti), s_i(:,:,jl1) ) 
     604 
     605         DO ji = 1, npti 
    608606 
    609607            !------------------------------------------------ 
     
    716714         ! special loop for e_i because of layers jk 
    717715         DO jk = 1, nlay_i 
    718             DO ji = 1, nidx 
     716            DO ji = 1, npti 
    719717               IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN 
    720718                  ! Compute ridging /rafting fractions 
     
    736734         DO jl2  = 1, jpl  
    737735            ! 
    738             DO ji = 1, nidx 
     736            DO ji = 1, npti 
    739737 
    740738               IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN 
     
    783781 
    784782            DO jk = 1, nlay_i 
    785                DO ji = 1, nidx 
     783               DO ji = 1, npti 
    786784                  IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp )   & 
    787785                     &   ze_i_2d(ji,jk,jl2) = ze_i_2d(ji,jk,jl2) + eirdg(ji,jk) * fvol(ji) + eirft(ji,jk) * zswitch(ji)                   
     
    794792      ! 
    795793      ! In case ridging/rafting lead to very small negative values (sometimes it happens) 
    796       WHERE( a_i_2d(1:nidx,:) < 0._wp )   a_i_2d(1:nidx,:) = 0._wp 
    797       WHERE( v_i_2d(1:nidx,:) < 0._wp )   v_i_2d(1:nidx,:) = 0._wp 
     794      WHERE( a_i_2d(1:npti,:) < 0._wp )   a_i_2d(1:npti,:) = 0._wp 
     795      WHERE( v_i_2d(1:npti,:) < 0._wp )   v_i_2d(1:npti,:) = 0._wp 
    798796      ! 
    799797   END SUBROUTINE rdgrft_shift 
     
    820818      !!---------------------------------------------------------------------- 
    821819      !                              !--------------------------------------------------! 
    822       IF( ln_str_R75 ) THEN          ! Ice strength => Rothrock (1975) method           ! 
    823       !                              !--------------------------------------------------! 
    824       !                              !--------------------------------------------------! 
    825       ELSEIF( ln_str_H79 ) THEN      ! Ice strength => Hibler (1979) method             ! 
     820      IF( ln_str_H79 ) THEN          ! Ice strength => Hibler (1979) method             ! 
    826821      !                              !--------------------------------------------------! 
    827822         strength(:,:) = rn_pstar * SUM( v_i(:,:,:), dim=3 ) * EXP( -rn_crhg * ( 1._wp - SUM( a_i(:,:,:), dim=3 ) ) ) 
    828          ! 
    829823         ismooth = 1 
    830          ! 
     824         !                           !--------------------------------------------------! 
     825      ELSE                           ! Zero strength                                    ! 
     826         !                           !--------------------------------------------------! 
     827         strength(:,:) = 0._wp 
     828         ismooth = 0 
    831829      ENDIF 
    832830      !                              !--------------------------------------------------! 
     
    896894      !! 
    897895      NAMELIST/namdyn_rdgrft/ ln_str_H79, rn_pstar, rn_crhg, & 
    898          &                    ln_str_R75, rn_perdg,          & 
    899896         &                    rn_csrdg  ,                    & 
    900897         &                    ln_partf_lin, rn_gstar,        & 
     
    921918         WRITE(numout,*) '            1st bulk-rheology parameter                        rn_pstar     = ', rn_pstar 
    922919         WRITE(numout,*) '            2nd bulk-rhelogy parameter                         rn_crhg      = ', rn_crhg 
    923          WRITE(numout,*) '      ice strength parameterization Rothrock (1975)            ln_str_R75   = ', ln_str_R75  
    924          WRITE(numout,*) '            Ratio of ridging work to PotEner change in ridging rn_perdg     = ', rn_perdg  
    925920         WRITE(numout,*) '      Fraction of shear energy contributing to ridging         rn_csrdg     = ', rn_csrdg  
    926921         WRITE(numout,*) '      linear ridging participation function                    ln_partf_lin = ', ln_partf_lin 
     
    940935      ENDIF 
    941936      ! 
    942       IF ( ( ln_str_H79 .AND. ln_str_R75 ) .OR. ( .NOT.ln_str_H79 .AND. .NOT.ln_str_R75 ) ) THEN 
    943          CALL ctl_stop( 'ice_dyn_rdgrft_init: choose one and only one formulation for ice strength (ln_str_H79 or ln_str_R75)' ) 
    944       ENDIF 
    945       ! 
    946937      IF ( ( ln_partf_lin .AND. ln_partf_exp ) .OR. ( .NOT.ln_partf_lin .AND. .NOT.ln_partf_exp ) ) THEN 
    947938         CALL ctl_stop( 'ice_dyn_rdgrft_init: choose one and only one participation function (ln_partf_lin or ln_partf_exp)' ) 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceitd.F90

    r8564 r8565  
    6060      ! 
    6161      INTEGER  ::   ji, jj, jl, jcat     ! dummy loop index 
    62       INTEGER  ::   nidx2                ! local integer 
     62      INTEGER  ::   ipti                 ! local integer 
    6363      REAL(wp) ::   zx1, zwk1, zdh0, zetamin, zdamax   ! local scalars 
    6464      REAL(wp) ::   zx2, zwk2, zda0, zetamax           !   -      - 
     
    6666      REAL(wp) ::   zslope          ! used to compute local thermodynamic "speeds" 
    6767 
    68       INTEGER , DIMENSION(jpij)       ::   idxice2         ! compute remapping or not 
     68      INTEGER , DIMENSION(jpij)       ::   iptidx          ! compute remapping or not 
    6969      INTEGER , DIMENSION(jpij,jpl-1) ::   jdonor          ! donor category index 
    7070      REAL(wp), DIMENSION(jpij,jpl)   ::   zdhice          ! ice thickness increment 
     
    8383      !  1) Identify grid cells with ice 
    8484      !----------------------------------------------------------------------------------------------- 
    85       nidx = 0   ;   idxice(:) = 0 
     85      npti = 0   ;   nptidx(:) = 0 
    8686      DO jj = 1, jpj 
    8787         DO ji = 1, jpi 
    8888            IF ( at_i(ji,jj) > epsi10 ) THEN 
    89                nidx         = nidx + 1 
    90                idxice( nidx ) = (jj - 1) * jpi + ji 
     89               npti         = npti + 1 
     90               nptidx( npti ) = (jj - 1) * jpi + ji 
    9191            ENDIF 
    9292         END DO 
     
    9696      !  2) Compute new category boundaries 
    9797      !----------------------------------------------------------------------------------------------- 
    98       IF( nidx > 0 ) THEN 
     98      IF( npti > 0 ) THEN 
    9999 
    100100         zdhice(:,:) = 0._wp 
    101101         zhbnew(:,:) = 0._wp 
    102102 
    103          CALL tab_3d_2d( nidx, idxice(1:nidx), h_i_2d (1:nidx,1:jpl), h_i   ) 
    104          CALL tab_3d_2d( nidx, idxice(1:nidx), h_ib_2d(1:nidx,1:jpl), h_i_b ) 
    105          CALL tab_3d_2d( nidx, idxice(1:nidx), a_i_2d  (1:nidx,1:jpl), a_i    ) 
    106          CALL tab_3d_2d( nidx, idxice(1:nidx), a_ib_2d (1:nidx,1:jpl), a_i_b  ) 
     103         CALL tab_3d_2d( npti, nptidx(1:npti), h_i_2d (1:npti,1:jpl), h_i   ) 
     104         CALL tab_3d_2d( npti, nptidx(1:npti), h_ib_2d(1:npti,1:jpl), h_i_b ) 
     105         CALL tab_3d_2d( npti, nptidx(1:npti), a_i_2d  (1:npti,1:jpl), a_i    ) 
     106         CALL tab_3d_2d( npti, nptidx(1:npti), a_ib_2d (1:npti,1:jpl), a_i_b  ) 
    107107 
    108108         DO jl = 1, jpl 
    109109            ! Compute thickness change in each ice category 
    110             DO ji = 1, nidx 
     110            DO ji = 1, npti 
    111111               zdhice(ji,jl) = h_i_2d(ji,jl) - h_ib_2d(ji,jl) 
    112112            END DO 
     
    116116         DO jl = 1, jpl - 1 
    117117            ! 
    118             DO ji = 1, nidx 
     118            DO ji = 1, npti 
    119119               ! 
    120120               ! --- New boundary: Hn* = Hn + Fn*dt --- ! 
     
    136136               !    Note: hn(t+1) must not be too close to either HR or HL otherwise a division by nearly 0 is possible  
    137137               !          in ice_itd_glinear in the case (HR-HL) = 3(Hice - HL) or = 3(HR - Hice) 
    138                IF( a_i_2d(ji,jl  ) > epsi10 .AND. h_i_2d(ji,jl  ) > ( zhbnew(ji,jl) - epsi10 ) )   idxice(ji) = 0 
    139                IF( a_i_2d(ji,jl+1) > epsi10 .AND. h_i_2d(ji,jl+1) < ( zhbnew(ji,jl) + epsi10 ) )   idxice(ji) = 0 
     138               IF( a_i_2d(ji,jl  ) > epsi10 .AND. h_i_2d(ji,jl  ) > ( zhbnew(ji,jl) - epsi10 ) )   nptidx(ji) = 0 
     139               IF( a_i_2d(ji,jl+1) > epsi10 .AND. h_i_2d(ji,jl+1) < ( zhbnew(ji,jl) + epsi10 ) )   nptidx(ji) = 0 
    140140                
    141141               ! 2) Hn-1 < Hn* < Hn+1   
    142                IF( zhbnew(ji,jl) < hi_max(jl-1) )   idxice(ji) = 0 
    143                IF( zhbnew(ji,jl) > hi_max(jl+1) )   idxice(ji) = 0 
     142               IF( zhbnew(ji,jl) < hi_max(jl-1) )   nptidx(ji) = 0 
     143               IF( zhbnew(ji,jl) > hi_max(jl+1) )   nptidx(ji) = 0 
    144144                
    145145            END DO 
     
    147147         ! 
    148148         ! --- New boundaries for category jpl --- ! 
    149          DO ji = 1, nidx 
     149         DO ji = 1, npti 
    150150            IF( a_i_2d(ji,jpl) > epsi10 ) THEN 
    151151               zhbnew(ji,jpl) = MAX( hi_max(jpl-1), 3._wp * h_i_2d(ji,jpl) - 2._wp * zhbnew(ji,jpl-1) ) 
     
    158158            !    h1(t) must not be too close to either HR or HL otherwise a division by nearly 0 is possible  
    159159            !    in ice_itd_glinear in the case (HR-HL) = 3(Hice - HL) or = 3(HR - Hice) 
    160             IF( h_ib_2d(ji,1) < ( hi_max(0) + epsi10 ) )   idxice(ji) = 0 
    161             IF( h_ib_2d(ji,1) > ( hi_max(1) - epsi10 ) )   idxice(ji) = 0 
     160            IF( h_ib_2d(ji,1) < ( hi_max(0) + epsi10 ) )   nptidx(ji) = 0 
     161            IF( h_ib_2d(ji,1) > ( hi_max(1) - epsi10 ) )   nptidx(ji) = 0 
    162162         END DO 
    163163         ! 
     
    165165         !  3) Identify cells where remapping 
    166166         !----------------------------------------------------------------------------------------------- 
    167          nidx2 = 0 ; idxice2(:) = 0 
    168          DO ji = 1, nidx 
    169             IF( idxice(ji) /= 0 ) THEN 
    170                nidx2 = nidx2 + 1 
    171                idxice2(nidx2) = idxice(ji) 
    172                zhbnew(nidx2,:) = zhbnew(ji,:)  ! adjust zhbnew to new indices 
     167         ipti = 0 ; iptidx(:) = 0 
     168         DO ji = 1, npti 
     169            IF( nptidx(ji) /= 0 ) THEN 
     170               ipti = ipti + 1 
     171               iptidx(ipti)   = nptidx(ji) 
     172               zhbnew(ipti,:) = zhbnew(ji,:)  ! adjust zhbnew to new indices 
    173173            ENDIF 
    174174         END DO 
    175          idxice(:) = idxice2(:) 
    176          nidx      = nidx2 
     175         nptidx(:) = iptidx(:) 
     176         npti      = ipti 
    177177         ! 
    178178      ENDIF 
     
    181181      !  4) Compute g(h)  
    182182      !----------------------------------------------------------------------------------------------- 
    183       IF( nidx > 0 ) THEN 
     183      IF( npti > 0 ) THEN 
    184184         ! 
    185185         zhb0(:) = hi_max(0)   ;   zhb1(:) = hi_max(1) 
     
    189189         DO jl = 1, jpl 
    190190            ! 
    191             CALL tab_2d_1d( nidx, idxice(1:nidx), h_ib_1d(1:nidx), h_i_b(:,:,jl) ) 
    192             CALL tab_2d_1d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,jl)   ) 
    193             CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,jl)    ) 
    194             CALL tab_2d_1d( nidx, idxice(1:nidx), v_i_1d  (1:nidx), v_i(:,:,jl)    ) 
     191            CALL tab_2d_1d( npti, nptidx(1:npti), h_ib_1d(1:npti), h_i_b(:,:,jl) ) 
     192            CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti), h_i(:,:,jl)   ) 
     193            CALL tab_2d_1d( npti, nptidx(1:npti), a_i_1d  (1:npti), a_i(:,:,jl)    ) 
     194            CALL tab_2d_1d( npti, nptidx(1:npti), v_i_1d  (1:npti), v_i(:,:,jl)    ) 
    195195            ! 
    196196            IF( jl == 1 ) THEN 
    197197               !   
    198198               ! --- g(h) for category 1 --- ! 
    199                CALL ice_itd_glinear( zhb0(1:nidx)  , zhb1(1:nidx)  , h_ib_1d(1:nidx)  , a_i_1d(1:nidx)  ,  &   ! in 
    200                   &                  g0  (1:nidx,1), g1  (1:nidx,1), hL      (1:nidx,1), hR    (1:nidx,1)   )   ! out 
     199               CALL ice_itd_glinear( zhb0(1:npti)  , zhb1(1:npti)  , h_ib_1d(1:npti)  , a_i_1d(1:npti)  ,  &   ! in 
     200                  &                  g0  (1:npti,1), g1  (1:npti,1), hL      (1:npti,1), hR    (1:npti,1)   )   ! out 
    201201                  ! 
    202202               ! Area lost due to melting of thin ice 
    203                DO ji = 1, nidx 
     203               DO ji = 1, npti 
    204204                  ! 
    205205                  IF( a_i_1d(ji) > epsi10 ) THEN 
     
    233233               END DO 
    234234               ! 
    235                CALL tab_1d_2d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,jl)   ) 
    236                CALL tab_1d_2d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,jl)    ) 
    237                CALL tab_1d_2d( nidx, idxice(1:nidx), v_i_1d  (1:nidx), v_i(:,:,jl)    ) 
     235               CALL tab_1d_2d( npti, nptidx(1:npti), h_i_1d (1:npti), h_i(:,:,jl)   ) 
     236               CALL tab_1d_2d( npti, nptidx(1:npti), a_i_1d  (1:npti), a_i(:,:,jl)    ) 
     237               CALL tab_1d_2d( npti, nptidx(1:npti), v_i_1d  (1:npti), v_i(:,:,jl)    ) 
    238238               ! 
    239239            ENDIF ! jl=1 
    240240            ! 
    241241            ! --- g(h) for each thickness category --- !   
    242             CALL ice_itd_glinear( zhbnew(1:nidx,jl-1), zhbnew(1:nidx,jl), h_i_1d(1:nidx)   , a_i_1d(1:nidx)   ,  &   ! in 
    243                &                  g0    (1:nidx,jl  ), g1    (1:nidx,jl), hL     (1:nidx,jl), hR    (1:nidx,jl)   )   ! out 
     242            CALL ice_itd_glinear( zhbnew(1:npti,jl-1), zhbnew(1:npti,jl), h_i_1d(1:npti)   , a_i_1d(1:npti)   ,  &   ! in 
     243               &                  g0    (1:npti,jl  ), g1    (1:npti,jl), hL     (1:npti,jl), hR    (1:npti,jl)   )   ! out 
    244244            ! 
    245245         END DO 
     
    250250         DO jl = 1, jpl - 1 
    251251            ! 
    252             DO ji = 1, nidx 
     252            DO ji = 1, npti 
    253253               ! 
    254254               ! left and right integration limits in eta space 
     
    281281         ! 6) Shift ice between categories 
    282282         !---------------------------------------------------------------------------------------------- 
    283          CALL ice_itd_shiftice ( jdonor(1:nidx,:), zdaice(1:nidx,:), zdvice(1:nidx,:) ) 
     283         CALL ice_itd_shiftice ( jdonor(1:npti,:), zdaice(1:npti,:), zdvice(1:npti,:) ) 
    284284          
    285285         !---------------------------------------------------------------------------------------------- 
    286286         ! 7) Make sure h_i >= minimum ice thickness hi_min 
    287287         !---------------------------------------------------------------------------------------------- 
    288          CALL tab_2d_1d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,1)   ) 
    289          CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,1)    ) 
    290          CALL tab_2d_1d( nidx, idxice(1:nidx), a_ip_1d (1:nidx), a_ip(:,:,1)   ) 
     288         CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti), h_i(:,:,1)   ) 
     289         CALL tab_2d_1d( npti, nptidx(1:npti), a_i_1d  (1:npti), a_i(:,:,1)    ) 
     290         CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_1d (1:npti), a_ip(:,:,1)   ) 
    291291          
    292          DO ji = 1, nidx 
     292         DO ji = 1, npti 
    293293            IF ( a_i_1d(ji) > epsi10 .AND. h_i_1d(ji) < rn_himin ) THEN 
    294294               a_i_1d (ji) = a_i_1d(ji) * h_i_1d(ji) / rn_himin  
     
    302302         END DO 
    303303         ! 
    304          CALL tab_1d_2d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,1)   ) 
    305          CALL tab_1d_2d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,1)    ) 
    306          CALL tab_1d_2d( nidx, idxice(1:nidx), a_ip_1d (1:nidx), a_ip(:,:,1)   ) 
     304         CALL tab_1d_2d( npti, nptidx(1:npti), h_i_1d (1:npti), h_i(:,:,1)   ) 
     305         CALL tab_1d_2d( npti, nptidx(1:npti), a_i_1d  (1:npti), a_i(:,:,1)    ) 
     306         CALL tab_1d_2d( npti, nptidx(1:npti), a_ip_1d (1:npti), a_ip(:,:,1)   ) 
    307307         ! 
    308308      ENDIF 
     
    338338      z2_3 = 2._wp / 3._wp 
    339339      ! 
    340       DO ji = 1, nidx 
     340      DO ji = 1, npti 
    341341         ! 
    342342         IF( paice(ji) > epsi10  .AND. phice(ji) > 0._wp )  THEN 
     
    392392      !!------------------------------------------------------------------ 
    393393          
    394       CALL tab_3d_2d( nidx, idxice(1:nidx), h_i_2d (1:nidx,1:jpl), h_i  ) 
    395       CALL tab_3d_2d( nidx, idxice(1:nidx), a_i_2d (1:nidx,1:jpl), a_i  ) 
    396       CALL tab_3d_2d( nidx, idxice(1:nidx), v_i_2d (1:nidx,1:jpl), v_i  ) 
    397       CALL tab_3d_2d( nidx, idxice(1:nidx), v_s_2d (1:nidx,1:jpl), v_s  ) 
    398       CALL tab_3d_2d( nidx, idxice(1:nidx), oa_i_2d(1:nidx,1:jpl), oa_i ) 
    399       CALL tab_3d_2d( nidx, idxice(1:nidx), sv_i_2d(1:nidx,1:jpl), sv_i ) 
    400       CALL tab_3d_2d( nidx, idxice(1:nidx), a_ip_2d(1:nidx,1:jpl), a_ip ) 
    401       CALL tab_3d_2d( nidx, idxice(1:nidx), v_ip_2d(1:nidx,1:jpl), v_ip ) 
    402       CALL tab_3d_2d( nidx, idxice(1:nidx), t_su_2d(1:nidx,1:jpl), t_su ) 
     394      CALL tab_3d_2d( npti, nptidx(1:npti), h_i_2d (1:npti,1:jpl), h_i  ) 
     395      CALL tab_3d_2d( npti, nptidx(1:npti), a_i_2d (1:npti,1:jpl), a_i  ) 
     396      CALL tab_3d_2d( npti, nptidx(1:npti), v_i_2d (1:npti,1:jpl), v_i  ) 
     397      CALL tab_3d_2d( npti, nptidx(1:npti), v_s_2d (1:npti,1:jpl), v_s  ) 
     398      CALL tab_3d_2d( npti, nptidx(1:npti), oa_i_2d(1:npti,1:jpl), oa_i ) 
     399      CALL tab_3d_2d( npti, nptidx(1:npti), sv_i_2d(1:npti,1:jpl), sv_i ) 
     400      CALL tab_3d_2d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip ) 
     401      CALL tab_3d_2d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip ) 
     402      CALL tab_3d_2d( npti, nptidx(1:npti), t_su_2d(1:npti,1:jpl), t_su ) 
    403403 
    404404      !---------------------------------------------------------------------------------------------- 
     
    406406      !---------------------------------------------------------------------------------------------- 
    407407      DO jl = 1, jpl 
    408          DO ji = 1, nidx 
     408         DO ji = 1, npti 
    409409            zaTsfn(ji,jl) = a_i_2d(ji,jl) * t_su_2d(ji,jl) 
    410410         END DO 
     
    415415      !------------------------------------------------------------------------------- 
    416416      DO jl = 1, jpl - 1 
    417          DO ji = 1, nidx 
     417         DO ji = 1, npti 
    418418            ! 
    419419            jl1 = kdonor(ji,jl) 
     
    475475         DO jk = 1, nlay_s         !--- Snow heat content 
    476476            ! 
    477             DO ji = 1, nidx 
    478                ii = MOD( idxice(ji) - 1, jpi ) + 1 
    479                ij = ( idxice(ji) - 1 ) / jpi + 1 
     477            DO ji = 1, npti 
     478               ii = MOD( nptidx(ji) - 1, jpi ) + 1 
     479               ij = ( nptidx(ji) - 1 ) / jpi + 1 
    480480               ! 
    481481               jl1 = kdonor(ji,jl) 
     
    494494 
    495495         DO jk = 1, nlay_i         !--- Ice heat content 
    496             DO ji = 1, nidx 
    497                ii = MOD( idxice(ji) - 1, jpi ) + 1 
    498                ij = ( idxice(ji) - 1 ) / jpi + 1 
     496            DO ji = 1, npti 
     497               ii = MOD( nptidx(ji) - 1, jpi ) + 1 
     498               ij = ( nptidx(ji) - 1 ) / jpi + 1 
    499499               ! 
    500500               jl1 = kdonor(ji,jl) 
     
    517517      ! 3) Update ice thickness and temperature 
    518518      !------------------------------------------------------------------------------- 
    519       WHERE( a_i_2d(1:nidx,:) >= epsi20 ) 
    520          h_i_2d(1:nidx,:)  =  v_i_2d(1:nidx,:) / a_i_2d(1:nidx,:)  
    521          t_su_2d(1:nidx,:)  =  zaTsfn(1:nidx,:) / a_i_2d(1:nidx,:)  
     519      WHERE( a_i_2d(1:npti,:) >= epsi20 ) 
     520         h_i_2d(1:npti,:)  =  v_i_2d(1:npti,:) / a_i_2d(1:npti,:)  
     521         t_su_2d(1:npti,:)  =  zaTsfn(1:npti,:) / a_i_2d(1:npti,:)  
    522522      ELSEWHERE 
    523          h_i_2d(1:nidx,:)  = 0._wp 
    524          t_su_2d(1:nidx,:)  = rt0 
     523         h_i_2d(1:npti,:)  = 0._wp 
     524         t_su_2d(1:npti,:)  = rt0 
    525525      END WHERE 
    526526      ! 
    527       CALL tab_2d_3d( nidx, idxice(1:nidx), h_i_2d (1:nidx,1:jpl), h_i  ) 
    528       CALL tab_2d_3d( nidx, idxice(1:nidx), a_i_2d (1:nidx,1:jpl), a_i  ) 
    529       CALL tab_2d_3d( nidx, idxice(1:nidx), v_i_2d (1:nidx,1:jpl), v_i  ) 
    530       CALL tab_2d_3d( nidx, idxice(1:nidx), v_s_2d (1:nidx,1:jpl), v_s  ) 
    531       CALL tab_2d_3d( nidx, idxice(1:nidx), oa_i_2d(1:nidx,1:jpl), oa_i ) 
    532       CALL tab_2d_3d( nidx, idxice(1:nidx), sv_i_2d(1:nidx,1:jpl), sv_i ) 
    533       CALL tab_2d_3d( nidx, idxice(1:nidx), a_ip_2d(1:nidx,1:jpl), a_ip ) 
    534       CALL tab_2d_3d( nidx, idxice(1:nidx), v_ip_2d(1:nidx,1:jpl), v_ip ) 
    535       CALL tab_2d_3d( nidx, idxice(1:nidx), t_su_2d(1:nidx,1:jpl), t_su ) 
     527      CALL tab_2d_3d( npti, nptidx(1:npti), h_i_2d (1:npti,1:jpl), h_i  ) 
     528      CALL tab_2d_3d( npti, nptidx(1:npti), a_i_2d (1:npti,1:jpl), a_i  ) 
     529      CALL tab_2d_3d( npti, nptidx(1:npti), v_i_2d (1:npti,1:jpl), v_i  ) 
     530      CALL tab_2d_3d( npti, nptidx(1:npti), v_s_2d (1:npti,1:jpl), v_s  ) 
     531      CALL tab_2d_3d( npti, nptidx(1:npti), oa_i_2d(1:npti,1:jpl), oa_i ) 
     532      CALL tab_2d_3d( npti, nptidx(1:npti), sv_i_2d(1:npti,1:jpl), sv_i ) 
     533      CALL tab_2d_3d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip ) 
     534      CALL tab_2d_3d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip ) 
     535      CALL tab_2d_3d( npti, nptidx(1:npti), t_su_2d(1:npti,1:jpl), t_su ) 
    536536      ! 
    537537   END SUBROUTINE ice_itd_shiftice 
     
    564564      DO jl = 1, jpl-1        ! identify thicknesses that are too big 
    565565         !                    !--------------------------------------- 
    566          nidx = 0   ;   idxice(:) = 0 
     566         npti = 0   ;   nptidx(:) = 0 
    567567         DO jj = 1, jpj 
    568568            DO ji = 1, jpi 
    569569               IF( a_i(ji,jj,jl) > epsi10 .AND. v_i(ji,jj,jl) > (a_i(ji,jj,jl) * hi_max(jl)) ) THEN 
    570                   nidx = nidx + 1 
    571                   idxice( nidx ) = (jj - 1) * jpi + ji                   
     570                  npti = npti + 1 
     571                  nptidx( npti ) = (jj - 1) * jpi + ji                   
    572572               ENDIF 
    573573            END DO 
    574574         END DO 
    575575         ! 
    576 !!clem   CALL tab_2d_1d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,jl)   ) 
    577          CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,jl)    ) 
    578          CALL tab_2d_1d( nidx, idxice(1:nidx), v_i_1d  (1:nidx), v_i(:,:,jl)    ) 
    579          ! 
    580          DO ji = 1, nidx 
     576!!clem   CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti), h_i(:,:,jl)   ) 
     577         CALL tab_2d_1d( npti, nptidx(1:npti), a_i_1d  (1:npti), a_i(:,:,jl)    ) 
     578         CALL tab_2d_1d( npti, nptidx(1:npti), v_i_1d  (1:npti), v_i(:,:,jl)    ) 
     579         ! 
     580         DO ji = 1, npti 
    581581            jdonor(ji,jl)  = jl  
    582582            ! how much of a_i you send in cat sup is somewhat arbitrary 
     
    592592         END DO 
    593593         ! 
    594          IF( nidx > 0 ) THEN 
    595             CALL ice_itd_shiftice( jdonor(1:nidx,:), zdaice(1:nidx,:), zdvice(1:nidx,:) )  ! Shift jl=>jl+1 
     594         IF( npti > 0 ) THEN 
     595            CALL ice_itd_shiftice( jdonor(1:npti,:), zdaice(1:npti,:), zdvice(1:npti,:) )  ! Shift jl=>jl+1 
    596596            ! Reset shift parameters 
    597             jdonor(1:nidx,jl) = 0 
    598             zdaice(1:nidx,jl) = 0._wp 
    599             zdvice(1:nidx,jl) = 0._wp 
     597            jdonor(1:npti,jl) = 0 
     598            zdaice(1:npti,jl) = 0._wp 
     599            zdvice(1:npti,jl) = 0._wp 
    600600         ENDIF 
    601601         ! 
     
    605605      DO jl = jpl-1, 1, -1    ! Identify thicknesses that are too small 
    606606         !                    !----------------------------------------- 
    607          nidx = 0 ; idxice(:) = 0 
     607         npti = 0 ; nptidx(:) = 0 
    608608         DO jj = 1, jpj 
    609609            DO ji = 1, jpi 
    610610               IF( a_i(ji,jj,jl+1) > epsi10 .AND. v_i(ji,jj,jl+1) <= (a_i(ji,jj,jl+1) * hi_max(jl)) ) THEN 
    611                   nidx = nidx + 1 
    612                   idxice( nidx ) = (jj - 1) * jpi + ji                   
     611                  npti = npti + 1 
     612                  nptidx( npti ) = (jj - 1) * jpi + ji                   
    613613               ENDIF 
    614614            END DO 
    615615         END DO 
    616616         ! 
    617          CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,jl+1)    ) ! jl+1 is ok 
    618          CALL tab_2d_1d( nidx, idxice(1:nidx), v_i_1d  (1:nidx), v_i(:,:,jl+1)    ) ! jl+1 is ok 
    619          DO ji = 1, nidx 
     617         CALL tab_2d_1d( npti, nptidx(1:npti), a_i_1d  (1:npti), a_i(:,:,jl+1)    ) ! jl+1 is ok 
     618         CALL tab_2d_1d( npti, nptidx(1:npti), v_i_1d  (1:npti), v_i(:,:,jl+1)    ) ! jl+1 is ok 
     619         DO ji = 1, npti 
    620620            jdonor(ji,jl) = jl + 1 
    621621            zdaice(ji,jl) = a_i_1d(ji)  
     
    623623         END DO 
    624624         ! 
    625          IF( nidx > 0 ) THEN 
    626             CALL ice_itd_shiftice( jdonor(1:nidx,:), zdaice(1:nidx,:), zdvice(1:nidx,:) )  ! Shift jl+1=>jl 
     625         IF( npti > 0 ) THEN 
     626            CALL ice_itd_shiftice( jdonor(1:npti,:), zdaice(1:npti,:), zdvice(1:npti,:) )  ! Shift jl+1=>jl 
    627627            ! Reset shift parameters 
    628             jdonor(1:nidx,jl) = 0 
    629             zdaice(1:nidx,jl) = 0._wp 
    630             zdvice(1:nidx,jl) = 0._wp 
     628            jdonor(1:npti,jl) = 0 
     629            zdaice(1:npti,jl) = 0._wp 
     630            zdvice(1:npti,jl) = 0._wp 
    631631         ENDIF 
    632632         ! 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd.F90

    r8564 r8565  
    211211 
    212212         ! select ice covered grid points 
    213          nidx = 0 ; idxice(:) = 0 
     213         npti = 0 ; nptidx(:) = 0 
    214214         DO jj = 1, jpj 
    215215            DO ji = 1, jpi 
    216216               IF ( a_i(ji,jj,jl) > epsi10 ) THEN      
    217                   nidx         = nidx  + 1 
    218                   idxice(nidx) = (jj - 1) * jpi + ji 
     217                  npti         = npti  + 1 
     218                  nptidx(npti) = (jj - 1) * jpi + ji 
    219219               ENDIF 
    220220            END DO 
    221221         END DO 
    222222 
    223          IF( lk_mpp )         CALL mpp_ini_ice( nidx , numout ) 
    224  
    225          IF( nidx > 0 ) THEN  ! If there is no ice, do nothing. 
     223         IF( lk_mpp )         CALL mpp_ini_ice( npti , numout ) 
     224 
     225         IF( npti > 0 ) THEN  ! If there is no ice, do nothing. 
    226226            !                                                                 
    227227                              CALL ice_thd_1d2d( jl, 1 )            ! --- Move to 1D arrays --- ! 
    228228            !                                                       ! --- & Change units of e_i, e_s from J/m2 to J/m3 --- ! 
    229229            ! 
    230             s_i_new   (1:nidx) = 0._wp ; dh_s_tot (1:nidx) = 0._wp  ! --- some init --- !  (important to have them here)  
    231             dh_i_surf (1:nidx) = 0._wp ; dh_i_bott(1:nidx) = 0._wp 
    232             dh_snowice(1:nidx) = 0._wp ; dh_i_sub (1:nidx) = 0._wp 
     230            s_i_new   (1:npti) = 0._wp ; dh_s_tot (1:npti) = 0._wp  ! --- some init --- !  (important to have them here)  
     231            dh_i_surf (1:npti) = 0._wp ; dh_i_bott(1:npti) = 0._wp 
     232            dh_snowice(1:npti) = 0._wp ; dh_i_sub (1:npti) = 0._wp 
    233233            ! 
    234234            IF( ln_icedH ) THEN                                     ! --- growing/melting --- ! 
    235235                              CALL ice_thd_zdf                             ! Ice/Snow Temperature profile 
    236236                              CALL ice_thd_dh                              ! Ice/Snow thickness    
    237                               CALL ice_thd_ent( e_i_1d(1:nidx,:) )         ! Ice enthalpy remapping 
     237                              CALL ice_thd_ent( e_i_1d(1:npti,:) )         ! Ice enthalpy remapping 
    238238            ENDIF 
    239239            ! 
     
    294294      ! Recover ice temperature 
    295295      DO jk = 1, nlay_i 
    296          DO ji = 1, nidx 
     296         DO ji = 1, npti 
    297297            ztmelts       = -tmut * sz_i_1d(ji,jk) 
    298298            ! Conversion q(S,T) -> T (second order equation) 
     
    323323      !!----------------------------------------------------------------------- 
    324324      ! 
    325       DO ji = 1, nidx 
     325      DO ji = 1, npti 
    326326         zdh_mel = MIN( 0._wp, dh_i_surf(ji) + dh_i_bott(ji) + dh_snowice(ji) + dh_i_sub(ji) ) 
    327327         IF( zdh_mel < 0._wp .AND. a_i_1d(ji) > 0._wp )  THEN 
     
    360360      CASE( 1 )            !==  from 2D to 1D  ==! 
    361361         !                 !---------------------! 
    362          CALL tab_2d_1d( nidx, idxice(1:nidx), at_i_1d(1:nidx), at_i             ) 
    363          CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d (1:nidx), a_i (:,:,kl)     ) 
    364          CALL tab_2d_1d( nidx, idxice(1:nidx), h_i_1d(1:nidx), h_i(:,:,kl)     ) 
    365          CALL tab_2d_1d( nidx, idxice(1:nidx), h_s_1d(1:nidx), h_s(:,:,kl)     ) 
    366          CALL tab_2d_1d( nidx, idxice(1:nidx), t_su_1d(1:nidx), t_su(:,:,kl)     ) 
    367          CALL tab_2d_1d( nidx, idxice(1:nidx), s_i_1d(1:nidx), s_i(:,:,kl)     ) 
     362         CALL tab_2d_1d( npti, nptidx(1:npti), at_i_1d(1:npti), at_i             ) 
     363         CALL tab_2d_1d( npti, nptidx(1:npti), a_i_1d (1:npti), a_i (:,:,kl)     ) 
     364         CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d(1:npti), h_i(:,:,kl)     ) 
     365         CALL tab_2d_1d( npti, nptidx(1:npti), h_s_1d(1:npti), h_s(:,:,kl)     ) 
     366         CALL tab_2d_1d( npti, nptidx(1:npti), t_su_1d(1:npti), t_su(:,:,kl)     ) 
     367         CALL tab_2d_1d( npti, nptidx(1:npti), s_i_1d(1:npti), s_i(:,:,kl)     ) 
    368368         DO jk = 1, nlay_s 
    369             CALL tab_2d_1d( nidx, idxice(1:nidx), t_s_1d(1:nidx,jk), t_s(:,:,jk,kl)   ) 
    370             CALL tab_2d_1d( nidx, idxice(1:nidx), e_s_1d(1:nidx,jk), e_s(:,:,jk,kl)   ) 
     369            CALL tab_2d_1d( npti, nptidx(1:npti), t_s_1d(1:npti,jk), t_s(:,:,jk,kl)   ) 
     370            CALL tab_2d_1d( npti, nptidx(1:npti), e_s_1d(1:npti,jk), e_s(:,:,jk,kl)   ) 
    371371         END DO 
    372372         DO jk = 1, nlay_i 
    373             CALL tab_2d_1d( nidx, idxice(1:nidx), t_i_1d(1:nidx,jk), t_i(:,:,jk,kl)   ) 
    374             CALL tab_2d_1d( nidx, idxice(1:nidx), e_i_1d(1:nidx,jk), e_i(:,:,jk,kl)   ) 
    375             CALL tab_2d_1d( nidx, idxice(1:nidx), sz_i_1d(1:nidx,jk), sz_i(:,:,jk,kl)   ) 
    376          END DO 
    377          ! 
    378          CALL tab_2d_1d( nidx, idxice(1:nidx), qprec_ice_1d(1:nidx), qprec_ice        ) 
    379          CALL tab_2d_1d( nidx, idxice(1:nidx), qsr_ice_1d  (1:nidx), qsr_ice (:,:,kl) ) 
    380          CALL tab_2d_1d( nidx, idxice(1:nidx), fr1_i0_1d   (1:nidx), fr1_i0           ) 
    381          CALL tab_2d_1d( nidx, idxice(1:nidx), fr2_i0_1d   (1:nidx), fr2_i0           ) 
    382          CALL tab_2d_1d( nidx, idxice(1:nidx), qns_ice_1d  (1:nidx), qns_ice (:,:,kl) ) 
    383          CALL tab_2d_1d( nidx, idxice(1:nidx), ftr_ice_1d  (1:nidx), ftr_ice (:,:,kl) ) 
    384          CALL tab_2d_1d( nidx, idxice(1:nidx), evap_ice_1d (1:nidx), evap_ice(:,:,kl) ) 
    385          CALL tab_2d_1d( nidx, idxice(1:nidx), dqns_ice_1d (1:nidx), dqns_ice(:,:,kl) ) 
    386          CALL tab_2d_1d( nidx, idxice(1:nidx), t_bo_1d     (1:nidx), t_bo             ) 
    387          CALL tab_2d_1d( nidx, idxice(1:nidx), sprecip_1d  (1:nidx), sprecip          )  
    388          CALL tab_2d_1d( nidx, idxice(1:nidx), fhtur_1d    (1:nidx), fhtur            ) 
    389          CALL tab_2d_1d( nidx, idxice(1:nidx), fhld_1d     (1:nidx), fhld             ) 
    390          ! 
    391          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_snw_sni_1d(1:nidx), wfx_snw_sni   ) 
    392          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_snw_sum_1d(1:nidx), wfx_snw_sum   ) 
    393          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_sub_1d    (1:nidx), wfx_sub       ) 
    394          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_snw_sub_1d(1:nidx), wfx_snw_sub   ) 
    395          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_ice_sub_1d(1:nidx), wfx_ice_sub   ) 
    396          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_err_sub_1d(1:nidx), wfx_err_sub   ) 
    397          ! 
    398          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_bog_1d (1:nidx), wfx_bog          ) 
    399          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_bom_1d (1:nidx), wfx_bom          ) 
    400          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_sum_1d (1:nidx), wfx_sum          ) 
    401          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_sni_1d (1:nidx), wfx_sni          ) 
    402          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_res_1d (1:nidx), wfx_res          ) 
    403          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_spr_1d (1:nidx), wfx_spr          ) 
    404          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_lam_1d (1:nidx), wfx_lam          ) 
    405          ! 
    406          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_bog_1d (1:nidx), sfx_bog          ) 
    407          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_bom_1d (1:nidx), sfx_bom          ) 
    408          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_sum_1d (1:nidx), sfx_sum          ) 
    409          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_sni_1d (1:nidx), sfx_sni          ) 
    410          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_bri_1d (1:nidx), sfx_bri          ) 
    411          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_res_1d (1:nidx), sfx_res          ) 
    412          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_sub_1d (1:nidx), sfx_sub          ) 
    413          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_lam_1d (1:nidx), sfx_lam          ) 
    414          ! 
    415          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_thd_1d (1:nidx), hfx_thd          ) 
    416          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_spr_1d (1:nidx), hfx_spr          ) 
    417          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_sum_1d (1:nidx), hfx_sum          ) 
    418          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_bom_1d (1:nidx), hfx_bom          ) 
    419          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_bog_1d (1:nidx), hfx_bog          ) 
    420          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_dif_1d (1:nidx), hfx_dif          ) 
    421          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_opw_1d (1:nidx), hfx_opw          ) 
    422          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_snw_1d (1:nidx), hfx_snw          ) 
    423          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_sub_1d (1:nidx), hfx_sub          ) 
    424          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_res_1d (1:nidx), hfx_res          ) 
    425          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_err_dif_1d(1:nidx), hfx_err_dif   ) 
    426          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_err_rem_1d(1:nidx), hfx_err_rem   ) 
    427          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_out_1d (1:nidx), hfx_out          ) 
     373            CALL tab_2d_1d( npti, nptidx(1:npti), t_i_1d(1:npti,jk), t_i(:,:,jk,kl)   ) 
     374            CALL tab_2d_1d( npti, nptidx(1:npti), e_i_1d(1:npti,jk), e_i(:,:,jk,kl)   ) 
     375            CALL tab_2d_1d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,kl)   ) 
     376         END DO 
     377         ! 
     378         CALL tab_2d_1d( npti, nptidx(1:npti), qprec_ice_1d(1:npti), qprec_ice        ) 
     379         CALL tab_2d_1d( npti, nptidx(1:npti), qsr_ice_1d  (1:npti), qsr_ice (:,:,kl) ) 
     380         CALL tab_2d_1d( npti, nptidx(1:npti), fr1_i0_1d   (1:npti), fr1_i0           ) 
     381         CALL tab_2d_1d( npti, nptidx(1:npti), fr2_i0_1d   (1:npti), fr2_i0           ) 
     382         CALL tab_2d_1d( npti, nptidx(1:npti), qns_ice_1d  (1:npti), qns_ice (:,:,kl) ) 
     383         CALL tab_2d_1d( npti, nptidx(1:npti), ftr_ice_1d  (1:npti), ftr_ice (:,:,kl) ) 
     384         CALL tab_2d_1d( npti, nptidx(1:npti), evap_ice_1d (1:npti), evap_ice(:,:,kl) ) 
     385         CALL tab_2d_1d( npti, nptidx(1:npti), dqns_ice_1d (1:npti), dqns_ice(:,:,kl) ) 
     386         CALL tab_2d_1d( npti, nptidx(1:npti), t_bo_1d     (1:npti), t_bo             ) 
     387         CALL tab_2d_1d( npti, nptidx(1:npti), sprecip_1d  (1:npti), sprecip          )  
     388         CALL tab_2d_1d( npti, nptidx(1:npti), fhtur_1d    (1:npti), fhtur            ) 
     389         CALL tab_2d_1d( npti, nptidx(1:npti), fhld_1d     (1:npti), fhld             ) 
     390         ! 
     391         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_snw_sni_1d(1:npti), wfx_snw_sni   ) 
     392         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_snw_sum_1d(1:npti), wfx_snw_sum   ) 
     393         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_sub_1d    (1:npti), wfx_sub       ) 
     394         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_snw_sub_1d(1:npti), wfx_snw_sub   ) 
     395         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_ice_sub_1d(1:npti), wfx_ice_sub   ) 
     396         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_err_sub_1d(1:npti), wfx_err_sub   ) 
     397         ! 
     398         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_bog_1d (1:npti), wfx_bog          ) 
     399         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_bom_1d (1:npti), wfx_bom          ) 
     400         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_sum_1d (1:npti), wfx_sum          ) 
     401         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_sni_1d (1:npti), wfx_sni          ) 
     402         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_res_1d (1:npti), wfx_res          ) 
     403         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_spr_1d (1:npti), wfx_spr          ) 
     404         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_lam_1d (1:npti), wfx_lam          ) 
     405         ! 
     406         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_bog_1d (1:npti), sfx_bog          ) 
     407         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_bom_1d (1:npti), sfx_bom          ) 
     408         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_sum_1d (1:npti), sfx_sum          ) 
     409         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_sni_1d (1:npti), sfx_sni          ) 
     410         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_bri_1d (1:npti), sfx_bri          ) 
     411         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_res_1d (1:npti), sfx_res          ) 
     412         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_sub_1d (1:npti), sfx_sub          ) 
     413         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_lam_1d (1:npti), sfx_lam          ) 
     414         ! 
     415         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_thd_1d (1:npti), hfx_thd          ) 
     416         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_spr_1d (1:npti), hfx_spr          ) 
     417         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_sum_1d (1:npti), hfx_sum          ) 
     418         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_bom_1d (1:npti), hfx_bom          ) 
     419         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_bog_1d (1:npti), hfx_bog          ) 
     420         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_dif_1d (1:npti), hfx_dif          ) 
     421         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_opw_1d (1:npti), hfx_opw          ) 
     422         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_snw_1d (1:npti), hfx_snw          ) 
     423         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_sub_1d (1:npti), hfx_sub          ) 
     424         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_res_1d (1:npti), hfx_res          ) 
     425         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_err_dif_1d(1:npti), hfx_err_dif   ) 
     426         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_err_rem_1d(1:npti), hfx_err_rem   ) 
     427         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_out_1d (1:npti), hfx_out          ) 
    428428         ! 
    429429         ! SIMIP diagnostics 
    430          CALL tab_2d_1d( nidx, idxice(1:nidx), diag_fc_bo_1d(1:nidx), diag_fc_bo   ) 
    431          CALL tab_2d_1d( nidx, idxice(1:nidx), diag_fc_su_1d(1:nidx), diag_fc_su   ) 
     430         CALL tab_2d_1d( npti, nptidx(1:npti), diag_fc_bo_1d(1:npti), diag_fc_bo   ) 
     431         CALL tab_2d_1d( npti, nptidx(1:npti), diag_fc_su_1d(1:npti), diag_fc_su   ) 
    432432         ! ocean surface fields 
    433          CALL tab_2d_1d( nidx, idxice(1:nidx), sst_1d(1:nidx), sst_m ) 
    434          CALL tab_2d_1d( nidx, idxice(1:nidx), sss_1d(1:nidx), sss_m ) 
     433         CALL tab_2d_1d( npti, nptidx(1:npti), sst_1d(1:npti), sst_m ) 
     434         CALL tab_2d_1d( npti, nptidx(1:npti), sss_1d(1:npti), sss_m ) 
    435435 
    436436         ! --- Change units of e_i, e_s from J/m2 to J/m3 --- ! 
    437437         DO jk = 1, nlay_i 
    438             WHERE( h_i_1d(1:nidx)>0._wp ) e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) / (h_i_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_i 
     438            WHERE( h_i_1d(1:npti)>0._wp ) e_i_1d(1:npti,jk) = e_i_1d(1:npti,jk) / (h_i_1d(1:npti) * a_i_1d(1:npti)) * nlay_i 
    439439         END DO 
    440440         DO jk = 1, nlay_s 
    441             WHERE( h_s_1d(1:nidx)>0._wp ) e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) / (h_s_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_s 
     441            WHERE( h_s_1d(1:npti)>0._wp ) e_s_1d(1:npti,jk) = e_s_1d(1:npti,jk) / (h_s_1d(1:npti) * a_i_1d(1:npti)) * nlay_s 
    442442         END DO 
    443443         ! 
     
    447447         ! --- Change units of e_i, e_s from J/m3 to J/m2 --- ! 
    448448         DO jk = 1, nlay_i 
    449             e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) * h_i_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_i 
     449            e_i_1d(1:npti,jk) = e_i_1d(1:npti,jk) * h_i_1d(1:npti) * a_i_1d(1:npti) * r1_nlay_i 
    450450         END DO 
    451451         DO jk = 1, nlay_s 
    452             e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) * h_s_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_s 
     452            e_s_1d(1:npti,jk) = e_s_1d(1:npti,jk) * h_s_1d(1:npti) * a_i_1d(1:npti) * r1_nlay_s 
    453453         END DO 
    454454         ! 
    455455         ! Change thickness to volume (replaces routine ice_var_eqv2glo) 
    456          v_i_1d(1:nidx)  = h_i_1d(1:nidx) * a_i_1d(1:nidx) 
    457          v_s_1d(1:nidx)  = h_s_1d(1:nidx) * a_i_1d(1:nidx) 
    458          sv_i_1d(1:nidx) = s_i_1d(1:nidx) * v_i_1d(1:nidx) 
     456         v_i_1d(1:npti)  = h_i_1d(1:npti) * a_i_1d(1:npti) 
     457         v_s_1d(1:npti)  = h_s_1d(1:npti) * a_i_1d(1:npti) 
     458         sv_i_1d(1:npti) = s_i_1d(1:npti) * v_i_1d(1:npti) 
    459459          
    460          CALL tab_1d_2d( nidx, idxice(1:nidx), at_i_1d(1:nidx), at_i             ) 
    461          CALL tab_1d_2d( nidx, idxice(1:nidx), a_i_1d (1:nidx), a_i (:,:,kl)     ) 
    462          CALL tab_1d_2d( nidx, idxice(1:nidx), h_i_1d(1:nidx), h_i(:,:,kl)     ) 
    463          CALL tab_1d_2d( nidx, idxice(1:nidx), h_s_1d(1:nidx), h_s(:,:,kl)     ) 
    464          CALL tab_1d_2d( nidx, idxice(1:nidx), t_su_1d(1:nidx), t_su(:,:,kl)     ) 
    465          CALL tab_1d_2d( nidx, idxice(1:nidx), s_i_1d(1:nidx), s_i(:,:,kl)     ) 
     460         CALL tab_1d_2d( npti, nptidx(1:npti), at_i_1d(1:npti), at_i             ) 
     461         CALL tab_1d_2d( npti, nptidx(1:npti), a_i_1d (1:npti), a_i (:,:,kl)     ) 
     462         CALL tab_1d_2d( npti, nptidx(1:npti), h_i_1d(1:npti), h_i(:,:,kl)     ) 
     463         CALL tab_1d_2d( npti, nptidx(1:npti), h_s_1d(1:npti), h_s(:,:,kl)     ) 
     464         CALL tab_1d_2d( npti, nptidx(1:npti), t_su_1d(1:npti), t_su(:,:,kl)     ) 
     465         CALL tab_1d_2d( npti, nptidx(1:npti), s_i_1d(1:npti), s_i(:,:,kl)     ) 
    466466         DO jk = 1, nlay_s 
    467             CALL tab_1d_2d( nidx, idxice(1:nidx), t_s_1d(1:nidx,jk), t_s(:,:,jk,kl) ) 
    468             CALL tab_1d_2d( nidx, idxice(1:nidx), e_s_1d(1:nidx,jk), e_s(:,:,jk,kl) ) 
     467            CALL tab_1d_2d( npti, nptidx(1:npti), t_s_1d(1:npti,jk), t_s(:,:,jk,kl) ) 
     468            CALL tab_1d_2d( npti, nptidx(1:npti), e_s_1d(1:npti,jk), e_s(:,:,jk,kl) ) 
    469469         END DO 
    470470         DO jk = 1, nlay_i 
    471             CALL tab_1d_2d( nidx, idxice(1:nidx), t_i_1d(1:nidx,jk), t_i(:,:,jk,kl) ) 
    472             CALL tab_1d_2d( nidx, idxice(1:nidx), e_i_1d(1:nidx,jk), e_i(:,:,jk,kl) ) 
    473             CALL tab_1d_2d( nidx, idxice(1:nidx), sz_i_1d(1:nidx,jk), sz_i(:,:,jk,kl) ) 
    474          END DO 
    475          ! 
    476          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_snw_sni_1d(1:nidx), wfx_snw_sni ) 
    477          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_snw_sum_1d(1:nidx), wfx_snw_sum ) 
    478          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_sub_1d    (1:nidx), wfx_sub     ) 
    479          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_snw_sub_1d(1:nidx), wfx_snw_sub ) 
    480          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_ice_sub_1d(1:nidx), wfx_ice_sub ) 
    481          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_err_sub_1d(1:nidx), wfx_err_sub ) 
    482          ! 
    483          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_bog_1d (1:nidx), wfx_bog        ) 
    484          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_bom_1d (1:nidx), wfx_bom        ) 
    485          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_sum_1d (1:nidx), wfx_sum        ) 
    486          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_sni_1d (1:nidx), wfx_sni        ) 
    487          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_res_1d (1:nidx), wfx_res        ) 
    488          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_spr_1d (1:nidx), wfx_spr        ) 
    489          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_lam_1d (1:nidx), wfx_lam        ) 
    490          ! 
    491          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_bog_1d (1:nidx), sfx_bog        ) 
    492          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_bom_1d (1:nidx), sfx_bom        ) 
    493          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_sum_1d (1:nidx), sfx_sum        ) 
    494          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_sni_1d (1:nidx), sfx_sni        ) 
    495          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_bri_1d (1:nidx), sfx_bri        ) 
    496          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_res_1d (1:nidx), sfx_res        ) 
    497          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_sub_1d (1:nidx), sfx_sub        ) 
    498          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_lam_1d (1:nidx), sfx_lam        ) 
    499          ! 
    500          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_thd_1d (1:nidx), hfx_thd        ) 
    501          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_spr_1d (1:nidx), hfx_spr        ) 
    502          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_sum_1d (1:nidx), hfx_sum        ) 
    503          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_bom_1d (1:nidx), hfx_bom        ) 
    504          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_bog_1d (1:nidx), hfx_bog        ) 
    505          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_dif_1d (1:nidx), hfx_dif        ) 
    506          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_opw_1d (1:nidx), hfx_opw        ) 
    507          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_snw_1d (1:nidx), hfx_snw        ) 
    508          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_sub_1d (1:nidx), hfx_sub        ) 
    509          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_res_1d (1:nidx), hfx_res        ) 
    510          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_err_dif_1d(1:nidx), hfx_err_dif ) 
    511          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_err_rem_1d(1:nidx), hfx_err_rem ) 
    512          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_out_1d (1:nidx), hfx_out        ) 
    513          ! 
    514          CALL tab_1d_2d( nidx, idxice(1:nidx), qns_ice_1d(1:nidx), qns_ice(:,:,kl) ) 
    515          CALL tab_1d_2d( nidx, idxice(1:nidx), ftr_ice_1d(1:nidx), ftr_ice(:,:,kl) ) 
     471            CALL tab_1d_2d( npti, nptidx(1:npti), t_i_1d(1:npti,jk), t_i(:,:,jk,kl) ) 
     472            CALL tab_1d_2d( npti, nptidx(1:npti), e_i_1d(1:npti,jk), e_i(:,:,jk,kl) ) 
     473            CALL tab_1d_2d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,kl) ) 
     474         END DO 
     475         ! 
     476         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_sni_1d(1:npti), wfx_snw_sni ) 
     477         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_sum_1d(1:npti), wfx_snw_sum ) 
     478         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_sub_1d    (1:npti), wfx_sub     ) 
     479         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_sub_1d(1:npti), wfx_snw_sub ) 
     480         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_ice_sub_1d(1:npti), wfx_ice_sub ) 
     481         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_err_sub_1d(1:npti), wfx_err_sub ) 
     482         ! 
     483         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_bog_1d (1:npti), wfx_bog        ) 
     484         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_bom_1d (1:npti), wfx_bom        ) 
     485         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_sum_1d (1:npti), wfx_sum        ) 
     486         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_sni_1d (1:npti), wfx_sni        ) 
     487         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_res_1d (1:npti), wfx_res        ) 
     488         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_spr_1d (1:npti), wfx_spr        ) 
     489         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_lam_1d (1:npti), wfx_lam        ) 
     490         ! 
     491         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_bog_1d (1:npti), sfx_bog        ) 
     492         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_bom_1d (1:npti), sfx_bom        ) 
     493         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_sum_1d (1:npti), sfx_sum        ) 
     494         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_sni_1d (1:npti), sfx_sni        ) 
     495         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_bri_1d (1:npti), sfx_bri        ) 
     496         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_res_1d (1:npti), sfx_res        ) 
     497         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_sub_1d (1:npti), sfx_sub        ) 
     498         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_lam_1d (1:npti), sfx_lam        ) 
     499         ! 
     500         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_thd_1d (1:npti), hfx_thd        ) 
     501         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_spr_1d (1:npti), hfx_spr        ) 
     502         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_sum_1d (1:npti), hfx_sum        ) 
     503         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_bom_1d (1:npti), hfx_bom        ) 
     504         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_bog_1d (1:npti), hfx_bog        ) 
     505         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_dif_1d (1:npti), hfx_dif        ) 
     506         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_opw_1d (1:npti), hfx_opw        ) 
     507         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_snw_1d (1:npti), hfx_snw        ) 
     508         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_sub_1d (1:npti), hfx_sub        ) 
     509         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_res_1d (1:npti), hfx_res        ) 
     510         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_err_dif_1d(1:npti), hfx_err_dif ) 
     511         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_err_rem_1d(1:npti), hfx_err_rem ) 
     512         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_out_1d (1:npti), hfx_out        ) 
     513         ! 
     514         CALL tab_1d_2d( npti, nptidx(1:npti), qns_ice_1d(1:npti), qns_ice(:,:,kl) ) 
     515         CALL tab_1d_2d( npti, nptidx(1:npti), ftr_ice_1d(1:npti), ftr_ice(:,:,kl) ) 
    516516         ! 
    517517         ! SIMIP diagnostics          
    518          CALL tab_1d_2d( nidx, idxice(1:nidx), t_si_1d      (1:nidx), t_si(:,:,kl) ) 
    519          CALL tab_1d_2d( nidx, idxice(1:nidx), diag_fc_bo_1d(1:nidx), diag_fc_bo   ) 
    520          CALL tab_1d_2d( nidx, idxice(1:nidx), diag_fc_su_1d(1:nidx), diag_fc_su   ) 
     518         CALL tab_1d_2d( npti, nptidx(1:npti), t_si_1d      (1:npti), t_si(:,:,kl) ) 
     519         CALL tab_1d_2d( npti, nptidx(1:npti), diag_fc_bo_1d(1:npti), diag_fc_bo   ) 
     520         CALL tab_1d_2d( npti, nptidx(1:npti), diag_fc_su_1d(1:npti), diag_fc_su   ) 
    521521         ! extensive variables 
    522          CALL tab_1d_2d( nidx, idxice(1:nidx), v_i_1d (1:nidx), v_i (:,:,kl) ) 
    523          CALL tab_1d_2d( nidx, idxice(1:nidx), v_s_1d (1:nidx), v_s (:,:,kl) ) 
    524          CALL tab_1d_2d( nidx, idxice(1:nidx), sv_i_1d(1:nidx), sv_i(:,:,kl) ) 
     522         CALL tab_1d_2d( npti, nptidx(1:npti), v_i_1d (1:npti), v_i (:,:,kl) ) 
     523         CALL tab_1d_2d( npti, nptidx(1:npti), v_s_1d (1:npti), v_s (:,:,kl) ) 
     524         CALL tab_1d_2d( npti, nptidx(1:npti), sv_i_1d(1:npti), sv_i(:,:,kl) ) 
    525525         ! 
    526526      END SELECT 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_da.F90

    r8564 r8565  
    118118      ! 
    119119      zastar = 1._wp / ( 1._wp - (rn_dmin / zdmax)**(1._wp/rn_beta) ) 
    120       DO ji = 1, nidx    
     120      DO ji = 1, npti    
    121121         ! --- Calculate reduction of total sea ice concentration --- ! 
    122122         zdfloe = rn_dmin * ( zastar / ( zastar - at_i_1d(ji) ) )**rn_beta         ! Mean floe caliper diameter [m] 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_dh.F90

    r8564 r8565  
    117117 
    118118      ! initialize layer thicknesses and enthalpies 
    119       h_i_old (1:nidx,0:nlay_i+1) = 0._wp 
    120       eh_i_old(1:nidx,0:nlay_i+1) = 0._wp 
     119      h_i_old (1:npti,0:nlay_i+1) = 0._wp 
     120      eh_i_old(1:npti,0:nlay_i+1) = 0._wp 
    121121      DO jk = 1, nlay_i 
    122          DO ji = 1, nidx 
     122         DO ji = 1, npti 
    123123            h_i_old (ji,jk) = h_i_1d(ji) * r1_nlay_i 
    124124            eh_i_old(ji,jk) = e_i_1d(ji,jk) * h_i_old(ji,jk) 
     
    130130      !------------------------------------------------------------------------------! 
    131131      ! 
    132       DO ji = 1, nidx 
     132      DO ji = 1, npti 
    133133         zdum       = qns_ice_1d(ji) + ( 1._wp - i0(ji) ) * qsr_ice_1d(ji) - fc_su(ji)  
    134134         zf_tt(ji)  = fc_bo_i(ji) + fhtur_1d(ji) + fhld_1d(ji)  
     
    142142      ! (should not happen but sometimes it does) 
    143143      !------------------------------------------------------------------------------! 
    144       DO ji = 1, nidx 
     144      DO ji = 1, npti 
    145145         IF( t_s_1d(ji,1) > rt0 ) THEN !!! Internal melting 
    146146            ! Contribution to heat flux to the ocean [W.m-2], < 0   
     
    159159      !------------------------------------------------------------! 
    160160      DO jk = 1, nlay_i 
    161          DO ji = 1, nidx 
     161         DO ji = 1, npti 
    162162            zh_i(ji,jk) = h_i_1d(ji) * r1_nlay_i 
    163163            zeh_i(ji)   = zeh_i(ji) + e_i_1d(ji,jk) * zh_i(ji,jk) 
     
    183183      ! Martin Vancoppenolle, December 2006 
    184184 
    185       CALL ice_thd_snwblow( 1. - at_i_1d(1:nidx), zsnw(1:nidx) ) ! snow distribution over ice after wind blowing 
    186  
    187       zdeltah(1:nidx,:) = 0._wp 
    188       DO ji = 1, nidx 
     185      CALL ice_thd_snwblow( 1. - at_i_1d(1:npti), zsnw(1:npti) ) ! snow distribution over ice after wind blowing 
     186 
     187      zdeltah(1:npti,:) = 0._wp 
     188      DO ji = 1, npti 
    189189         !----------- 
    190190         ! Snow fall 
     
    220220 
    221221      ! If heat still available (zq_su > 0), then melt more snow 
    222       zdeltah(1:nidx,:) = 0._wp 
    223       zdh_s_mel(1:nidx) = 0._wp 
     222      zdeltah(1:npti,:) = 0._wp 
     223      zdh_s_mel(1:npti) = 0._wp 
    224224      DO jk = 1, nlay_s 
    225          DO ji = 1, nidx 
     225         DO ji = 1, npti 
    226226            ! thickness change 
    227227            rswitch          = 1._wp - MAX( 0._wp, SIGN( 1._wp, - h_s_1d(ji) ) )  
     
    245245      ! qla_ice is always >=0 (upwards), heat goes to the atmosphere, therefore snow sublimates 
    246246      ! clem comment: not counted in mass/heat exchange in iceupdate.F90 since this is an exchange with atm. (not ocean) 
    247       zdeltah(1:nidx,:) = 0._wp 
    248       DO ji = 1, nidx 
     247      zdeltah(1:npti,:) = 0._wp 
     248      DO ji = 1, npti 
    249249         zdh_s_sub(ji)  = MAX( - h_s_1d(ji) , - evap_ice_1d(ji) * r1_rhosn * rdt_ice ) 
    250250         ! remaining evap in kg.m-2 (used for ice melting later on) 
     
    265265       
    266266      ! --- Update snow diags --- ! 
    267       DO ji = 1, nidx 
     267      DO ji = 1, npti 
    268268         dh_s_tot(ji) = zdh_s_mel(ji) + zdh_s_pre(ji) + zdh_s_sub(ji) 
    269269      END DO 
     
    274274      ! new temp and enthalpy of the snow (remaining snow precip + remaining pre-existing snow) 
    275275      DO jk = 1, nlay_s 
    276          DO ji = 1,nidx 
     276         DO ji = 1,npti 
    277277            rswitch       = MAX( 0._wp , SIGN( 1._wp, h_s_1d(ji) - epsi20 ) ) 
    278278            e_s_1d(ji,jk) = rswitch / MAX( h_s_1d(ji), epsi20 ) *           & 
     
    285285      ! 3.4 Surface ice ablation  
    286286      !-------------------------- 
    287       zdeltah(1:nidx,:) = 0._wp ! important 
     287      zdeltah(1:npti,:) = 0._wp ! important 
    288288      DO jk = 1, nlay_i 
    289          DO ji = 1, nidx 
     289         DO ji = 1, npti 
    290290            ztmelts           = - tmut * sz_i_1d(ji,jk)          ! Melting point of layer k [C] 
    291291             
     
    375375      END DO 
    376376      ! update ice thickness 
    377       DO ji = 1, nidx 
     377      DO ji = 1, npti 
    378378         h_i_1d(ji) =  MAX( 0._wp , h_i_1d(ji) + dh_i_surf(ji) + dh_i_sub(ji) ) 
    379379      END DO 
    380380 
    381381      ! remaining "potential" evap is sent to ocean 
    382       DO ji = 1, nidx 
     382      DO ji = 1, npti 
    383383         wfx_err_sub_1d(ji) = wfx_err_sub_1d(ji) - zevap_rema(ji) * a_i_1d(ji) * r1_rdtice  ! <=0 (net evap for the ocean in kg.m-2.s-1) 
    384384      END DO 
     
    407407       
    408408      ! Iterative procedure 
    409       DO ji = 1, nidx 
     409      DO ji = 1, npti 
    410410         IF(  zf_tt(ji) < 0._wp  ) THEN 
    411411            DO iter = 1, num_iter_max 
     
    478478      ! 4.2 Basal melt 
    479479      !---------------- 
    480       zdeltah(1:nidx,:) = 0._wp ! important 
     480      zdeltah(1:npti,:) = 0._wp ! important 
    481481      DO jk = nlay_i, 1, -1 
    482          DO ji = 1, nidx 
     482         DO ji = 1, npti 
    483483            IF(  zf_tt(ji)  >  0._wp  .AND. jk > icount(ji,jk) ) THEN   ! do not calculate where layer has already disappeared by surface melting  
    484484 
     
    554554      ! Update temperature, energy 
    555555      !------------------------------------------- 
    556       DO ji = 1, nidx 
     556      DO ji = 1, npti 
    557557         h_i_1d(ji) =  MAX( 0._wp , h_i_1d(ji) + dh_i_bott(ji) ) 
    558558      END DO   
     
    563563      ! If heat still available for melting and snow remains, then melt more snow 
    564564      !------------------------------------------- 
    565       zdeltah(1:nidx,:) = 0._wp ! important 
    566       DO ji = 1, nidx 
     565      zdeltah(1:npti,:) = 0._wp ! important 
     566      DO ji = 1, npti 
    567567         zq_rema(ji)     = zq_su(ji) + zq_bo(ji)  
    568568         rswitch         = 1._wp - MAX( 0._wp, SIGN( 1._wp, - h_s_1d(ji) ) )   ! =1 if snow 
     
    592592      ! flooding of seawater transforms snow into ice dh_snowice is positive for the ice 
    593593      z1_rho = 1._wp / ( rhosn+rau0-rhoic ) 
    594       DO ji = 1, nidx 
     594      DO ji = 1, npti 
    595595         ! 
    596596         dh_snowice(ji) = MAX(  0._wp , ( rhosn * h_s_1d(ji) + (rhoic-rau0) * h_i_1d(ji) ) * z1_rho ) 
     
    631631      ! Update temperature, energy 
    632632      !------------------------------------------- 
    633       DO ji = 1, nidx 
     633      DO ji = 1, npti 
    634634         rswitch     =  1.0 - MAX( 0._wp , SIGN( 1._wp , - h_i_1d(ji) ) )  
    635635         t_su_1d(ji) =  rswitch * t_su_1d(ji) + ( 1.0 - rswitch ) * rt0 
     
    637637 
    638638      DO jk = 1, nlay_s 
    639          DO ji = 1,nidx 
     639         DO ji = 1,npti 
    640640            ! mask enthalpy 
    641641            rswitch       = 1._wp - MAX(  0._wp , SIGN( 1._wp, - h_s_1d(ji) )  ) 
     
    647647 
    648648      ! --- ensure that a_i = 0 where h_i = 0 --- 
    649       WHERE( h_i_1d(1:nidx) == 0._wp )   a_i_1d(1:nidx) = 0._wp 
     649      WHERE( h_i_1d(1:npti) == 0._wp )   a_i_1d(1:npti) = 0._wp 
    650650      ! 
    651651   END SUBROUTINE ice_thd_dh 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_do.F90

    r8564 r8565  
    222222      !------------------------------------- 
    223223      ! This occurs if open water energy budget is negative (cooling) and there is no landfast ice 
    224       nidx = 0 ; idxice(:) = 0 
     224      npti = 0 ; nptidx(:) = 0 
    225225      DO jj = 1, jpj 
    226226         DO ji = 1, jpi 
    227227            IF ( qlead(ji,jj)  <  0._wp .AND. tau_icebfr(ji,jj) == 0._wp ) THEN 
    228                nidx = nidx + 1 
    229                idxice( nidx ) = (jj - 1) * jpi + ji 
     228               npti = npti + 1 
     229               nptidx( npti ) = (jj - 1) * jpi + ji 
    230230            ENDIF 
    231231         END DO 
     
    237237      ! If ocean gains heat do nothing. Otherwise compute new ice formation 
    238238 
    239       IF ( nidx > 0 ) THEN 
    240  
    241          CALL tab_2d_1d( nidx, idxice(1:nidx), at_i_1d(1:nidx)      , at_i        ) 
    242          CALL tab_3d_2d( nidx, idxice(1:nidx), a_i_2d (1:nidx,1:jpl), a_i (:,:,:) ) 
    243          CALL tab_3d_2d( nidx, idxice(1:nidx), v_i_2d (1:nidx,1:jpl), v_i (:,:,:) ) 
    244          CALL tab_3d_2d( nidx, idxice(1:nidx), sv_i_2d(1:nidx,1:jpl), sv_i(:,:,:) ) 
     239      IF ( npti > 0 ) THEN 
     240 
     241         CALL tab_2d_1d( npti, nptidx(1:npti), at_i_1d(1:npti)      , at_i        ) 
     242         CALL tab_3d_2d( npti, nptidx(1:npti), a_i_2d (1:npti,1:jpl), a_i (:,:,:) ) 
     243         CALL tab_3d_2d( npti, nptidx(1:npti), v_i_2d (1:npti,1:jpl), v_i (:,:,:) ) 
     244         CALL tab_3d_2d( npti, nptidx(1:npti), sv_i_2d(1:npti,1:jpl), sv_i(:,:,:) ) 
    245245         DO jl = 1, jpl 
    246246            DO jk = 1, nlay_i 
    247                CALL tab_2d_1d( nidx, idxice(1:nidx), ze_i_2d(1:nidx,jk,jl), e_i(:,:,jk,jl) ) 
    248             END DO 
    249          END DO 
    250          CALL tab_2d_1d( nidx, idxice(1:nidx), qlead_1d  (1:nidx) , qlead       ) 
    251          CALL tab_2d_1d( nidx, idxice(1:nidx), t_bo_1d   (1:nidx) , t_bo        ) 
    252          CALL tab_2d_1d( nidx, idxice(1:nidx), sfx_opw_1d(1:nidx) , sfx_opw     ) 
    253          CALL tab_2d_1d( nidx, idxice(1:nidx), wfx_opw_1d(1:nidx) , wfx_opw     ) 
    254          CALL tab_2d_1d( nidx, idxice(1:nidx), zh_newice (1:nidx) , ht_i_new    ) 
    255          CALL tab_2d_1d( nidx, idxice(1:nidx), zvrel_1d  (1:nidx) , zvrel       ) 
    256  
    257          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_thd_1d(1:nidx) , hfx_thd     ) 
    258          CALL tab_2d_1d( nidx, idxice(1:nidx), hfx_opw_1d(1:nidx) , hfx_opw     ) 
    259          CALL tab_2d_1d( nidx, idxice(1:nidx), rn_amax_1d(1:nidx) , rn_amax_2d  ) 
    260          CALL tab_2d_1d( nidx, idxice(1:nidx), sss_1d    (1:nidx) , sss_m       ) 
     247               CALL tab_2d_1d( npti, nptidx(1:npti), ze_i_2d(1:npti,jk,jl), e_i(:,:,jk,jl) ) 
     248            END DO 
     249         END DO 
     250         CALL tab_2d_1d( npti, nptidx(1:npti), qlead_1d  (1:npti) , qlead       ) 
     251         CALL tab_2d_1d( npti, nptidx(1:npti), t_bo_1d   (1:npti) , t_bo        ) 
     252         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_opw_1d(1:npti) , sfx_opw     ) 
     253         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_opw_1d(1:npti) , wfx_opw     ) 
     254         CALL tab_2d_1d( npti, nptidx(1:npti), zh_newice (1:npti) , ht_i_new    ) 
     255         CALL tab_2d_1d( npti, nptidx(1:npti), zvrel_1d  (1:npti) , zvrel       ) 
     256 
     257         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_thd_1d(1:npti) , hfx_thd     ) 
     258         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_opw_1d(1:npti) , hfx_opw     ) 
     259         CALL tab_2d_1d( npti, nptidx(1:npti), rn_amax_1d(1:npti) , rn_amax_2d  ) 
     260         CALL tab_2d_1d( npti, nptidx(1:npti), sss_1d    (1:npti) , sss_m       ) 
    261261 
    262262         !------------------------------------------------------------------------------| 
     
    265265         DO jl = 1, jpl 
    266266            DO jk = 1, nlay_i                
    267                WHERE( v_i_2d(1:nidx,jl) > 0._wp ) 
    268                   ze_i_2d(1:nidx,jk,jl) = ze_i_2d(1:nidx,jk,jl) / v_i_2d(1:nidx,jl) * REAL( nlay_i ) 
     267               WHERE( v_i_2d(1:npti,jl) > 0._wp ) 
     268                  ze_i_2d(1:npti,jk,jl) = ze_i_2d(1:npti,jk,jl) / v_i_2d(1:npti,jl) * REAL( nlay_i ) 
    269269               ELSEWHERE 
    270                   ze_i_2d(1:nidx,jk,jl) = 0._wp 
     270                  ze_i_2d(1:npti,jk,jl) = 0._wp 
    271271               END WHERE 
    272272            END DO 
     
    279279         ! Keep old ice areas and volume in memory 
    280280         !----------------------------------------- 
    281          zv_b(1:nidx,:) = v_i_2d(1:nidx,:)  
    282          za_b(1:nidx,:) = a_i_2d(1:nidx,:) 
     281         zv_b(1:npti,:) = v_i_2d(1:npti,:)  
     282         za_b(1:npti,:) = a_i_2d(1:npti,:) 
    283283 
    284284         !---------------------- 
     
    287287         SELECT CASE ( nn_icesal ) 
    288288         CASE ( 1 )                    ! Sice = constant  
    289             zs_newice(1:nidx) = rn_icesal 
     289            zs_newice(1:npti) = rn_icesal 
    290290         CASE ( 2 )                    ! Sice = F(z,t) [Vancoppenolle et al (2005)] 
    291             DO ji = 1, nidx 
     291            DO ji = 1, npti 
    292292               zs_newice(ji) = MIN(  4.606 + 0.91 / zh_newice(ji) , rn_simax , 0.5 * sss_1d(ji) ) 
    293293            END DO 
    294294         CASE ( 3 )                    ! Sice = F(z) [multiyear ice] 
    295             zs_newice(1:nidx) =   2.3 
     295            zs_newice(1:npti) =   2.3 
    296296         END SELECT 
    297297 
     
    300300         !------------------------- 
    301301         ! We assume that new ice is formed at the seawater freezing point 
    302          DO ji = 1, nidx 
     302         DO ji = 1, npti 
    303303            ztmelts       = - tmut * zs_newice(ji)                  ! Melting point (C) 
    304304            ze_newice(ji) =   rhoic * (  cpic * ( ztmelts - ( t_bo_1d(ji) - rt0 ) )                     & 
     
    310310         ! Age of new ice 
    311311         !---------------- 
    312          zo_newice(1:nidx) = 0._wp 
     312         zo_newice(1:npti) = 0._wp 
    313313 
    314314         !------------------- 
    315315         ! Volume of new ice 
    316316         !------------------- 
    317          DO ji = 1, nidx 
     317         DO ji = 1, npti 
    318318 
    319319            zEi           = - ze_newice(ji) * r1_rhoic             ! specific enthalpy of forming ice [J/kg] 
     
    340340         END DO 
    341341          
    342          zv_frazb(1:nidx) = 0._wp 
     342         zv_frazb(1:npti) = 0._wp 
    343343         IF( ln_frazil ) THEN 
    344344            ! A fraction zfrazb of frazil ice is accreted at the ice bottom 
    345             DO ji = 1, nidx 
     345            DO ji = 1, npti 
    346346               rswitch       = 1._wp - MAX( 0._wp, SIGN( 1._wp , - at_i_1d(ji) ) ) 
    347347               zfrazb        = rswitch * ( TANH( rn_Cfraz * ( zvrel_1d(ji) - rn_vfraz ) ) + 1.0 ) * 0.5 * rn_maxfraz 
     
    354354         ! Area of new ice 
    355355         !----------------- 
    356          DO ji = 1, nidx 
     356         DO ji = 1, npti 
    357357            za_newice(ji) = zv_newice(ji) / zh_newice(ji) 
    358358         END DO 
     
    367367         ! If lateral ice growth gives an ice concentration gt 1, then 
    368368         ! we keep the excessive volume in memory and attribute it later to bottom accretion 
    369          DO ji = 1, nidx 
     369         DO ji = 1, npti 
    370370            IF ( za_newice(ji) >  ( rn_amax_1d(ji) - at_i_1d(ji) ) ) THEN 
    371371               zda_res(ji)   = za_newice(ji) - ( rn_amax_1d(ji) - at_i_1d(ji) ) 
     
    381381         ! find which category to fill 
    382382         DO jl = 1, jpl 
    383             DO ji = 1, nidx 
     383            DO ji = 1, npti 
    384384               IF( zh_newice(ji) > hi_max(jl-1) .AND. zh_newice(ji) <= hi_max(jl) ) THEN 
    385385                  a_i_2d(ji,jl) = a_i_2d(ji,jl) + za_newice(ji) 
     
    389389            END DO 
    390390         END DO 
    391          at_i_1d(1:nidx) = SUM( a_i_2d(1:nidx,:), dim=2 ) 
     391         at_i_1d(1:npti) = SUM( a_i_2d(1:npti,:), dim=2 ) 
    392392 
    393393         ! Heat content 
    394          DO ji = 1, nidx 
     394         DO ji = 1, npti 
    395395            jl = jcat(ji)                                                    ! categroy in which new ice is put 
    396396            zswinew  (ji) = MAX( 0._wp , SIGN( 1._wp , - za_b(ji,jl) ) )   ! 0 if old ice 
     
    398398 
    399399         DO jk = 1, nlay_i 
    400             DO ji = 1, nidx 
     400            DO ji = 1, npti 
    401401               jl = jcat(ji) 
    402402               rswitch = MAX( 0._wp, SIGN( 1._wp , v_i_2d(ji,jl) - epsi20 ) ) 
     
    413413 
    414414            ! for remapping 
    415             h_i_old (1:nidx,0:nlay_i+1) = 0._wp 
    416             eh_i_old(1:nidx,0:nlay_i+1) = 0._wp 
     415            h_i_old (1:npti,0:nlay_i+1) = 0._wp 
     416            eh_i_old(1:npti,0:nlay_i+1) = 0._wp 
    417417            DO jk = 1, nlay_i 
    418                DO ji = 1, nidx 
     418               DO ji = 1, npti 
    419419                  h_i_old (ji,jk) = v_i_2d(ji,jl) * r1_nlay_i 
    420420                  eh_i_old(ji,jk) = ze_i_2d(ji,jk,jl) * h_i_old(ji,jk) 
     
    423423 
    424424            ! new volumes including lateral/bottom accretion + residual 
    425             DO ji = 1, nidx 
     425            DO ji = 1, npti 
    426426               rswitch        = MAX( 0._wp, SIGN( 1._wp , at_i_1d(ji) - epsi20 ) ) 
    427427               zv_newfra     = rswitch * ( zdv_res(ji) + zv_frazb(ji) ) * a_i_2d(ji,jl) / MAX( at_i_1d(ji) , epsi20 ) 
     
    433433            ENDDO 
    434434            ! --- Ice enthalpy remapping --- ! 
    435             CALL ice_thd_ent( ze_i_2d(1:nidx,:,jl) )  
     435            CALL ice_thd_ent( ze_i_2d(1:npti,:,jl) )  
    436436         ENDDO 
    437437 
     
    440440         !----------------- 
    441441         DO jl = 1, jpl 
    442             DO ji = 1, nidx 
     442            DO ji = 1, npti 
    443443               sv_i_2d(ji,jl) = sv_i_2d(ji,jl) + zs_newice(ji) * ( v_i_2d(ji,jl) - zv_b(ji,jl) ) 
    444444            END DO 
     
    450450         DO jl = 1, jpl 
    451451            DO jk = 1, nlay_i 
    452                ze_i_2d(1:nidx,jk,jl) = ze_i_2d(1:nidx,jk,jl) * v_i_2d(1:nidx,jl) * r1_nlay_i  
     452               ze_i_2d(1:npti,jk,jl) = ze_i_2d(1:npti,jk,jl) * v_i_2d(1:npti,jl) * r1_nlay_i  
    453453            END DO 
    454454         END DO 
     
    456456         ! 7) Change 2D vectors to 1D vectors  
    457457         !------------------------------------------------------------------------------! 
    458          CALL tab_2d_3d( nidx, idxice(1:nidx), a_i_2d (1:nidx,1:jpl), a_i (:,:,:) ) 
    459          CALL tab_2d_3d( nidx, idxice(1:nidx), v_i_2d (1:nidx,1:jpl), v_i (:,:,:) ) 
    460          CALL tab_2d_3d( nidx, idxice(1:nidx), sv_i_2d(1:nidx,1:jpl), sv_i(:,:,:) ) 
     458         CALL tab_2d_3d( npti, nptidx(1:npti), a_i_2d (1:npti,1:jpl), a_i (:,:,:) ) 
     459         CALL tab_2d_3d( npti, nptidx(1:npti), v_i_2d (1:npti,1:jpl), v_i (:,:,:) ) 
     460         CALL tab_2d_3d( npti, nptidx(1:npti), sv_i_2d(1:npti,1:jpl), sv_i(:,:,:) ) 
    461461          DO jl = 1, jpl 
    462462            DO jk = 1, nlay_i 
    463                CALL tab_1d_2d( nidx, idxice(1:nidx), ze_i_2d(1:nidx,jk,jl), e_i(:,:,jk,jl) ) 
    464             END DO 
    465          END DO 
    466          CALL tab_1d_2d( nidx, idxice(1:nidx), sfx_opw_1d(1:nidx), sfx_opw ) 
    467          CALL tab_1d_2d( nidx, idxice(1:nidx), wfx_opw_1d(1:nidx), wfx_opw ) 
    468          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_thd_1d(1:nidx), hfx_thd ) 
    469          CALL tab_1d_2d( nidx, idxice(1:nidx), hfx_opw_1d(1:nidx), hfx_opw ) 
     463               CALL tab_1d_2d( npti, nptidx(1:npti), ze_i_2d(1:npti,jk,jl), e_i(:,:,jk,jl) ) 
     464            END DO 
     465         END DO 
     466         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_opw_1d(1:npti), sfx_opw ) 
     467         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_opw_1d(1:npti), wfx_opw ) 
     468         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_thd_1d(1:npti), hfx_thd ) 
     469         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_opw_1d(1:npti), hfx_opw ) 
    470470         ! 
    471       ENDIF ! nidx > 0 
     471      ENDIF ! npti > 0 
    472472      ! 
    473473      IF( ln_icediachk )   CALL ice_cons_hsm(1, 'icethd_do', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_ent.F90

    r8562 r8565  
    7878      !  1) Cumulative integral of old enthalpy * thickness and layers interfaces 
    7979      !-------------------------------------------------------------------------- 
    80       zeh_cum0(1:nidx,0) = 0._wp  
    81       zh_cum0 (1:nidx,0) = 0._wp 
     80      zeh_cum0(1:npti,0) = 0._wp  
     81      zh_cum0 (1:npti,0) = 0._wp 
    8282      DO jk0 = 1, nlay_i+2 
    83          DO ji = 1, nidx 
     83         DO ji = 1, npti 
    8484            zeh_cum0(ji,jk0) = zeh_cum0(ji,jk0-1) + eh_i_old(ji,jk0-1) 
    8585            zh_cum0 (ji,jk0) = zh_cum0 (ji,jk0-1) + h_i_old (ji,jk0-1) 
     
    9191      !------------------------------------ 
    9292      ! new layer thickesses 
    93       DO ji = 1, nidx 
     93      DO ji = 1, npti 
    9494         zhnew(ji) = SUM( h_i_old(ji,0:nlay_i+1) ) * r1_nlay_i   
    9595      END DO 
    9696 
    9797      ! new layers interfaces 
    98       zh_cum1(1:nidx,0) = 0._wp 
     98      zh_cum1(1:npti,0) = 0._wp 
    9999      DO jk1 = 1, nlay_i 
    100          DO ji = 1, nidx 
     100         DO ji = 1, npti 
    101101            zh_cum1(ji,jk1) = zh_cum1(ji,jk1-1) + zhnew(ji) 
    102102         END DO 
    103103      END DO 
    104104 
    105       zeh_cum1(1:nidx,0) = 0._wp  
     105      zeh_cum1(1:npti,0) = 0._wp  
    106106      ! new cumulative q*h => linear interpolation 
    107107      DO jk0 = 1, nlay_i+1 
    108108         DO jk1 = 1, nlay_i-1 
    109             DO ji = 1, nidx 
     109            DO ji = 1, npti 
    110110               IF( zh_cum1(ji,jk1) <= zh_cum0(ji,jk0) .AND. zh_cum1(ji,jk1) > zh_cum0(ji,jk0-1) ) THEN 
    111111                  zeh_cum1(ji,jk1) = ( zeh_cum0(ji,jk0-1) * ( zh_cum0(ji,jk0) - zh_cum1(ji,jk1  ) ) +  & 
     
    117117      END DO 
    118118      ! to ensure that total heat content is strictly conserved, set: 
    119       zeh_cum1(1:nidx,nlay_i) = zeh_cum0(1:nidx,nlay_i+2)  
     119      zeh_cum1(1:npti,nlay_i) = zeh_cum0(1:npti,nlay_i+2)  
    120120 
    121121      ! new enthalpies 
    122122      DO jk1 = 1, nlay_i 
    123          DO ji = 1, nidx 
     123         DO ji = 1, npti 
    124124            rswitch      = MAX( 0._wp , SIGN( 1._wp , zhnew(ji) - epsi20 ) )  
    125125            qnew(ji,jk1) = rswitch * ( zeh_cum1(ji,jk1) - zeh_cum1(ji,jk1-1) ) / MAX( zhnew(ji), epsi20 ) 
     
    130130      ! comment: if input h_i_old and eh_i_old are already multiplied by a_i (as in icethd_do),  
    131131      ! then we should not (* a_i) again but not important since this is just to check that remap error is ~0 
    132       DO ji = 1, nidx 
     132      DO ji = 1, npti 
    133133         hfx_err_rem_1d(ji) = hfx_err_rem_1d(ji) + a_i_1d(ji) * r1_rdtice *  & 
    134134            &               ( SUM( qnew(ji,1:nlay_i) ) * zhnew(ji) - SUM( eh_i_old(ji,0:nlay_i+1) ) )  
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_sal.F90

    r8564 r8565  
    7070         z1_time_fl = 1._wp / rn_time_fl * rdt_ice 
    7171         ! 
    72          DO ji = 1, nidx 
     72         DO ji = 1, npti 
    7373 
    7474            !--------------------------------------------------------- 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_zdf.F90

    r8564 r8565  
    144144 
    145145      ! --- diag error on heat diffusion - PART 1 --- ! 
    146       DO ji = 1, nidx 
     146      DO ji = 1, npti 
    147147         zq_ini(ji) = ( SUM( e_i_1d(ji,1:nlay_i) ) * h_i_1d(ji) * r1_nlay_i +  & 
    148148            &           SUM( e_s_1d(ji,1:nlay_s) ) * h_s_1d(ji) * r1_nlay_s )  
     
    152152      ! 1) Initialization 
    153153      !------------------ 
    154       DO ji = 1, nidx 
     154      DO ji = 1, npti 
    155155         isnow(ji)= 1._wp - MAX( 0._wp , SIGN(1._wp, - h_s_1d(ji) ) )  ! is there snow or not 
    156156         ! layer thickness 
     
    159159      END DO 
    160160      ! 
    161       WHERE( zh_i(1:nidx) >= epsi10 )   ;   z1_h_i(1:nidx) = 1._wp / zh_i(1:nidx) 
    162       ELSEWHERE                         ;   z1_h_i(1:nidx) = 0._wp 
     161      WHERE( zh_i(1:npti) >= epsi10 )   ;   z1_h_i(1:npti) = 1._wp / zh_i(1:npti) 
     162      ELSEWHERE                         ;   z1_h_i(1:npti) = 0._wp 
    163163      END WHERE 
    164164 
    165       WHERE( zh_s(1:nidx) >= epsi10 )   ;   z1_h_s(1:nidx) = 1._wp / zh_s(1:nidx) 
    166       ELSEWHERE                         ;   z1_h_s(1:nidx) = 0._wp 
     165      WHERE( zh_s(1:npti) >= epsi10 )   ;   z1_h_s(1:npti) = 1._wp / zh_s(1:npti) 
     166      ELSEWHERE                         ;   z1_h_s(1:npti) = 0._wp 
    167167      END WHERE 
    168168      ! 
    169169      ! temperatures 
    170       ztsub  (1:nidx)   = t_su_1d(1:nidx)  ! temperature at the previous iteration 
    171       ztsold (1:nidx,:) = t_s_1d(1:nidx,:)   ! Old snow temperature 
    172       ztiold (1:nidx,:) = t_i_1d(1:nidx,:)   ! Old ice temperature 
    173       t_su_1d(1:nidx)   = MIN( t_su_1d(1:nidx), rt0 - ztsu_err )   ! necessary 
     170      ztsub  (1:npti)   = t_su_1d(1:npti)  ! temperature at the previous iteration 
     171      ztsold (1:npti,:) = t_s_1d(1:npti,:)   ! Old snow temperature 
     172      ztiold (1:npti,:) = t_i_1d(1:npti,:)   ! Old ice temperature 
     173      t_su_1d(1:npti)   = MIN( t_su_1d(1:npti), rt0 - ztsu_err )   ! necessary 
    174174      ! 
    175175      !------------- 
     
    177177      !------------- 
    178178      z1_hsu = 1._wp / 0.1_wp ! threshold for the computation of i0 
    179       DO ji = 1, nidx 
     179      DO ji = 1, npti 
    180180         ! --- Computation of i0 --- ! 
    181181         ! i0 describes the fraction of solar radiation which does not contribute 
     
    199199 
    200200      ! --- Transmission/absorption of solar radiation in the ice --- ! 
    201       zradtr_s(1:nidx,0) = zftrice(1:nidx) 
     201      zradtr_s(1:npti,0) = zftrice(1:npti) 
    202202      DO jk = 1, nlay_s 
    203          DO ji = 1, nidx 
     203         DO ji = 1, npti 
    204204            !                             ! radiation transmitted below the layer-th snow layer 
    205205            zradtr_s(ji,jk) = zradtr_s(ji,0) * EXP( - zraext_s * zh_s(ji) * REAL(jk) ) 
     
    209209      END DO 
    210210          
    211       zradtr_i(1:nidx,0) = zradtr_s(1:nidx,nlay_s) * isnow(1:nidx) + zftrice(1:nidx) * ( 1._wp - isnow(1:nidx) ) 
     211      zradtr_i(1:npti,0) = zradtr_s(1:npti,nlay_s) * isnow(1:npti) + zftrice(1:npti) * ( 1._wp - isnow(1:npti) ) 
    212212      DO jk = 1, nlay_i  
    213          DO ji = 1, nidx 
     213         DO ji = 1, npti 
    214214            !                             ! radiation transmitted below the layer-th ice layer 
    215215            zradtr_i(ji,jk) = zradtr_i(ji,0) * EXP( - rn_kappa_i * zh_i(ji) * REAL(jk) ) 
     
    219219      END DO 
    220220 
    221       ftr_ice_1d(1:nidx) = zradtr_i(1:nidx,nlay_i)   ! record radiation transmitted below the ice 
     221      ftr_ice_1d(1:npti) = zradtr_i(1:npti,nlay_i)   ! record radiation transmitted below the ice 
    222222      ! 
    223223      iconv    =  0          ! number of iterations 
     
    228228         iconv = iconv + 1 
    229229         ! 
    230          ztib(1:nidx,:) = t_i_1d(1:nidx,:) 
    231          ztsb(1:nidx,:) = t_s_1d(1:nidx,:) 
     230         ztib(1:npti,:) = t_i_1d(1:npti,:) 
     231         ztsb(1:npti,:) = t_s_1d(1:npti,:) 
    232232         ! 
    233233         !-------------------------------- 
     
    236236         IF( ln_cndi_U64 ) THEN         !-- Untersteiner (1964) formula: k = k0 + beta.S/T 
    237237            ! 
    238             DO ji = 1, nidx 
     238            DO ji = 1, npti 
    239239               ztcond_i(ji,0)      = rcdic + zbeta * sz_i_1d(ji,1)      / MIN( -epsi10, t_i_1d(ji,1) - rt0 ) 
    240240               ztcond_i(ji,nlay_i) = rcdic + zbeta * sz_i_1d(ji,nlay_i) / MIN( -epsi10, t_bo_1d(ji)  - rt0 ) 
    241241            END DO 
    242242            DO jk = 1, nlay_i-1 
    243                DO ji = 1, nidx 
     243               DO ji = 1, npti 
    244244                  ztcond_i(ji,jk) = rcdic + zbeta * 0.5_wp * ( sz_i_1d(ji,jk) + sz_i_1d(ji,jk+1) ) /  & 
    245245                     &                      MIN( -epsi10, 0.5_wp * (t_i_1d(ji,jk) + t_i_1d(ji,jk+1)) - rt0 ) 
     
    249249         ELSEIF( ln_cndi_P07 ) THEN     !-- Pringle et al formula: k = k0 + beta1.S/T - beta2.T 
    250250            ! 
    251             DO ji = 1, nidx 
     251            DO ji = 1, npti 
    252252               ztcond_i(ji,0)      = rcdic + 0.09_wp * sz_i_1d(ji,1) / MIN( -epsi10, t_i_1d(ji,1) - rt0 )  & 
    253253                  &                        - 0.011_wp * ( t_i_1d(ji,1) - rt0 ) 
     
    256256            END DO 
    257257            DO jk = 1, nlay_i-1 
    258                DO ji = 1, nidx 
     258               DO ji = 1, npti 
    259259                  ztcond_i(ji,jk) = rcdic + 0.09_wp * 0.5_wp * ( sz_i_1d(ji,jk) + sz_i_1d(ji,jk+1) ) /          & 
    260260                     &                     MIN( -epsi10, 0.5_wp * (t_i_1d(ji,jk) + t_i_1d(ji,jk+1)) - rt0 )   & 
     
    264264            ! 
    265265         ENDIF 
    266          ztcond_i(1:nidx,:) = MAX( zkimin, ztcond_i(1:nidx,:) )         
     266         ztcond_i(1:npti,:) = MAX( zkimin, ztcond_i(1:npti,:) )         
    267267         ! 
    268268         !--- G(he) : enhancement of thermal conductivity in mono-category case 
     
    270270         ! Used in mono-category case only to simulate an ITD implicitly 
    271271         ! Fichefet and Morales Maqueda, JGR 1997 
    272          zghe(1:nidx) = 1._wp 
     272         zghe(1:npti) = 1._wp 
    273273         ! 
    274274         SELECT CASE ( nn_monocat ) 
     
    277277 
    278278            zepsilon = 0.1_wp 
    279             DO ji = 1, nidx 
     279            DO ji = 1, npti 
    280280               zcnd_i = SUM( ztcond_i(ji,:) ) / REAL( nlay_i+1, wp )                             ! Mean sea ice thermal conductivity 
    281281               zhe = ( rn_cnd_s * h_i_1d(ji) + zcnd_i * h_s_1d(ji) ) / ( rn_cnd_s + zcnd_i )   ! Effective thickness he (zhe) 
     
    292292         !--- Snow 
    293293         DO jk = 0, nlay_s-1 
    294             DO ji = 1, nidx 
     294            DO ji = 1, npti 
    295295               zkappa_s(ji,jk) = zghe(ji) * rn_cnd_s * z1_h_s(ji) 
    296296            END DO 
    297297         END DO 
    298          DO ji = 1, nidx   ! Snow-ice interface 
     298         DO ji = 1, npti   ! Snow-ice interface 
    299299            zfac = 0.5_wp * ( ztcond_i(ji,0) * zh_s(ji) + rn_cnd_s * zh_i(ji) ) 
    300300            IF( zfac > epsi10 ) THEN 
     
    307307         !--- Ice 
    308308         DO jk = 0, nlay_i 
    309             DO ji = 1, nidx 
     309            DO ji = 1, npti 
    310310               zkappa_i(ji,jk) = zghe(ji) * ztcond_i(ji,jk) * z1_h_i(ji) 
    311311            END DO 
    312312         END DO 
    313          DO ji = 1, nidx   ! Snow-ice interface 
     313         DO ji = 1, npti   ! Snow-ice interface 
    314314            zkappa_i(ji,0) = zkappa_s(ji,nlay_s) * isnow(ji) + zkappa_i(ji,0) * ( 1._wp - isnow(ji) ) 
    315315         END DO 
     
    319319         !-------------------------------------- 
    320320         DO jk = 1, nlay_i 
    321             DO ji = 1, nidx 
     321            DO ji = 1, npti 
    322322               zcpi = cpic + zgamma * sz_i_1d(ji,jk) / MAX( ( t_i_1d(ji,jk) - rt0 ) * ( ztiold(ji,jk) - rt0 ), epsi10 ) 
    323323               zeta_i(ji,jk) = rdt_ice * r1_rhoic * z1_h_i(ji) / MAX( epsi10, zcpi )  
     
    326326 
    327327         DO jk = 1, nlay_s 
    328             DO ji = 1, nidx 
     328            DO ji = 1, npti 
    329329               zeta_s(ji,jk) = rdt_ice * r1_rhosn * r1_cpic * z1_h_s(ji) 
    330330            END DO 
     
    335335         !---------------------------- 
    336336         IF ( ln_dqns_i ) THEN  
    337             DO ji = 1, nidx 
     337            DO ji = 1, npti 
    338338               ! update of the non solar flux according to the update in T_su 
    339339               qns_ice_1d(ji) = qns_ice_1d(ji) + dqns_ice_1d(ji) * ( t_su_1d(ji) - ztsub(ji) ) 
     
    341341         ENDIF 
    342342 
    343          DO ji = 1, nidx 
     343         DO ji = 1, npti 
    344344            zfnet(ji) = zfsw(ji) + qns_ice_1d(ji) ! incoming = net absorbed solar radiation + non solar total flux (LWup, LWdw, SH, LH) 
    345345         END DO 
     
    354354 
    355355         !!ice interior terms (top equation has the same form as the others) 
    356          ztrid   (1:nidx,:,:) = 0._wp 
    357          zindterm(1:nidx,:)   = 0._wp 
    358          zindtbis(1:nidx,:)   = 0._wp 
    359          zdiagbis(1:nidx,:)   = 0._wp 
     356         ztrid   (1:npti,:,:) = 0._wp 
     357         zindterm(1:npti,:)   = 0._wp 
     358         zindtbis(1:npti,:)   = 0._wp 
     359         zdiagbis(1:npti,:)   = 0._wp 
    360360 
    361361         DO jm = nlay_s + 2, nlay_s + nlay_i  
    362             DO ji = 1, nidx 
     362            DO ji = 1, npti 
    363363               jk = jm - nlay_s - 1 
    364364               ztrid(ji,jm,1)  = - zeta_i(ji,jk) * zkappa_i(ji,jk-1) 
     
    370370 
    371371         jm =  nlay_s + nlay_i + 1 
    372          DO ji = 1, nidx 
     372         DO ji = 1, npti 
    373373            !!ice bottom term 
    374374            ztrid(ji,jm,1)  = - zeta_i(ji,nlay_i)*zkappa_i(ji,nlay_i-1)    
     
    380380 
    381381 
    382          DO ji = 1, nidx 
     382         DO ji = 1, npti 
    383383            !                               !---------------------! 
    384384            IF ( h_s_1d(ji) > 0.0 ) THEN   !  snow-covered cells ! 
     
    497497         jm_maxt = 0 
    498498         jm_mint = nlay_i+5 
    499          DO ji = 1, nidx 
     499         DO ji = 1, npti 
    500500            jm_mint = MIN(jm_min(ji),jm_mint) 
    501501            jm_maxt = MAX(jm_max(ji),jm_maxt) 
     
    503503 
    504504         DO jk = jm_mint+1, jm_maxt 
    505             DO ji = 1, nidx 
     505            DO ji = 1, npti 
    506506               jm = min(max(jm_min(ji)+1,jk),jm_max(ji)) 
    507507               zdiagbis(ji,jm) = ztrid(ji,jm,2)  - ztrid(ji,jm,1) * ztrid(ji,jm-1,3)  / zdiagbis(ji,jm-1) 
     
    510510         END DO 
    511511 
    512          DO ji = 1, nidx 
     512         DO ji = 1, npti 
    513513            ! ice temperatures 
    514514            t_i_1d(ji,nlay_i) = zindtbis(ji,jm_max(ji)) / zdiagbis(ji,jm_max(ji)) 
     
    516516 
    517517         DO jm = nlay_i + nlay_s, nlay_s + 2, -1 
    518             DO ji = 1, nidx 
     518            DO ji = 1, npti 
    519519               jk    =  jm - nlay_s - 1 
    520520               t_i_1d(ji,jk) = ( zindtbis(ji,jm) - ztrid(ji,jm,3) * t_i_1d(ji,jk+1) ) / zdiagbis(ji,jm) 
     
    522522         END DO 
    523523 
    524          DO ji = 1, nidx 
     524         DO ji = 1, npti 
    525525            ! snow temperatures       
    526526            IF( h_s_1d(ji) > 0._wp ) THEN 
     
    542542         ! zdti_max is a measure of error, it has to be under zdti_bnd 
    543543         zdti_max = 0._wp 
    544          DO ji = 1, nidx 
     544         DO ji = 1, npti 
    545545            t_su_1d(ji) = MAX( MIN( t_su_1d(ji) , rt0 ) , rt0 - 100._wp ) 
    546546            zdti_max = MAX( zdti_max, ABS( t_su_1d(ji) - ztsub(ji) ) )      
     
    548548 
    549549         DO jk  =  1, nlay_s 
    550             DO ji = 1, nidx 
     550            DO ji = 1, npti 
    551551               t_s_1d(ji,jk) = MAX( MIN( t_s_1d(ji,jk), rt0 ), rt0 - 100._wp ) 
    552552               zdti_max = MAX( zdti_max, ABS( t_s_1d(ji,jk) - ztsb(ji,jk) ) ) 
     
    555555 
    556556         DO jk  =  1, nlay_i 
    557             DO ji = 1, nidx 
     557            DO ji = 1, npti 
    558558               ztmelt_i      = -tmut * sz_i_1d(ji,jk) + rt0  
    559559               t_i_1d(ji,jk) =  MAX( MIN( t_i_1d(ji,jk), ztmelt_i ), rt0 - 100._wp ) 
     
    576576      ! 10) Fluxes at the interfaces 
    577577      !----------------------------- 
    578       DO ji = 1, nidx 
     578      DO ji = 1, npti 
    579579         !                                ! surface ice conduction flux 
    580580         fc_su(ji)   =  -           isnow(ji)   * zkappa_s(ji,0) * zg1s * (t_s_1d(ji,1) - t_su_1d(ji))   & 
     
    589589      ! --- diagnose the change in non-solar flux due to surface temperature change --- ! 
    590590      IF ( ln_dqns_i ) THEN 
    591          DO ji = 1, nidx 
     591         DO ji = 1, npti 
    592592            hfx_err_dif_1d(ji) = hfx_err_dif_1d(ji) - ( qns_ice_1d(ji)  - zqns_ice_b(ji) ) * a_i_1d(ji)  
    593593         END DO 
     
    597597      !     hfx_dif = Heat flux used to warm/cool ice in W.m-2 
    598598      !     zhfx_err = correction on the diagnosed heat flux due to non-convergence of the algorithm used to solve heat equation 
    599       DO ji = 1, nidx 
     599      DO ji = 1, npti 
    600600         zdq = - zq_ini(ji) + ( SUM( e_i_1d(ji,1:nlay_i) ) * h_i_1d(ji) * r1_nlay_i +  & 
    601601            &                   SUM( e_s_1d(ji,1:nlay_s) ) * h_s_1d(ji) * r1_nlay_s ) 
     
    614614 
    615615      ! --- Diagnostics SIMIP --- ! 
    616       DO ji = 1, nidx 
     616      DO ji = 1, npti 
    617617         !--- Conduction fluxes (positive downwards) 
    618618         diag_fc_bo_1d(ji) = diag_fc_bo_1d(ji) + fc_bo_i(ji) * a_i_1d(ji) / at_i_1d(ji) 
     
    645645      ! 
    646646      DO jk = 1, nlay_i             ! Sea ice energy of melting 
    647          DO ji = 1, nidx 
     647         DO ji = 1, npti 
    648648            ztmelts      = - tmut  * sz_i_1d(ji,jk) 
    649649            t_i_1d(ji,jk) = MIN( t_i_1d(ji,jk), ztmelts + rt0 ) ! Force t_i_1d to be lower than melting point 
     
    655655      END DO 
    656656      DO jk = 1, nlay_s             ! Snow energy of melting 
    657          DO ji = 1, nidx 
     657         DO ji = 1, npti 
    658658            e_s_1d(ji,jk) = rhosn * ( cpic * ( rt0 - t_s_1d(ji,jk) ) + lfus ) 
    659659         END DO 
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90

    r8564 r8565  
    378378      CASE( 1 )       !  constant salinity in time and space  ! 
    379379         !            !---------------------------------------! 
    380          sz_i_1d(1:nidx,:) = rn_icesal 
     380         sz_i_1d(1:npti,:) = rn_icesal 
    381381         ! 
    382382         !            !---------------------------------------------! 
     
    387387         ! 
    388388         !                                      ! Slope of the linear profile  
    389          WHERE( h_i_1d(1:nidx) > epsi20 )   ;   z_slope_s(1:nidx) = 2._wp * s_i_1d(1:nidx) / h_i_1d(1:nidx) 
    390          ELSEWHERE                           ;   z_slope_s(1:nidx) = 0._wp 
     389         WHERE( h_i_1d(1:npti) > epsi20 )   ;   z_slope_s(1:npti) = 2._wp * s_i_1d(1:npti) / h_i_1d(1:npti) 
     390         ELSEWHERE                           ;   z_slope_s(1:npti) = 0._wp 
    391391         END WHERE 
    392392          
    393393         z1_dS = 1._wp / ( zsi1 - zsi0 ) 
    394          DO ji = 1, nidx 
     394         DO ji = 1, npti 
    395395            zalpha(ji) = MAX(  0._wp , MIN(  ( zsi1 - s_i_1d(ji) ) * z1_dS , 1._wp  )  ) 
    396396            !                             ! force a constant profile when SSS too low (Baltic Sea) 
     
    400400         ! Computation of the profile 
    401401         DO jk = 1, nlay_i 
    402             DO ji = 1, nidx 
     402            DO ji = 1, npti 
    403403               !                          ! linear profile with 0 surface value 
    404404               zs0 = z_slope_s(ji) * ( REAL(jk,wp) - 0.5_wp ) * h_i_1d(ji) * r1_nlay_i 
     
    414414         !            !-------------------------------------------!                                   (mean = 2.30) 
    415415         ! 
    416          s_i_1d(1:nidx) = 2.30_wp 
     416         s_i_1d(1:npti) = 2.30_wp 
    417417         ! 
    418418!!gm cf remark in ice_var_salprof routine, CASE( 3 ) 
     
    420420            zargtemp  = ( REAL(jk,wp) - 0.5_wp ) * r1_nlay_i 
    421421            zsal =  1.6_wp * ( 1._wp - COS( rpi * zargtemp**( 0.407_wp / ( 0.573_wp + zargtemp ) ) ) ) 
    422             DO ji = 1, nidx 
     422            DO ji = 1, npti 
    423423               sz_i_1d(ji,jk) = zsal 
    424424            END DO 
Note: See TracChangeset for help on using the changeset viewer.