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 921 for trunk/NEMO/LIM_SRC_3/limdia.F90 – NEMO

Ignore:
Timestamp:
2008-05-13T10:28:52+02:00 (16 years ago)
Author:
rblod
Message:

Correct indentation and print for debug in LIM3, see ticket #134, step I

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMO/LIM_SRC_3/limdia.F90

    r895 r921  
    9595      !!------------------------------------------------------------------- 
    9696      !! * Local variables 
    97        INTEGER  ::   jv,ji,jj,jl ! dummy loop indices 
    98        REAL(wp), DIMENSION(jpinfmx) ::  &  
    99           vinfor           ! temporary working space  
    100        REAL(wp) ::    & 
    101           zshift_date   , & ! date from the minimum ice extent 
    102           zday, zday_min, & ! current day, day of minimum extent 
    103           zafy, zamy,     & ! temporary area of fy and my ice 
    104           zindb 
    105        !!------------------------------------------------------------------- 
    106  
    107        ! 0) date from the minimum of ice extent 
    108        !--------------------------------------- 
    109        zday_min = 273.0        ! zday_min = date of minimum extent, here September 30th 
    110        zday = FLOAT(numit-nit000) * rdt_ice / ( 86400.0 * FLOAT(nn_fsbc) ) 
    111        IF (zday.GT.zday_min) THEN  
    112           zshift_date  =  zday - zday_min 
    113        ELSE 
    114           zshift_date  =  zday - (365.0 - zday_min) 
    115        ENDIF 
    116  
    117        IF( numit == nstart )   CALL lim_dia_init   ! initialisation of ice_evolu file       
    118  
    119        ! temporal diagnostics  
    120        vinfor(1) = REAL(numit) 
    121        vinfor(2) = nyear 
    122   
    123        ! put everything to zero 
    124        DO jv = nbvt + 1, nvinfo 
    125           vinfor(jv) = 0.0 
    126        END DO 
    127  
    128        !!------------------------------------------------------------------- 
    129        !! 1) Northern hemisphere 
    130        !!------------------------------------------------------------------- 
    131        !! 1.1) Diagnostics independent on age 
    132        !!------------------------------------ 
    133        DO jj = njeq, jpjm1 
    134           DO ji = fs_2, fs_jpim1   ! vector opt. 
    135              IF( tms(ji,jj) == 1 ) THEN 
    136                 vinfor(3)  = vinfor(3)  + at_i(ji,jj)*aire(ji,jj) / 1.0e12 !ice area 
    137                 IF (at_i(ji,jj).GT.0.15) vinfor(5) = vinfor(5) + aire(ji,jj) / 1.0e12 !ice extent 
    138                 vinfor(7)  = vinfor(7)  + vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !ice volume 
    139                 vinfor(9)  = vinfor(9)  + vt_s(ji,jj)*aire(ji,jj) / 1.0e12 !snow volume 
    140                 vinfor(15) = vinfor(15) + ot_i(ji,jj) *vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !mean age 
    141                 vinfor(29) = vinfor(29) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !mean salinity 
    142                 ! the computation of this diagnostic is not reliable 
    143                 vinfor(31) = vinfor(31) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + &  
    144                                                         v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12  
    145                 vinfor(53) = vinfor(53) + emps(ji,jj)*aire(ji,jj) / 1.0e12 !salt flux 
    146                 vinfor(55) = vinfor(55) + fsbri(ji,jj)*aire(ji,jj) / 1.0e12 !brine drainage flux 
    147                 vinfor(57) = vinfor(57) + fseqv(ji,jj)*aire(ji,jj) / 1.0e12 !equivalent salt flux 
    148                 vinfor(59) = vinfor(59) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SST 
    149                 vinfor(61) = vinfor(61) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SSS 
    150                 vinfor(65) = vinfor(65) + et_s(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12  ! snow temperature 
    151                 vinfor(67) = vinfor(67) + et_i(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12       ! ice heat content 
    152                 vinfor(69) = vinfor(69) + v_i(ji,jj,1)*aire(ji,jj) / 1.0e12 !ice volume 
    153                 vinfor(71) = vinfor(71) + v_i(ji,jj,2)*aire(ji,jj) / 1.0e12 !ice volume 
    154                 vinfor(73) = vinfor(73) + v_i(ji,jj,3)*aire(ji,jj) / 1.0e12 !ice volume 
    155                 vinfor(75) = vinfor(75) + v_i(ji,jj,4)*aire(ji,jj) / 1.0e12 !ice volume 
    156                 vinfor(77) = vinfor(77) + v_i(ji,jj,5)*aire(ji,jj) / 1.0e12 !ice volume 
    157                 vinfor(79) = 0.0 
    158                 vinfor(81) = vinfor(81) + emp(ji,jj)*aire(ji,jj) / 1.0e12 ! mass flux 
    159              ENDIF 
    160           END DO 
    161        END DO 
    162  
    163        DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 
    164           DO jj = njeq, jpjm1 
    165              DO ji = fs_2, fs_jpim1   ! vector opt. 
    166                 IF( tms(ji,jj) == 1 ) THEN 
    167                    vinfor(11) = vinfor(11) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !undef def ice volume 
    168                 ENDIF 
    169              END DO 
    170           END DO 
    171        END DO 
    172  
    173        vinfor(13) = 0.0 
    174  
    175        vinfor(15) = vinfor(15) / MAX(vinfor(7),epsi06) ! these have to be divided by total ice volume to have the 
    176        vinfor(29) = vinfor(29) / MAX(vinfor(7),epsi06) ! right value 
    177        vinfor(31) = SQRT( vinfor(31) / MAX( vinfor(7) , epsi06 ) ) 
    178        vinfor(67) = vinfor(67) / MAX(vinfor(7),epsi06) 
    179  
    180        vinfor(53) = vinfor(53) / MAX(vinfor(5),epsi06) ! these have to be divided by total ice extent to have the 
    181        vinfor(55) = vinfor(55) / MAX(vinfor(5),epsi06) ! right value  
    182        vinfor(57) = vinfor(57) / MAX(vinfor(5),epsi06) !  
    183        vinfor(79) = vinfor(79) / MAX(vinfor(5),epsi06) ! 
    184  
    185        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(3))) ! 
    186        vinfor(59) = zindb*vinfor(59) / MAX(vinfor(3),epsi06) ! divide by ice area 
    187        vinfor(61) = zindb*vinfor(61) / MAX(vinfor(3),epsi06) ! 
    188  
    189        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(9))) ! 
    190        vinfor(65) = zindb*vinfor(65) / MAX(vinfor(9),epsi06) ! divide it by snow volume 
    191  
    192  
    193        DO jl = 1, jpl 
    194           DO jj = njeq, jpjm1 
    195              DO ji = fs_2, fs_jpim1   ! vector opt. 
    196                 IF( tms(ji,jj) == 1 ) THEN 
    197                    vinfor(33) = vinfor(33) + d_v_i_trp(ji,jj,jl)*aire(ji,jj) / 1.0e12 !ice volume 
    198                    vinfor(35) = vinfor(35) + d_v_i_thd(ji,jj,jl)*aire(ji,jj) / 1.0e12 !ice volume 
    199                 ENDIF 
    200              END DO 
    201           END DO 
    202        END DO 
    203  
    204        DO jj = njeq, jpjm1 
    205           DO ji = fs_2, fs_jpim1   ! vector opt. 
    206                 IF( tms(ji,jj) == 1 ) THEN 
    207                    vinfor(37) = vinfor(37) + diag_sni_gr(ji,jj)*aire(ji,jj) / 1.0e12 !th growth rates 
    208                    vinfor(39) = vinfor(39) + diag_lat_gr(ji,jj)*aire(ji,jj) / 1.0e12  
    209                    vinfor(41) = vinfor(41) + diag_bot_gr(ji,jj)*aire(ji,jj) / 1.0e12 
    210                    vinfor(43) = vinfor(43) + diag_dyn_gr(ji,jj)*aire(ji,jj) / 1.0e12  
    211                    vinfor(45) = vinfor(45) + dv_dt_thd(ji,jj,5)*aire(ji,jj) / 1.0e12 
    212                    vinfor(47) = vinfor(47) + v_newice(ji,jj) *aire(ji,jj) / 1.0e12 / rdt_ice ! volume acc in OW 
    213                 ENDIF 
    214           END DO 
    215        END DO 
    216  
    217        DO jl = 1, jpl 
    218           DO jj = njeq, jpjm1 
    219              DO ji = fs_2, fs_jpim1   ! vector opt. 
    220                 IF( tms(ji,jj) == 1 ) THEN 
    221                    vinfor(63) = vinfor(63) + t_su(ji,jj,jl)*a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 
    222                 ENDIF 
    223              END DO 
    224           END DO 
    225        END DO 
    226        vinfor(63) = vinfor(63) / MAX(vinfor(3),epsi06) ! these have to be divided by total ice area 
    227  
    228        !! 1.2) Diagnostics dependent on age 
    229        !!------------------------------------ 
    230        DO jj = njeq, jpjm1 
    231           DO ji = fs_2, fs_jpim1   ! vector opt. 
    232              IF( tms(ji,jj) == 1 ) THEN 
    233                 zafy = 0.0 
    234                 zamy = 0.0 
    235                 DO jl = 1, jpl 
    236                    IF ((o_i(ji,jj,jl) - zshift_date).LT.0.0) THEN 
    237                       vinfor(17) = vinfor(17) + a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! FY ice area 
    238                       vinfor(25) = vinfor(25) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! FY ice volume 
    239                       vinfor(49) = vinfor(49) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !FY ice salinity 
    240                       zafy = zafy + a_i(ji,jj,jl) 
    241                    ENDIF 
    242                    IF ((o_i(ji,jj,jl) - zshift_date).GT.0.0) THEN 
    243                       vinfor(19) = vinfor(19) + a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12    ! MY ice area 
    244                       vinfor(27) = vinfor(27) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! MY ice volume 
    245                       vinfor(51) = vinfor(51) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !MY ice salinity 
    246                       zamy = zamy + a_i(ji,jj,jl) 
    247                    ENDIF 
    248                 END DO 
    249                 IF ((at_i(ji,jj).GT.0.15).AND.(zafy.GT.zamy)) THEN 
    250                    vinfor(21) = vinfor(21) + aire(ji,jj) / 1.0e12 ! Seasonal ice extent 
    251                 ENDIF 
    252                 IF ((at_i(ji,jj).GT.0.15).AND.(zafy.LE.zamy)) THEN 
    253                    vinfor(23) = vinfor(23) + aire(ji,jj) / 1.0e12 ! Perennial ice extent 
    254                 ENDIF 
    255              ENDIF 
    256           END DO 
    257        END DO 
    258        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(25))) !=0 if no multiyear ice 1 if yes 
    259        vinfor(49) = zindb*vinfor(49) / MAX(vinfor(25),epsi06) 
    260        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(27))) !=0 if no multiyear ice 1 if yes 
    261        vinfor(51) = zindb*vinfor(51) / MAX(vinfor(27),epsi06) 
    262         
    263        !! Fram Strait Export 
    264        !! 83 = area export 
    265        !! 84 = volume export 
    266        !! Fram strait in ORCA2 = 5 points 
    267        !! export = -v_ice*e1t*ddtb*at_i or -v_ice*e1t*ddtb*at_i*h_i 
    268        jj = 136 ! C grid 
    269        vinfor(83) = 0.0 
    270        vinfor(84) = 0.0 
    271        DO ji = 134, 138 
    272           vinfor(83) = vinfor(83) - v_ice(ji,jj) * &  
    273                                       e1t(ji,jj)*at_i(ji,jj)*rdt_ice / 1.0e12 
    274           vinfor(84) = vinfor(84) - v_ice(ji,jj) * &  
    275                                       e1t(ji,jj)*vt_i(ji,jj)*rdt_ice / 1.0e12 
    276        END DO 
    277  
    278        !!------------------------------------------------------------------- 
    279        !! 2) Southern hemisphere 
    280        !!------------------------------------------------------------------- 
    281        !! 2.1) Diagnostics independent on age 
    282        !!------------------------------------ 
    283        DO jj = 2, njeqm1 
    284           DO ji = fs_2, fs_jpim1   ! vector opt. 
    285              IF( tms(ji,jj) == 1 ) THEN 
    286                 vinfor(4)  = vinfor(4)  + at_i(ji,jj)*aire(ji,jj) / 1.0e12 !ice area 
    287                 IF (at_i(ji,jj).GT.0.15) vinfor(6) = vinfor(6) + aire(ji,jj) / 1.0e12 !ice extent 
    288                 vinfor(8)  = vinfor(8)  + vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !ice volume 
    289                 vinfor(10) = vinfor(10) + vt_s(ji,jj)*aire(ji,jj) / 1.0e12 !snow volume 
    290                 vinfor(16) = vinfor(16) + ot_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !mean age 
    291                 vinfor(30) = vinfor(30) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !mean salinity 
    292                 ! this diagnostic is not well computed (weighted by vol instead 
    293                 ! of area) 
    294                 vinfor(32) = vinfor(32) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + &  
    295                                                         v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12 !ice vel 
    296                 vinfor(54) = vinfor(54) + at_i(ji,jj)*emps(ji,jj)*aire(ji,jj) / 1.0e12 ! Total salt flux 
    297                 vinfor(56) = vinfor(56) + at_i(ji,jj)*fsbri(ji,jj)*aire(ji,jj) / 1.0e12 ! Brine drainage salt flux 
    298                 vinfor(58) = vinfor(58) + at_i(ji,jj)*fseqv(ji,jj)*aire(ji,jj) / 1.0e12 ! Equivalent salt flux 
    299                 vinfor(60) = vinfor(60) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SST 
    300                 vinfor(62) = vinfor(62) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SSS 
    301                 vinfor(66) = vinfor(66) + et_s(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12 ! snow temperature 
    302                 vinfor(68) = vinfor(68) + et_i(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12 ! ice enthalpy 
    303                 vinfor(70) = vinfor(70) + v_i(ji,jj,1)*aire(ji,jj) / 1.0e12 !ice volume 
    304                 vinfor(72) = vinfor(72) + v_i(ji,jj,2)*aire(ji,jj) / 1.0e12 !ice volume 
    305                 vinfor(74) = vinfor(74) + v_i(ji,jj,3)*aire(ji,jj) / 1.0e12 !ice volume 
    306                 vinfor(76) = vinfor(76) + v_i(ji,jj,4)*aire(ji,jj) / 1.0e12 !ice volume 
    307                 vinfor(78) = vinfor(78) + v_i(ji,jj,5)*aire(ji,jj) / 1.0e12 !ice volume 
    308                 vinfor(80) = 0.0 
    309                 vinfor(82) = vinfor(82) + emp(ji,jj)*aire(ji,jj) / 1.0e12 ! mass flux 
    310              ENDIF 
    311           END DO 
    312        END DO 
    313  
    314        DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 
    315           DO jj = 2, njeqm1 
    316              DO ji = fs_2, fs_jpim1   ! vector opt. 
    317                 vinfor(12) = vinfor(12) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !undef def ice volume 
    318              END DO 
    319           END DO 
    320        END DO 
    321  
    322        vinfor(14) = 0.0 
    323  
    324        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(8)))  
    325        vinfor(16) = zindb * vinfor(16) / MAX(vinfor(8),epsi06) ! these have to be divided by ice vol 
    326        vinfor(30) = zindb * vinfor(30) / MAX(vinfor(8),epsi06) !  
    327        vinfor(32) = zindb * SQRT( vinfor(32) / MAX( vinfor(8) , epsi06 ) ) 
    328        vinfor(68) = zindb * vinfor(68) / MAX(vinfor(8),epsi06) !  
    329  
    330        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(6)))  
    331        vinfor(54) = zindb * vinfor(54) / MAX(vinfor(6),epsi06) ! these have to be divided by ice extt 
    332        vinfor(56) = zindb * vinfor(56) / MAX(vinfor(6),epsi06) !  
    333        vinfor(58) = zindb * vinfor(58) / MAX(vinfor(6),epsi06) !  
    334        vinfor(80) = zindb * vinfor(80) / MAX(vinfor(6),epsi06) ! 
    335 !      vinfor(84) = vinfor(84) / vinfor(6) ! 
    336   
    337        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(4))) ! 
    338        vinfor(60) = zindb*vinfor(60) / ( MAX(vinfor(4), epsi06) ) ! divide by ice area 
    339        vinfor(62) = zindb*vinfor(62) / ( MAX(vinfor(4), epsi06) ) ! 
    340  
    341        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(10))) ! 
    342        vinfor(66) = zindb*vinfor(66) / MAX(vinfor(10),epsi06) ! divide it by snow volume 
    343  
    344        DO jl = 1, jpl 
    345           DO jj = 2, njeqm1 
    346              DO ji = fs_2, fs_jpim1   ! vector opt. 
    347                 IF( tms(ji,jj) == 1 ) THEN 
    348                    vinfor(34) = vinfor(34) + d_v_i_trp(ji,jj,jl)*aire(ji,jj) / 1.0e12 !ice volume 
    349                    vinfor(36) = vinfor(36) + d_v_i_thd(ji,jj,jl)*aire(ji,jj) / 1.0e12 !ice volume 
    350                 ENDIF 
    351              END DO 
    352           END DO 
    353        END DO 
    354  
    355        DO jj = 2, njeqm1 
    356           DO ji = fs_2, fs_jpim1   ! vector opt. 
    357                 IF( tms(ji,jj) == 1 ) THEN 
    358                    vinfor(38) = vinfor(38) + diag_sni_gr(ji,jj)*aire(ji,jj) / 1.0e12 !th growth rates 
    359                    vinfor(40) = vinfor(40) + diag_lat_gr(ji,jj)*aire(ji,jj) / 1.0e12  
    360                    vinfor(42) = vinfor(42) + diag_bot_gr(ji,jj)*aire(ji,jj) / 1.0e12 
    361                    vinfor(44) = vinfor(44) + diag_dyn_gr(ji,jj)*aire(ji,jj) / 1.0e12  
    362                    vinfor(46) = vinfor(46) + dv_dt_thd(ji,jj,5)*aire(ji,jj) / 1.0e12 
    363                    vinfor(48) = vinfor(48) + v_newice(ji,jj) *aire(ji,jj) / 1.0e12 / rdt_ice ! volume acc in OW 
    364                 ENDIF 
    365           END DO 
    366        END DO 
    367  
    368  
    369        DO jl = 1, jpl 
    370           DO jj = 2, njeqm1 
    371              DO ji = fs_2, fs_jpim1   ! vector opt. 
    372                 IF( tms(ji,jj) == 1 ) THEN 
    373                    vinfor(64) = vinfor(64) + t_su(ji,jj,jl)*a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 
    374                 ENDIF 
    375              END DO 
    376           END DO 
    377        END DO 
    378        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(4))) ! 
    379        vinfor(64) = zindb * vinfor(64) / MAX(vinfor(4),epsi06) ! divide by ice extt 
    380        !! 2.2) Diagnostics dependent on age 
    381        !!------------------------------------ 
    382        DO jj = 2, njeqm1 
    383           DO ji = fs_2, fs_jpim1   ! vector opt. 
    384              IF( tms(ji,jj) == 1 ) THEN 
    385                 zafy = 0.0 
    386                 zamy = 0.0 
    387                 DO jl = 1, jpl 
    388                    IF ((o_i(ji,jj,jl) - zshift_date).LT.0.0) THEN 
    389                       vinfor(18) = vinfor(18) + a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! FY ice area 
    390                       vinfor(26) = vinfor(26) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! FY ice volume 
    391                       zafy = zafy + a_i(ji,jj,jl) 
    392                       vinfor(50) = vinfor(50) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !FY ice salinity 
    393                    ENDIF 
    394                    IF ((o_i(ji,jj,jl) - zshift_date).GT.0.0) THEN 
    395                       vinfor(20) = vinfor(20) + a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12    ! MY ice area 
    396                       vinfor(28) = vinfor(28) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 
    397                       vinfor(52) = vinfor(52) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !FY ice salinity 
    398                       zamy = zamy + a_i(ji,jj,jl) 
    399                    ENDIF 
    400                 END DO ! jl 
    401                 IF ((at_i(ji,jj).GT.0.15).AND.(zafy.GT.zamy)) THEN 
    402                    vinfor(22) = vinfor(22) + aire(ji,jj) / 1.0e12 ! Seasonal ice extent 
    403                 ENDIF 
    404                 IF ((at_i(ji,jj).GT.0.15).AND.(zafy.LE.zamy)) THEN 
    405                    vinfor(24) = vinfor(24) + aire(ji,jj) / 1.0e12 ! Perennial ice extent 
    406                 ENDIF 
    407              ENDIF ! tms 
    408           END DO ! jj 
    409        END DO ! ji 
    410        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(26))) !=0 if no multiyear ice 1 if yes 
    411        vinfor(50) = zindb*vinfor(50) / MAX(vinfor(26),epsi06) 
    412        zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(28))) !=0 if no multiyear ice 1 if yes 
    413        vinfor(52) = zindb*vinfor(52) / MAX(vinfor(28),epsi06) 
    414  
    415        !  Accumulation before averaging  
    416        DO jv = 1, nvinfo 
    417           vinfom(jv) = vinfom(jv) + vinfor(jv) 
    418        END DO 
    419        naveg = naveg + 1   
    420      
    421        ! oututs on file ice_evolu     
    422 !MV      IF( MOD( numit , ninfo ) == 0 ) THEN 
    423           WRITE(numevo_ice,fmtw) ( titvar(jv), vinfom(jv)/naveg, jv = 1, nvinfo ) 
    424           naveg = 0 
    425           DO jv = 1, nvinfo 
    426              vinfom(jv)=0.0 
    427           END DO 
    428 !MV      ENDIF 
    429    
    430     END SUBROUTINE lim_dia 
    431   
    432     SUBROUTINE lim_dia_init 
    433        !!------------------------------------------------------------------- 
    434        !!                  ***  ROUTINE lim_dia_init  *** 
    435        !!              
    436        !! ** Purpose : Preparation of the file ice_evolu for the output of 
    437        !!      the temporal evolution of key variables 
    438        !! 
    439        !! ** input   : Namelist namicedia 
    440        !! 
    441        !! history : 
    442        !!  8.5  ! 03-08 (C. Ethe) original code 
    443        !!  9.0  ! 08-03 (M. Vancoppenolle) LIM3 
    444        !!------------------------------------------------------------------- 
    445        NAMELIST/namicedia/fmtinf, nfrinf, ninfo, ntmoy 
    446  
    447        INTEGER  ::   jv   ,     &  ! dummy loop indice 
    448           &          ntot ,     & 
    449           &          ndeb ,     & 
    450           &          irecl 
    451  
    452        REAL(wp) ::   zxx0, zxx1    ! temporary scalars 
    453  
    454        CHARACTER(len=jpchinf) ::   titinf 
    455        CHARACTER(len=50)      ::   clname 
    456        !!------------------------------------------------------------------- 
    457  
    458  
    459        ! Read Namelist namicedia 
    460        REWIND ( numnam_ice ) 
    461        READ   ( numnam_ice  , namicedia ) 
    462        IF(lwp) THEN 
    463           WRITE(numout,*) 
    464           WRITE(numout,*) 'lim_dia_init : ice parameters for ice diagnostics ' 
    465           WRITE(numout,*) '~~~~~~~~~~~~' 
    466           WRITE(numout,*) '   format of the output values                                 fmtinf = ', fmtinf 
    467           WRITE(numout,*) '   number of variables written in one line                     nfrinf = ', nfrinf  
    468           WRITE(numout,*) '   Instantaneous values of ice evolution or averaging          ntmoy  = ', ntmoy 
    469           WRITE(numout,*) '   frequency of ouputs on file ice_evolu in case of averaging  ninfo  = ', ninfo 
    470        ENDIF 
    471  
    472        ! masked grid cell area 
    473        aire(:,:) = area(:,:) * tms(:,:) 
    474  
    475        ! Titles of ice key variables : 
    476        titvar(1) = 'NoIt'  ! iteration number 
    477        titvar(2) = 'T yr'  ! time step in years 
    478        nbvt = 2            ! number of time variables 
    479  
    480        titvar(3) = 'AI_N'  ! sea ice area in the northern Hemisp.(10^12 km2) 
    481        titvar(4) = 'AI_S'  ! sea ice area in the southern Hemisp.(10^12 km2) 
    482        titvar(5) = 'EI_N'  ! sea ice extent (15%) in the northern Hemisp.(10^12 km2) 
    483        titvar(6) = 'EI_S'  ! sea ice extent (15%) in the southern Hemisp.(10^12 km2) 
    484        titvar(7) = 'VI_N'  ! sea ice volume in the northern Hemisp.(10^3 km3) 
    485        titvar(8) = 'VI_S'  ! sea ice volume in the southern Hemisp.(10^3 km3) 
    486        titvar(9) = 'VS_N'  ! snow volume over sea ice in the northern Hemisp.(10^3 km3) 
    487        titvar(10)= 'VS_S'  ! snow volume over sea ice in the northern Hemisp.(10^3 km3) 
    488        titvar(11)= 'VuIN'  ! undeformed sea ice volume in the northern Hemisp.(10^3 km3) 
    489        titvar(12)= 'VuIS'  ! undeformed sea ice volume in the southern Hemisp.(10^3 km3) 
    490        titvar(13)= 'VdIN'  ! deformed sea ice volume in the northern Hemisp.(10^3 km3) 
    491        titvar(14)= 'VdIS'  ! deformed sea ice volume in the southern Hemisp.(10^3 km3) 
    492        titvar(15)= 'OI_N'  ! sea ice mean age in the northern Hemisp.(years) 
    493        titvar(16)= 'OI_S'  ! sea ice mean age in the southern Hemisp.(years) 
    494        titvar(17)= 'AFYN'  ! total FY ice area northern Hemisp.(10^12 km2) 
    495        titvar(18)= 'AFYS'  ! total FY ice area southern Hemisp.(10^12 km2) 
    496        titvar(19)= 'AMYN'  ! total MY ice area northern Hemisp.(10^12 km2) 
    497        titvar(20)= 'AMYS'  ! total MY ice area southern Hemisp.(10^12 km2) 
    498        titvar(21)= 'EFYN'  ! total FY ice extent northern Hemisp.(10^12 km2) (with more 50% FY ice) 
    499        titvar(22)= 'EFYS'  ! total FY ice extent southern Hemisp.(10^12 km2) (with more 50% FY ice) 
    500        titvar(23)= 'EMYN'  ! total MY ice extent northern Hemisp.(10^12 km2) (with more 50% MY ice) 
    501        titvar(24)= 'EMYS'  ! total MY ice extent southern Hemisp.(10^12 km2) (with more 50% MY ice) 
    502        titvar(25)= 'VFYN'  ! total undeformed FY ice volume northern Hemisp.(10^3 km3)  
    503        titvar(26)= 'VFYS'  ! total undeformed FY ice volume southern Hemisp.(10^3 km3) 
    504        titvar(27)= 'VMYN'  ! total undeformed MY ice volume northern Hemisp.(10^3 km3)  
    505        titvar(28)= 'VMYS'  ! total undeformed MY ice volume southern Hemisp.(10^3 km3)  
    506        titvar(29)= 'IS_N'  ! sea ice mean salinity in the northern hemisphere (ppt)   
    507        titvar(30)= 'IS_S'  ! sea ice mean salinity in the southern hemisphere (ppt)   
    508        titvar(31)= 'IVeN'  ! sea ice mean velocity in the northern hemisphere (m/s)  
    509        titvar(32)= 'IVeS'  ! sea ice mean velocity in the southern hemisphere (m/s)  
    510        titvar(33)= 'DVDN'  ! variation of sea ice volume due to dynamics in the northern hemisphere 
    511        titvar(34)= 'DVDS'  ! variation of sea ice volume due to dynamics in the southern hemisphere 
    512        titvar(35)= 'DVTN'  ! variation of sea ice volume due to thermo in the   northern hemisphere 
    513        titvar(36)= 'DVTS'  ! variation of sea ice volume due to thermo in the   southern hemisphere 
    514        titvar(37)= 'TG1N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 1   
    515        titvar(38)= 'TG1S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 1   
    516        titvar(39)= 'TG2N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 2   
    517        titvar(40)= 'TG2S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 2   
    518        titvar(41)= 'TG3N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 3   
    519        titvar(42)= 'TG3S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 3   
    520        titvar(43)= 'TG4N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 4   
    521        titvar(44)= 'TG4S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 4   
    522        titvar(45)= 'TG5N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 5   
    523        titvar(46)= 'TG5S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 5   
    524        titvar(47)= 'LA_N'  ! lateral accretion growth rate, northern hemisphere 
    525        titvar(48)= 'LA_S'  ! lateral accretion growth rate, southern hemisphere  
    526        titvar(49)= 'SF_N'  ! Salinity FY, NH  
    527        titvar(50)= 'SF_S'  ! Salinity FY, SH  
    528        titvar(51)= 'SF_N'  ! Salinity MY, NH  
    529        titvar(52)= 'SF_S'  ! Salinity MY, SH  
    530        titvar(53)= 'Fs_N'  ! Total salt flux NH 
    531        titvar(54)= 'Fs_S'  ! Total salt flux SH 
    532        titvar(55)= 'FsbN'  ! Salt - brine drainage flux NH 
    533        titvar(56)= 'FsbS'  ! Salt - brine drainage flux SH 
    534        titvar(57)= 'FseN'  ! Salt - Equivalent salt flux NH 
    535        titvar(58)= 'FseS'  ! Salt - Equivalent salt flux SH 
    536        titvar(59)= 'SSTN'  ! SST, NH 
    537        titvar(60)= 'SSTS'  ! SST, SH 
    538        titvar(61)= 'SSSN'  ! SSS, NH 
    539        titvar(62)= 'SSSS'  ! SSS, SH 
    540        titvar(63)= 'TsuN'  ! Tsu, NH 
    541        titvar(64)= 'TsuS'  ! Tsu, SH 
    542        titvar(65)= 'TsnN'  ! Tsn, NH 
    543        titvar(66)= 'TsnS'  ! Tsn, SH 
    544        titvar(67)= 'ei_N'  ! ei, NH 
    545        titvar(68)= 'ei_S'  ! ei, SH 
    546        titvar(69)= 'vi1N'  ! vi1, NH 
    547        titvar(70)= 'vi1S'  ! vi1, SH 
    548        titvar(71)= 'vi2N'  ! vi2, NH 
    549        titvar(72)= 'vi2S'  ! vi2, SH 
    550        titvar(73)= 'vi3N'  ! vi3, NH 
    551        titvar(74)= 'vi3S'  ! vi3, SH 
    552        titvar(75)= 'vi4N'  ! vi4, NH 
    553        titvar(76)= 'vi4S'  ! vi4, SH 
    554        titvar(77)= 'vi5N'  ! vi5, NH 
    555        titvar(78)= 'vi5S'  ! vi5, SH 
    556        titvar(79)= 'vi6N'  ! vi6, NH 
    557        titvar(80)= 'vi6S'  ! vi6, SH 
    558        titvar(81)= 'fmaN'  ! mass flux in the ocean, NH 
    559        titvar(82)= 'fmaS'  ! mass flux in the ocean, SH 
    560        titvar(83)= 'AFSE'  ! Fram Strait Area export 
    561        titvar(84)= 'VFSE'  ! Fram Strait Volume export 
    562        nvinfo = 84 
    563  
    564        ! Definition et Ecriture de l'entete : nombre d'enregistrements  
    565        ndeb   = ( nstart - 1 ) / ninfo 
    566        IF( nstart == 1 ) ndeb = -1 
    567  
    568        nferme = ( nstart - 1 + nitrun) / ninfo 
    569        ntot   = nferme - ndeb 
    570        ndeb   = ninfo * ( 1 + ndeb ) 
    571        nferme = ninfo * nferme 
    572  
    573        ! definition of formats  
    574        WRITE( fmtw  , '(A,I3,A2,I1,A)' )  '(', nfrinf, '(A', jpchsep, ','//fmtinf//'))' 
    575        WRITE( fmtr  , '(A,I3,A,I1,A)'  )  '(', nfrinf, '(', jpchsep, 'X,'//fmtinf//'))' 
    576        WRITE( fmtitr, '(A,I3,A,I1,A)'  )  '(', nvinfo, 'A', jpchinf, ')' 
    577  
    578        ! opening  "ice_evolu" file 
    579        clname = 'ice.evolu' 
    580        irecl = ( jpchinf + 1 ) * nvinfo  
    581        CALL ctlopn( numevo_ice, clname, 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL',    & 
    582           &         irecl, numout, lwp, 1 ) 
    583  
    584        !- ecriture de 2 lignes d''entete : 
    585        WRITE(numevo_ice,1000) fmtr, fmtw, fmtitr, nvinfo, ntot, 0, nfrinf 
    586        zxx0 = 0.001 * REAL(ninfo) 
    587        zxx1 = 0.001 * REAL(ndeb) 
    588        WRITE(numevo_ice,1111) REAL(jpchinf), 0., zxx1, zxx0, 0., 0., 0 
    589  
    590        !- ecriture de 2 lignes de titre : 
    591        WRITE(numevo_ice,'(A,I8,A,I8,A,I5)')                                      & 
    592           'Evolution chronologique - Experience '//cexper   & 
    593           //'   de', ndeb, ' a', nferme, ' pas', ninfo 
    594        WRITE(numevo_ice,fmtitr) ( titvar(jv), jv = 1, nvinfo ) 
    595  
    596  
    597        !--preparation de "titvar" pour l''ecriture parmi les valeurs numeriques : 
    598        DO  jv = 2 , nvinfo 
    599           titinf     = titvar(jv)(:jpchinf) 
    600           titvar(jv) = '  '//titinf 
    601        END DO 
    602  
    603        !--Initialisation of the arrays for the accumulation 
    604        DO  jv = 1, nvinfo 
    605           vinfom(jv) = 0. 
    606        END DO 
    607        naveg = 0 
    608  
    609 1000   FORMAT( 3(A20),4(1x,I6) ) 
    610 1111   FORMAT( 3(F7.1,1X,F7.3,1X),I3,A )   
    611  
    612     END SUBROUTINE lim_dia_init 
     97      INTEGER  ::   jv,ji,jj,jl ! dummy loop indices 
     98      REAL(wp), DIMENSION(jpinfmx) ::  &  
     99         vinfor           ! temporary working space  
     100      REAL(wp) ::    & 
     101         zshift_date   , & ! date from the minimum ice extent 
     102         zday, zday_min, & ! current day, day of minimum extent 
     103         zafy, zamy,     & ! temporary area of fy and my ice 
     104         zindb 
     105      !!------------------------------------------------------------------- 
     106 
     107      ! 0) date from the minimum of ice extent 
     108      !--------------------------------------- 
     109      zday_min = 273.0        ! zday_min = date of minimum extent, here September 30th 
     110      zday = FLOAT(numit-nit000) * rdt_ice / ( 86400.0 * FLOAT(nn_fsbc) ) 
     111      IF (zday.GT.zday_min) THEN  
     112         zshift_date  =  zday - zday_min 
     113      ELSE 
     114         zshift_date  =  zday - (365.0 - zday_min) 
     115      ENDIF 
     116 
     117      IF( numit == nstart )   CALL lim_dia_init   ! initialisation of ice_evolu file       
     118 
     119      ! temporal diagnostics  
     120      vinfor(1) = REAL(numit) 
     121      vinfor(2) = nyear 
     122 
     123      ! put everything to zero 
     124      DO jv = nbvt + 1, nvinfo 
     125         vinfor(jv) = 0.0 
     126      END DO 
     127 
     128      !!------------------------------------------------------------------- 
     129      !! 1) Northern hemisphere 
     130      !!------------------------------------------------------------------- 
     131      !! 1.1) Diagnostics independent on age 
     132      !!------------------------------------ 
     133      DO jj = njeq, jpjm1 
     134         DO ji = fs_2, fs_jpim1   ! vector opt. 
     135            IF( tms(ji,jj) == 1 ) THEN 
     136               vinfor(3)  = vinfor(3)  + at_i(ji,jj)*aire(ji,jj) / 1.0e12 !ice area 
     137               IF (at_i(ji,jj).GT.0.15) vinfor(5) = vinfor(5) + aire(ji,jj) / 1.0e12 !ice extent 
     138               vinfor(7)  = vinfor(7)  + vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !ice volume 
     139               vinfor(9)  = vinfor(9)  + vt_s(ji,jj)*aire(ji,jj) / 1.0e12 !snow volume 
     140               vinfor(15) = vinfor(15) + ot_i(ji,jj) *vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !mean age 
     141               vinfor(29) = vinfor(29) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !mean salinity 
     142               ! the computation of this diagnostic is not reliable 
     143               vinfor(31) = vinfor(31) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + &  
     144                  v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12  
     145               vinfor(53) = vinfor(53) + emps(ji,jj)*aire(ji,jj) / 1.0e12 !salt flux 
     146               vinfor(55) = vinfor(55) + fsbri(ji,jj)*aire(ji,jj) / 1.0e12 !brine drainage flux 
     147               vinfor(57) = vinfor(57) + fseqv(ji,jj)*aire(ji,jj) / 1.0e12 !equivalent salt flux 
     148               vinfor(59) = vinfor(59) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SST 
     149               vinfor(61) = vinfor(61) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SSS 
     150               vinfor(65) = vinfor(65) + et_s(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12  ! snow temperature 
     151               vinfor(67) = vinfor(67) + et_i(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12       ! ice heat content 
     152               vinfor(69) = vinfor(69) + v_i(ji,jj,1)*aire(ji,jj) / 1.0e12 !ice volume 
     153               vinfor(71) = vinfor(71) + v_i(ji,jj,2)*aire(ji,jj) / 1.0e12 !ice volume 
     154               vinfor(73) = vinfor(73) + v_i(ji,jj,3)*aire(ji,jj) / 1.0e12 !ice volume 
     155               vinfor(75) = vinfor(75) + v_i(ji,jj,4)*aire(ji,jj) / 1.0e12 !ice volume 
     156               vinfor(77) = vinfor(77) + v_i(ji,jj,5)*aire(ji,jj) / 1.0e12 !ice volume 
     157               vinfor(79) = 0.0 
     158               vinfor(81) = vinfor(81) + emp(ji,jj)*aire(ji,jj) / 1.0e12 ! mass flux 
     159            ENDIF 
     160         END DO 
     161      END DO 
     162 
     163      DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 
     164         DO jj = njeq, jpjm1 
     165            DO ji = fs_2, fs_jpim1   ! vector opt. 
     166               IF( tms(ji,jj) == 1 ) THEN 
     167                  vinfor(11) = vinfor(11) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !undef def ice volume 
     168               ENDIF 
     169            END DO 
     170         END DO 
     171      END DO 
     172 
     173      vinfor(13) = 0.0 
     174 
     175      vinfor(15) = vinfor(15) / MAX(vinfor(7),epsi06) ! these have to be divided by total ice volume to have the 
     176      vinfor(29) = vinfor(29) / MAX(vinfor(7),epsi06) ! right value 
     177      vinfor(31) = SQRT( vinfor(31) / MAX( vinfor(7) , epsi06 ) ) 
     178      vinfor(67) = vinfor(67) / MAX(vinfor(7),epsi06) 
     179 
     180      vinfor(53) = vinfor(53) / MAX(vinfor(5),epsi06) ! these have to be divided by total ice extent to have the 
     181      vinfor(55) = vinfor(55) / MAX(vinfor(5),epsi06) ! right value  
     182      vinfor(57) = vinfor(57) / MAX(vinfor(5),epsi06) !  
     183      vinfor(79) = vinfor(79) / MAX(vinfor(5),epsi06) ! 
     184 
     185      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(3))) ! 
     186      vinfor(59) = zindb*vinfor(59) / MAX(vinfor(3),epsi06) ! divide by ice area 
     187      vinfor(61) = zindb*vinfor(61) / MAX(vinfor(3),epsi06) ! 
     188 
     189      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(9))) ! 
     190      vinfor(65) = zindb*vinfor(65) / MAX(vinfor(9),epsi06) ! divide it by snow volume 
     191 
     192 
     193      DO jl = 1, jpl 
     194         DO jj = njeq, jpjm1 
     195            DO ji = fs_2, fs_jpim1   ! vector opt. 
     196               IF( tms(ji,jj) == 1 ) THEN 
     197                  vinfor(33) = vinfor(33) + d_v_i_trp(ji,jj,jl)*aire(ji,jj) / 1.0e12 !ice volume 
     198                  vinfor(35) = vinfor(35) + d_v_i_thd(ji,jj,jl)*aire(ji,jj) / 1.0e12 !ice volume 
     199               ENDIF 
     200            END DO 
     201         END DO 
     202      END DO 
     203 
     204      DO jj = njeq, jpjm1 
     205         DO ji = fs_2, fs_jpim1   ! vector opt. 
     206            IF( tms(ji,jj) == 1 ) THEN 
     207               vinfor(37) = vinfor(37) + diag_sni_gr(ji,jj)*aire(ji,jj) / 1.0e12 !th growth rates 
     208               vinfor(39) = vinfor(39) + diag_lat_gr(ji,jj)*aire(ji,jj) / 1.0e12  
     209               vinfor(41) = vinfor(41) + diag_bot_gr(ji,jj)*aire(ji,jj) / 1.0e12 
     210               vinfor(43) = vinfor(43) + diag_dyn_gr(ji,jj)*aire(ji,jj) / 1.0e12  
     211               vinfor(45) = vinfor(45) + dv_dt_thd(ji,jj,5)*aire(ji,jj) / 1.0e12 
     212               vinfor(47) = vinfor(47) + v_newice(ji,jj) *aire(ji,jj) / 1.0e12 / rdt_ice ! volume acc in OW 
     213            ENDIF 
     214         END DO 
     215      END DO 
     216 
     217      DO jl = 1, jpl 
     218         DO jj = njeq, jpjm1 
     219            DO ji = fs_2, fs_jpim1   ! vector opt. 
     220               IF( tms(ji,jj) == 1 ) THEN 
     221                  vinfor(63) = vinfor(63) + t_su(ji,jj,jl)*a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 
     222               ENDIF 
     223            END DO 
     224         END DO 
     225      END DO 
     226      vinfor(63) = vinfor(63) / MAX(vinfor(3),epsi06) ! these have to be divided by total ice area 
     227 
     228      !! 1.2) Diagnostics dependent on age 
     229      !!------------------------------------ 
     230      DO jj = njeq, jpjm1 
     231         DO ji = fs_2, fs_jpim1   ! vector opt. 
     232            IF( tms(ji,jj) == 1 ) THEN 
     233               zafy = 0.0 
     234               zamy = 0.0 
     235               DO jl = 1, jpl 
     236                  IF ((o_i(ji,jj,jl) - zshift_date).LT.0.0) THEN 
     237                     vinfor(17) = vinfor(17) + a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! FY ice area 
     238                     vinfor(25) = vinfor(25) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! FY ice volume 
     239                     vinfor(49) = vinfor(49) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !FY ice salinity 
     240                     zafy = zafy + a_i(ji,jj,jl) 
     241                  ENDIF 
     242                  IF ((o_i(ji,jj,jl) - zshift_date).GT.0.0) THEN 
     243                     vinfor(19) = vinfor(19) + a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12    ! MY ice area 
     244                     vinfor(27) = vinfor(27) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! MY ice volume 
     245                     vinfor(51) = vinfor(51) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !MY ice salinity 
     246                     zamy = zamy + a_i(ji,jj,jl) 
     247                  ENDIF 
     248               END DO 
     249               IF ((at_i(ji,jj).GT.0.15).AND.(zafy.GT.zamy)) THEN 
     250                  vinfor(21) = vinfor(21) + aire(ji,jj) / 1.0e12 ! Seasonal ice extent 
     251               ENDIF 
     252               IF ((at_i(ji,jj).GT.0.15).AND.(zafy.LE.zamy)) THEN 
     253                  vinfor(23) = vinfor(23) + aire(ji,jj) / 1.0e12 ! Perennial ice extent 
     254               ENDIF 
     255            ENDIF 
     256         END DO 
     257      END DO 
     258      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(25))) !=0 if no multiyear ice 1 if yes 
     259      vinfor(49) = zindb*vinfor(49) / MAX(vinfor(25),epsi06) 
     260      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(27))) !=0 if no multiyear ice 1 if yes 
     261      vinfor(51) = zindb*vinfor(51) / MAX(vinfor(27),epsi06) 
     262 
     263      !! Fram Strait Export 
     264      !! 83 = area export 
     265      !! 84 = volume export 
     266      !! Fram strait in ORCA2 = 5 points 
     267      !! export = -v_ice*e1t*ddtb*at_i or -v_ice*e1t*ddtb*at_i*h_i 
     268      jj = 136 ! C grid 
     269      vinfor(83) = 0.0 
     270      vinfor(84) = 0.0 
     271      DO ji = 134, 138 
     272         vinfor(83) = vinfor(83) - v_ice(ji,jj) * &  
     273            e1t(ji,jj)*at_i(ji,jj)*rdt_ice / 1.0e12 
     274         vinfor(84) = vinfor(84) - v_ice(ji,jj) * &  
     275            e1t(ji,jj)*vt_i(ji,jj)*rdt_ice / 1.0e12 
     276      END DO 
     277 
     278      !!------------------------------------------------------------------- 
     279      !! 2) Southern hemisphere 
     280      !!------------------------------------------------------------------- 
     281      !! 2.1) Diagnostics independent on age 
     282      !!------------------------------------ 
     283      DO jj = 2, njeqm1 
     284         DO ji = fs_2, fs_jpim1   ! vector opt. 
     285            IF( tms(ji,jj) == 1 ) THEN 
     286               vinfor(4)  = vinfor(4)  + at_i(ji,jj)*aire(ji,jj) / 1.0e12 !ice area 
     287               IF (at_i(ji,jj).GT.0.15) vinfor(6) = vinfor(6) + aire(ji,jj) / 1.0e12 !ice extent 
     288               vinfor(8)  = vinfor(8)  + vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !ice volume 
     289               vinfor(10) = vinfor(10) + vt_s(ji,jj)*aire(ji,jj) / 1.0e12 !snow volume 
     290               vinfor(16) = vinfor(16) + ot_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !mean age 
     291               vinfor(30) = vinfor(30) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) / 1.0e12 !mean salinity 
     292               ! this diagnostic is not well computed (weighted by vol instead 
     293               ! of area) 
     294               vinfor(32) = vinfor(32) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + &  
     295                  v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12 !ice vel 
     296               vinfor(54) = vinfor(54) + at_i(ji,jj)*emps(ji,jj)*aire(ji,jj) / 1.0e12 ! Total salt flux 
     297               vinfor(56) = vinfor(56) + at_i(ji,jj)*fsbri(ji,jj)*aire(ji,jj) / 1.0e12 ! Brine drainage salt flux 
     298               vinfor(58) = vinfor(58) + at_i(ji,jj)*fseqv(ji,jj)*aire(ji,jj) / 1.0e12 ! Equivalent salt flux 
     299               vinfor(60) = vinfor(60) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SST 
     300               vinfor(62) = vinfor(62) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SSS 
     301               vinfor(66) = vinfor(66) + et_s(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12 ! snow temperature 
     302               vinfor(68) = vinfor(68) + et_i(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12 ! ice enthalpy 
     303               vinfor(70) = vinfor(70) + v_i(ji,jj,1)*aire(ji,jj) / 1.0e12 !ice volume 
     304               vinfor(72) = vinfor(72) + v_i(ji,jj,2)*aire(ji,jj) / 1.0e12 !ice volume 
     305               vinfor(74) = vinfor(74) + v_i(ji,jj,3)*aire(ji,jj) / 1.0e12 !ice volume 
     306               vinfor(76) = vinfor(76) + v_i(ji,jj,4)*aire(ji,jj) / 1.0e12 !ice volume 
     307               vinfor(78) = vinfor(78) + v_i(ji,jj,5)*aire(ji,jj) / 1.0e12 !ice volume 
     308               vinfor(80) = 0.0 
     309               vinfor(82) = vinfor(82) + emp(ji,jj)*aire(ji,jj) / 1.0e12 ! mass flux 
     310            ENDIF 
     311         END DO 
     312      END DO 
     313 
     314      DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 
     315         DO jj = 2, njeqm1 
     316            DO ji = fs_2, fs_jpim1   ! vector opt. 
     317               vinfor(12) = vinfor(12) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !undef def ice volume 
     318            END DO 
     319         END DO 
     320      END DO 
     321 
     322      vinfor(14) = 0.0 
     323 
     324      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(8)))  
     325      vinfor(16) = zindb * vinfor(16) / MAX(vinfor(8),epsi06) ! these have to be divided by ice vol 
     326      vinfor(30) = zindb * vinfor(30) / MAX(vinfor(8),epsi06) !  
     327      vinfor(32) = zindb * SQRT( vinfor(32) / MAX( vinfor(8) , epsi06 ) ) 
     328      vinfor(68) = zindb * vinfor(68) / MAX(vinfor(8),epsi06) !  
     329 
     330      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(6)))  
     331      vinfor(54) = zindb * vinfor(54) / MAX(vinfor(6),epsi06) ! these have to be divided by ice extt 
     332      vinfor(56) = zindb * vinfor(56) / MAX(vinfor(6),epsi06) !  
     333      vinfor(58) = zindb * vinfor(58) / MAX(vinfor(6),epsi06) !  
     334      vinfor(80) = zindb * vinfor(80) / MAX(vinfor(6),epsi06) ! 
     335      !      vinfor(84) = vinfor(84) / vinfor(6) ! 
     336 
     337      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(4))) ! 
     338      vinfor(60) = zindb*vinfor(60) / ( MAX(vinfor(4), epsi06) ) ! divide by ice area 
     339      vinfor(62) = zindb*vinfor(62) / ( MAX(vinfor(4), epsi06) ) ! 
     340 
     341      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(10))) ! 
     342      vinfor(66) = zindb*vinfor(66) / MAX(vinfor(10),epsi06) ! divide it by snow volume 
     343 
     344      DO jl = 1, jpl 
     345         DO jj = 2, njeqm1 
     346            DO ji = fs_2, fs_jpim1   ! vector opt. 
     347               IF( tms(ji,jj) == 1 ) THEN 
     348                  vinfor(34) = vinfor(34) + d_v_i_trp(ji,jj,jl)*aire(ji,jj) / 1.0e12 !ice volume 
     349                  vinfor(36) = vinfor(36) + d_v_i_thd(ji,jj,jl)*aire(ji,jj) / 1.0e12 !ice volume 
     350               ENDIF 
     351            END DO 
     352         END DO 
     353      END DO 
     354 
     355      DO jj = 2, njeqm1 
     356         DO ji = fs_2, fs_jpim1   ! vector opt. 
     357            IF( tms(ji,jj) == 1 ) THEN 
     358               vinfor(38) = vinfor(38) + diag_sni_gr(ji,jj)*aire(ji,jj) / 1.0e12 !th growth rates 
     359               vinfor(40) = vinfor(40) + diag_lat_gr(ji,jj)*aire(ji,jj) / 1.0e12  
     360               vinfor(42) = vinfor(42) + diag_bot_gr(ji,jj)*aire(ji,jj) / 1.0e12 
     361               vinfor(44) = vinfor(44) + diag_dyn_gr(ji,jj)*aire(ji,jj) / 1.0e12  
     362               vinfor(46) = vinfor(46) + dv_dt_thd(ji,jj,5)*aire(ji,jj) / 1.0e12 
     363               vinfor(48) = vinfor(48) + v_newice(ji,jj) *aire(ji,jj) / 1.0e12 / rdt_ice ! volume acc in OW 
     364            ENDIF 
     365         END DO 
     366      END DO 
     367 
     368 
     369      DO jl = 1, jpl 
     370         DO jj = 2, njeqm1 
     371            DO ji = fs_2, fs_jpim1   ! vector opt. 
     372               IF( tms(ji,jj) == 1 ) THEN 
     373                  vinfor(64) = vinfor(64) + t_su(ji,jj,jl)*a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 
     374               ENDIF 
     375            END DO 
     376         END DO 
     377      END DO 
     378      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(4))) ! 
     379      vinfor(64) = zindb * vinfor(64) / MAX(vinfor(4),epsi06) ! divide by ice extt 
     380      !! 2.2) Diagnostics dependent on age 
     381      !!------------------------------------ 
     382      DO jj = 2, njeqm1 
     383         DO ji = fs_2, fs_jpim1   ! vector opt. 
     384            IF( tms(ji,jj) == 1 ) THEN 
     385               zafy = 0.0 
     386               zamy = 0.0 
     387               DO jl = 1, jpl 
     388                  IF ((o_i(ji,jj,jl) - zshift_date).LT.0.0) THEN 
     389                     vinfor(18) = vinfor(18) + a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! FY ice area 
     390                     vinfor(26) = vinfor(26) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 ! FY ice volume 
     391                     zafy = zafy + a_i(ji,jj,jl) 
     392                     vinfor(50) = vinfor(50) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !FY ice salinity 
     393                  ENDIF 
     394                  IF ((o_i(ji,jj,jl) - zshift_date).GT.0.0) THEN 
     395                     vinfor(20) = vinfor(20) + a_i(ji,jj,jl)*aire(ji,jj) / 1.0e12    ! MY ice area 
     396                     vinfor(28) = vinfor(28) + v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 
     397                     vinfor(52) = vinfor(52) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) / 1.0e12 !FY ice salinity 
     398                     zamy = zamy + a_i(ji,jj,jl) 
     399                  ENDIF 
     400               END DO ! jl 
     401               IF ((at_i(ji,jj).GT.0.15).AND.(zafy.GT.zamy)) THEN 
     402                  vinfor(22) = vinfor(22) + aire(ji,jj) / 1.0e12 ! Seasonal ice extent 
     403               ENDIF 
     404               IF ((at_i(ji,jj).GT.0.15).AND.(zafy.LE.zamy)) THEN 
     405                  vinfor(24) = vinfor(24) + aire(ji,jj) / 1.0e12 ! Perennial ice extent 
     406               ENDIF 
     407            ENDIF ! tms 
     408         END DO ! jj 
     409      END DO ! ji 
     410      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(26))) !=0 if no multiyear ice 1 if yes 
     411      vinfor(50) = zindb*vinfor(50) / MAX(vinfor(26),epsi06) 
     412      zindb      = 1.0 - MAX(0.0,SIGN(1.0,-vinfor(28))) !=0 if no multiyear ice 1 if yes 
     413      vinfor(52) = zindb*vinfor(52) / MAX(vinfor(28),epsi06) 
     414 
     415      !  Accumulation before averaging  
     416      DO jv = 1, nvinfo 
     417         vinfom(jv) = vinfom(jv) + vinfor(jv) 
     418      END DO 
     419      naveg = naveg + 1   
     420 
     421      ! oututs on file ice_evolu     
     422      !MV      IF( MOD( numit , ninfo ) == 0 ) THEN 
     423      WRITE(numevo_ice,fmtw) ( titvar(jv), vinfom(jv)/naveg, jv = 1, nvinfo ) 
     424      naveg = 0 
     425      DO jv = 1, nvinfo 
     426         vinfom(jv)=0.0 
     427      END DO 
     428      !MV      ENDIF 
     429 
     430   END SUBROUTINE lim_dia 
     431 
     432   SUBROUTINE lim_dia_init 
     433      !!------------------------------------------------------------------- 
     434      !!                  ***  ROUTINE lim_dia_init  *** 
     435      !!              
     436      !! ** Purpose : Preparation of the file ice_evolu for the output of 
     437      !!      the temporal evolution of key variables 
     438      !! 
     439      !! ** input   : Namelist namicedia 
     440      !! 
     441      !! history : 
     442      !!  8.5  ! 03-08 (C. Ethe) original code 
     443      !!  9.0  ! 08-03 (M. Vancoppenolle) LIM3 
     444      !!------------------------------------------------------------------- 
     445      NAMELIST/namicedia/fmtinf, nfrinf, ninfo, ntmoy 
     446 
     447      INTEGER  ::   jv   ,     &  ! dummy loop indice 
     448         &          ntot ,     & 
     449         &          ndeb ,     & 
     450         &          irecl 
     451 
     452      REAL(wp) ::   zxx0, zxx1    ! temporary scalars 
     453 
     454      CHARACTER(len=jpchinf) ::   titinf 
     455      CHARACTER(len=50)      ::   clname 
     456      !!------------------------------------------------------------------- 
     457 
     458 
     459      ! Read Namelist namicedia 
     460      REWIND ( numnam_ice ) 
     461      READ   ( numnam_ice  , namicedia ) 
     462      IF(lwp) THEN 
     463         WRITE(numout,*) 
     464         WRITE(numout,*) 'lim_dia_init : ice parameters for ice diagnostics ' 
     465         WRITE(numout,*) '~~~~~~~~~~~~' 
     466         WRITE(numout,*) '   format of the output values                                 fmtinf = ', fmtinf 
     467         WRITE(numout,*) '   number of variables written in one line                     nfrinf = ', nfrinf  
     468         WRITE(numout,*) '   Instantaneous values of ice evolution or averaging          ntmoy  = ', ntmoy 
     469         WRITE(numout,*) '   frequency of ouputs on file ice_evolu in case of averaging  ninfo  = ', ninfo 
     470      ENDIF 
     471 
     472      ! masked grid cell area 
     473      aire(:,:) = area(:,:) * tms(:,:) 
     474 
     475      ! Titles of ice key variables : 
     476      titvar(1) = 'NoIt'  ! iteration number 
     477      titvar(2) = 'T yr'  ! time step in years 
     478      nbvt = 2            ! number of time variables 
     479 
     480      titvar(3) = 'AI_N'  ! sea ice area in the northern Hemisp.(10^12 km2) 
     481      titvar(4) = 'AI_S'  ! sea ice area in the southern Hemisp.(10^12 km2) 
     482      titvar(5) = 'EI_N'  ! sea ice extent (15%) in the northern Hemisp.(10^12 km2) 
     483      titvar(6) = 'EI_S'  ! sea ice extent (15%) in the southern Hemisp.(10^12 km2) 
     484      titvar(7) = 'VI_N'  ! sea ice volume in the northern Hemisp.(10^3 km3) 
     485      titvar(8) = 'VI_S'  ! sea ice volume in the southern Hemisp.(10^3 km3) 
     486      titvar(9) = 'VS_N'  ! snow volume over sea ice in the northern Hemisp.(10^3 km3) 
     487      titvar(10)= 'VS_S'  ! snow volume over sea ice in the northern Hemisp.(10^3 km3) 
     488      titvar(11)= 'VuIN'  ! undeformed sea ice volume in the northern Hemisp.(10^3 km3) 
     489      titvar(12)= 'VuIS'  ! undeformed sea ice volume in the southern Hemisp.(10^3 km3) 
     490      titvar(13)= 'VdIN'  ! deformed sea ice volume in the northern Hemisp.(10^3 km3) 
     491      titvar(14)= 'VdIS'  ! deformed sea ice volume in the southern Hemisp.(10^3 km3) 
     492      titvar(15)= 'OI_N'  ! sea ice mean age in the northern Hemisp.(years) 
     493      titvar(16)= 'OI_S'  ! sea ice mean age in the southern Hemisp.(years) 
     494      titvar(17)= 'AFYN'  ! total FY ice area northern Hemisp.(10^12 km2) 
     495      titvar(18)= 'AFYS'  ! total FY ice area southern Hemisp.(10^12 km2) 
     496      titvar(19)= 'AMYN'  ! total MY ice area northern Hemisp.(10^12 km2) 
     497      titvar(20)= 'AMYS'  ! total MY ice area southern Hemisp.(10^12 km2) 
     498      titvar(21)= 'EFYN'  ! total FY ice extent northern Hemisp.(10^12 km2) (with more 50% FY ice) 
     499      titvar(22)= 'EFYS'  ! total FY ice extent southern Hemisp.(10^12 km2) (with more 50% FY ice) 
     500      titvar(23)= 'EMYN'  ! total MY ice extent northern Hemisp.(10^12 km2) (with more 50% MY ice) 
     501      titvar(24)= 'EMYS'  ! total MY ice extent southern Hemisp.(10^12 km2) (with more 50% MY ice) 
     502      titvar(25)= 'VFYN'  ! total undeformed FY ice volume northern Hemisp.(10^3 km3)  
     503      titvar(26)= 'VFYS'  ! total undeformed FY ice volume southern Hemisp.(10^3 km3) 
     504      titvar(27)= 'VMYN'  ! total undeformed MY ice volume northern Hemisp.(10^3 km3)  
     505      titvar(28)= 'VMYS'  ! total undeformed MY ice volume southern Hemisp.(10^3 km3)  
     506      titvar(29)= 'IS_N'  ! sea ice mean salinity in the northern hemisphere (ppt)   
     507      titvar(30)= 'IS_S'  ! sea ice mean salinity in the southern hemisphere (ppt)   
     508      titvar(31)= 'IVeN'  ! sea ice mean velocity in the northern hemisphere (m/s)  
     509      titvar(32)= 'IVeS'  ! sea ice mean velocity in the southern hemisphere (m/s)  
     510      titvar(33)= 'DVDN'  ! variation of sea ice volume due to dynamics in the northern hemisphere 
     511      titvar(34)= 'DVDS'  ! variation of sea ice volume due to dynamics in the southern hemisphere 
     512      titvar(35)= 'DVTN'  ! variation of sea ice volume due to thermo in the   northern hemisphere 
     513      titvar(36)= 'DVTS'  ! variation of sea ice volume due to thermo in the   southern hemisphere 
     514      titvar(37)= 'TG1N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 1   
     515      titvar(38)= 'TG1S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 1   
     516      titvar(39)= 'TG2N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 2   
     517      titvar(40)= 'TG2S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 2   
     518      titvar(41)= 'TG3N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 3   
     519      titvar(42)= 'TG3S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 3   
     520      titvar(43)= 'TG4N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 4   
     521      titvar(44)= 'TG4S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 4   
     522      titvar(45)= 'TG5N'  ! thermodynamic vertical growth rate in the northern hemisphere, cat 5   
     523      titvar(46)= 'TG5S'  ! thermodynamic vertical growth rate in the souhtern hemisphere, cat 5   
     524      titvar(47)= 'LA_N'  ! lateral accretion growth rate, northern hemisphere 
     525      titvar(48)= 'LA_S'  ! lateral accretion growth rate, southern hemisphere  
     526      titvar(49)= 'SF_N'  ! Salinity FY, NH  
     527      titvar(50)= 'SF_S'  ! Salinity FY, SH  
     528      titvar(51)= 'SF_N'  ! Salinity MY, NH  
     529      titvar(52)= 'SF_S'  ! Salinity MY, SH  
     530      titvar(53)= 'Fs_N'  ! Total salt flux NH 
     531      titvar(54)= 'Fs_S'  ! Total salt flux SH 
     532      titvar(55)= 'FsbN'  ! Salt - brine drainage flux NH 
     533      titvar(56)= 'FsbS'  ! Salt - brine drainage flux SH 
     534      titvar(57)= 'FseN'  ! Salt - Equivalent salt flux NH 
     535      titvar(58)= 'FseS'  ! Salt - Equivalent salt flux SH 
     536      titvar(59)= 'SSTN'  ! SST, NH 
     537      titvar(60)= 'SSTS'  ! SST, SH 
     538      titvar(61)= 'SSSN'  ! SSS, NH 
     539      titvar(62)= 'SSSS'  ! SSS, SH 
     540      titvar(63)= 'TsuN'  ! Tsu, NH 
     541      titvar(64)= 'TsuS'  ! Tsu, SH 
     542      titvar(65)= 'TsnN'  ! Tsn, NH 
     543      titvar(66)= 'TsnS'  ! Tsn, SH 
     544      titvar(67)= 'ei_N'  ! ei, NH 
     545      titvar(68)= 'ei_S'  ! ei, SH 
     546      titvar(69)= 'vi1N'  ! vi1, NH 
     547      titvar(70)= 'vi1S'  ! vi1, SH 
     548      titvar(71)= 'vi2N'  ! vi2, NH 
     549      titvar(72)= 'vi2S'  ! vi2, SH 
     550      titvar(73)= 'vi3N'  ! vi3, NH 
     551      titvar(74)= 'vi3S'  ! vi3, SH 
     552      titvar(75)= 'vi4N'  ! vi4, NH 
     553      titvar(76)= 'vi4S'  ! vi4, SH 
     554      titvar(77)= 'vi5N'  ! vi5, NH 
     555      titvar(78)= 'vi5S'  ! vi5, SH 
     556      titvar(79)= 'vi6N'  ! vi6, NH 
     557      titvar(80)= 'vi6S'  ! vi6, SH 
     558      titvar(81)= 'fmaN'  ! mass flux in the ocean, NH 
     559      titvar(82)= 'fmaS'  ! mass flux in the ocean, SH 
     560      titvar(83)= 'AFSE'  ! Fram Strait Area export 
     561      titvar(84)= 'VFSE'  ! Fram Strait Volume export 
     562      nvinfo = 84 
     563 
     564      ! Definition et Ecriture de l'entete : nombre d'enregistrements  
     565      ndeb   = ( nstart - 1 ) / ninfo 
     566      IF( nstart == 1 ) ndeb = -1 
     567 
     568      nferme = ( nstart - 1 + nitrun) / ninfo 
     569      ntot   = nferme - ndeb 
     570      ndeb   = ninfo * ( 1 + ndeb ) 
     571      nferme = ninfo * nferme 
     572 
     573      ! definition of formats  
     574      WRITE( fmtw  , '(A,I3,A2,I1,A)' )  '(', nfrinf, '(A', jpchsep, ','//fmtinf//'))' 
     575      WRITE( fmtr  , '(A,I3,A,I1,A)'  )  '(', nfrinf, '(', jpchsep, 'X,'//fmtinf//'))' 
     576      WRITE( fmtitr, '(A,I3,A,I1,A)'  )  '(', nvinfo, 'A', jpchinf, ')' 
     577 
     578      ! opening  "ice_evolu" file 
     579      clname = 'ice.evolu' 
     580      irecl = ( jpchinf + 1 ) * nvinfo  
     581      CALL ctlopn( numevo_ice, clname, 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL',    & 
     582         &         irecl, numout, lwp, 1 ) 
     583 
     584      !- ecriture de 2 lignes d''entete : 
     585      WRITE(numevo_ice,1000) fmtr, fmtw, fmtitr, nvinfo, ntot, 0, nfrinf 
     586      zxx0 = 0.001 * REAL(ninfo) 
     587      zxx1 = 0.001 * REAL(ndeb) 
     588      WRITE(numevo_ice,1111) REAL(jpchinf), 0., zxx1, zxx0, 0., 0., 0 
     589 
     590      !- ecriture de 2 lignes de titre : 
     591      WRITE(numevo_ice,'(A,I8,A,I8,A,I5)')                                      & 
     592         'Evolution chronologique - Experience '//cexper   & 
     593         //'   de', ndeb, ' a', nferme, ' pas', ninfo 
     594      WRITE(numevo_ice,fmtitr) ( titvar(jv), jv = 1, nvinfo ) 
     595 
     596 
     597      !--preparation de "titvar" pour l''ecriture parmi les valeurs numeriques : 
     598      DO  jv = 2 , nvinfo 
     599         titinf     = titvar(jv)(:jpchinf) 
     600         titvar(jv) = '  '//titinf 
     601      END DO 
     602 
     603      !--Initialisation of the arrays for the accumulation 
     604      DO  jv = 1, nvinfo 
     605         vinfom(jv) = 0. 
     606      END DO 
     607      naveg = 0 
     608 
     6091000  FORMAT( 3(A20),4(1x,I6) ) 
     6101111  FORMAT( 3(F7.1,1X,F7.3,1X),I3,A )   
     611 
     612   END SUBROUTINE lim_dia_init 
    613613 
    614614#else 
Note: See TracChangeset for help on using the changeset viewer.