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.
limwri.F90 in branches/2016/v3_6_CMIP6_ice_diagnostics/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2016/v3_6_CMIP6_ice_diagnostics/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90 @ 8150

Last change on this file since 8150 was 8150, checked in by vancop, 7 years ago

SIMIP outputs, phase 2

  • Property svn:keywords set to Id
File size: 39.1 KB
Line 
1MODULE limwri
2   !!======================================================================
3   !!                     ***  MODULE  limwri  ***
4   !!         Ice diagnostics :  write ice output files
5   !!======================================================================
6#if defined key_lim3
7   !!----------------------------------------------------------------------
8   !!   'key_lim3'                                      LIM3 sea-ice model
9   !!----------------------------------------------------------------------
10   !!   lim_wri      : write of the diagnostics variables in ouput file
11   !!   lim_wri_state : write for initial state or/and abandon
12   !!----------------------------------------------------------------------
13   USE ioipsl
14   USE dianam          ! build name of file (routine)
15   USE phycst
16   USE dom_oce
17   USE sbc_oce         ! Surface boundary condition: ocean fields
18   USE sbc_ice         ! Surface boundary condition: ice fields
19   USE ice
20   USE dom_ice
21   USE limvar
22   USE in_out_manager
23   USE lbclnk
24   USE lib_mpp         ! MPP library
25   USE wrk_nemo        ! work arrays
26   USE iom
27   USE timing          ! Timing
28   USE lib_fortran     ! Fortran utilities
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC lim_wri        ! routine called by lim_step.F90
34   PUBLIC lim_wri_state  ! called by dia_wri_state
35
36   !!----------------------------------------------------------------------
37   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
38   !! $Id$
39   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
40   !!----------------------------------------------------------------------
41CONTAINS
42
43   SUBROUTINE lim_wri( kindic )
44      !!-------------------------------------------------------------------
45      !!  This routine computes the average of some variables and write it
46      !!  on the ouput files.
47      !!  ATTENTION cette routine n'est valable que si le pas de temps est
48      !!  egale a une fraction entiere de 1 jours.
49      !!  Diff 1-D 3-D : suppress common also included in etat
50      !!                 suppress cmoymo 11-18
51      !!  modif : 03/06/98
52      !!-------------------------------------------------------------------
53      INTEGER, INTENT(in) ::   kindic   ! if kindic < 0 there has been an error somewhere
54      !
55      INTEGER  ::  ii, ji, jj, jk, jl  ! dummy loop indices
56      REAL(wp) ::  z1_365
57      REAL(wp) ::  z2da, z2db, ztmp, zrho1, zrho2
58      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zswi2
59      REAL(wp), POINTER, DIMENSION(:,:)   ::  z2d, zswi    ! 2D workspace
60
61      ! Global ice diagnostics (SIMIP)
62      REAL(wp) ::  zdiag_area_nh, &   ! area, extent, volume
63         &         zdiag_extt_nh, &
64         &         zdiag_area_sh, & 
65         &         zdiag_extt_sh, & 
66         &         zdiag_volu_nh, & 
67         &         zdiag_volu_sh 
68
69      ! Strait / passage fluxes (SIMIP)
70      REAL(wp), DIMENSION(4) ::  &                                    ! Strait fluxes for output
71         &         zdiag_area_strait  ,   &                                    ! 1=Fram Strait, 2=CAA, 3= Barents, 4 = Bering
72         &         zdiag_mice_strait  ,   &
73         &         zdiag_msno_strait
74
75      REAL(wp) :: zfarea_u, zfmice_u, zfmsno_u, zfarea_v, zfmice_v, zfmsno_v   ! dummy fluxes
76
77      REAL(wp), DIMENSION(11) :: &
78         &         zui, zuj, zvi, zvj                                          ! strait addresses
79
80      INTEGER  :: Nu, Nv                                                       ! passage size
81
82      INTEGER, DIMENSION(4)  :: ji0, ji1, jj0, jj1
83     
84      !!-------------------------------------------------------------------
85
86      IF( nn_timing == 1 )  CALL timing_start('limwri')
87
88      CALL wrk_alloc( jpi, jpj, jpl, zswi2 )
89      CALL wrk_alloc( jpi, jpj     , z2d, zswi )
90
91      !-----------------------------
92      ! Mean category values
93      !-----------------------------
94      z1_365 = 1._wp / 365._wp
95
96      ! brine volume
97      CALL lim_var_bv 
98
99      ! tresholds for outputs
100      DO jj = 1, jpj
101         DO ji = 1, jpi
102            zswi(ji,jj)  = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice, 0 if no ice
103         END DO
104      END DO
105      DO jl = 1, jpl
106         DO jj = 1, jpj
107            DO ji = 1, jpi
108               zswi2(ji,jj,jl)  = MAX( 0._wp , SIGN( 1._wp , a_i(ji,jj,jl) - epsi06 ) )
109            END DO
110         END DO
111      END DO
112      !
113      ! fluxes
114      ! pfrld is the lead fraction at the previous time step (actually between TRP and THD)
115      IF( iom_use('qsr_oce') )   CALL iom_put( "qsr_oce" , qsr_oce(:,:) * pfrld(:,:) )                                   !     solar flux at ocean surface
116      IF( iom_use('qns_oce') )   CALL iom_put( "qns_oce" , qns_oce(:,:) * pfrld(:,:) + qemp_oce(:,:) )                   ! non-solar flux at ocean surface
117      IF( iom_use('qsr_ice') )   CALL iom_put( "qsr_ice" , SUM( qsr_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) )                 !     solar flux at ice surface
118      IF( iom_use('qns_ice') )   CALL iom_put( "qns_ice" , SUM( qns_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) + qemp_ice(:,:) ) ! non-solar flux at ice surface
119      IF( iom_use('qtr_ice') )   CALL iom_put( "qtr_ice" , SUM( ftr_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) )                 !     solar flux transmitted thru ice
120      IF( iom_use('qt_oce' ) )   CALL iom_put( "qt_oce"  , ( qsr_oce(:,:) + qns_oce(:,:) ) * pfrld(:,:) + qemp_oce(:,:) ) 
121      IF( iom_use('qt_ice' ) )   CALL iom_put( "qt_ice"  , SUM( ( qns_ice(:,:,:) + qsr_ice(:,:,:) )   &
122         &                                                      * a_i_b(:,:,:),dim=3 ) + qemp_ice(:,:) )
123      IF( iom_use('qemp_oce') )  CALL iom_put( "qemp_oce" , qemp_oce(:,:) ) 
124      IF( iom_use('qemp_ice') )  CALL iom_put( "qemp_ice" , qemp_ice(:,:) ) 
125      IF( iom_use('emp_oce' ) )  CALL iom_put( "emp_oce"  , emp_oce(:,:) )   !emp over ocean (taking into account the snow blown away from the ice)
126      IF( iom_use('emp_ice' ) )  CALL iom_put( "emp_ice"  , emp_ice(:,:) )   !emp over ice   (taking into account the snow blown away from the ice)
127
128      ! velocity
129      IF ( iom_use( "uice_ipa" ) .OR. iom_use( "vice_ipa" ) .OR. iom_use( "icevel" ) ) THEN
130         DO jj = 2 , jpjm1
131            DO ji = 2 , jpim1
132               z2da  = ( u_ice(ji,jj) * umask(ji,jj,1) + u_ice(ji-1,jj) * umask(ji-1,jj,1) ) * 0.5_wp
133               z2db  = ( v_ice(ji,jj) * vmask(ji,jj,1) + v_ice(ji,jj-1) * vmask(ji,jj-1,1) ) * 0.5_wp
134               z2d(ji,jj) = SQRT( z2da * z2da + z2db * z2db )
135           END DO
136         END DO
137         CALL lbc_lnk( z2d, 'T', 1. )
138         CALL iom_put( "uice_ipa"     , u_ice      )       ! ice velocity u component
139         CALL iom_put( "vice_ipa"     , v_ice      )       ! ice velocity v component
140         CALL iom_put( "icevel"       , z2d        )       ! ice velocity module
141      ENDIF
142      !
143      IF ( iom_use( "miceage" ) )       CALL iom_put( "miceage"     , om_i * zswi * z1_365   )  ! mean ice age
144      IF ( iom_use( "micet" ) )         CALL iom_put( "micet"       , ( tm_i  - rt0 ) * zswi )  ! ice mean    temperature
145      IF ( iom_use( "icest" ) )         CALL iom_put( "icest"       , ( tm_su - rt0 ) * zswi )  ! ice surface temperature
146      IF ( iom_use( "icecolf" ) )       CALL iom_put( "icecolf"     , hicol                  )  ! frazil ice collection thickness
147      !
148      CALL iom_put( "isst"        , sst_m               )        ! sea surface temperature
149      CALL iom_put( "isss"        , sss_m               )        ! sea surface salinity
150      CALL iom_put( "iceconc"     , at_i  * zswi        )        ! ice concentration
151      CALL iom_put( "icevolu"     , vt_i  * zswi        )        ! ice volume = mean ice thickness over the cell
152      CALL iom_put( "icehc"       , et_i  * zswi        )        ! ice total heat content
153      CALL iom_put( "isnowhc"     , et_s  * zswi        )        ! snow total heat content
154      CALL iom_put( "ibrinv"      , bvm_i * zswi * 100. )        ! brine volume
155      CALL iom_put( "utau_ice"    , utau_ice*zswi       )        ! wind stress over ice along i-axis at I-point
156      CALL iom_put( "vtau_ice"    , vtau_ice*zswi       )        ! wind stress over ice along j-axis at I-point
157      CALL iom_put( "snowpre"     , sprecip * 86400.    )        ! snow precipitation
158      CALL iom_put( "micesalt"    , smt_i   * zswi      )        ! mean ice salinity
159
160      CALL iom_put( "icestr"      , strength * zswi )    ! ice strength
161      CALL iom_put( "idive"       , divu_i              )    ! divergence
162      CALL iom_put( "ishear"      , shear_i             )    ! shear
163      CALL iom_put( "snowvol"     , vt_s   * zswi       )        ! snow volume
164     
165      CALL iom_put( "icetrp"      , diag_trp_vi * rday  )        ! ice volume transport
166      CALL iom_put( "snwtrp"      , diag_trp_vs * rday  )        ! snw volume transport
167      CALL iom_put( "saltrp"      , diag_trp_smv * rday * rhoic ) ! salt content transport
168      CALL iom_put( "deitrp"      , diag_trp_ei         )        ! advected ice enthalpy (W/m2)
169      CALL iom_put( "destrp"      , diag_trp_es         )        ! advected snw enthalpy (W/m2)
170
171      CALL iom_put( "sfxbog"      , sfx_bog * rday      )        ! salt flux from bottom growth
172      CALL iom_put( "sfxbom"      , sfx_bom * rday      )        ! salt flux from bottom melting
173      CALL iom_put( "sfxsum"      , sfx_sum * rday      )        ! salt flux from surface melting
174      CALL iom_put( "sfxsni"      , sfx_sni * rday      )        ! salt flux from snow ice formation
175      CALL iom_put( "sfxopw"      , sfx_opw * rday      )        ! salt flux from open water formation
176      CALL iom_put( "sfxdyn"      , sfx_dyn * rday      )        ! salt flux from ridging rafting
177      CALL iom_put( "sfxres"      , sfx_res * rday      )        ! salt flux from limupdate (resultant)
178      CALL iom_put( "sfxbri"      , sfx_bri * rday      )        ! salt flux from brines
179      CALL iom_put( "sfxsub"      , sfx_sub * rday      )        ! salt flux from sublimation
180      CALL iom_put( "sfx"         , sfx     * rday      )        ! total salt flux
181
182      ztmp = rday / rhoic
183      CALL iom_put( "vfxres"     , wfx_res * ztmp       )        ! daily prod./melting due to limupdate
184      CALL iom_put( "vfxopw"     , wfx_opw * ztmp       )        ! daily lateral thermodynamic ice production
185      CALL iom_put( "vfxsni"     , wfx_sni * ztmp       )        ! daily snowice ice production
186      CALL iom_put( "vfxbog"     , wfx_bog * ztmp       )        ! daily bottom thermodynamic ice production
187      CALL iom_put( "vfxdyn"     , wfx_dyn * ztmp       )        ! daily dynamic ice production (rid/raft)
188      CALL iom_put( "vfxsum"     , wfx_sum * ztmp       )        ! surface melt
189      CALL iom_put( "vfxbom"     , wfx_bom * ztmp       )        ! bottom melt
190      CALL iom_put( "vfxice"     , wfx_ice * ztmp       )        ! total ice growth/melt
191
192      IF ( iom_use( "vfxthin" ) ) THEN   ! ice production for open water + thin ice (<20cm) => comparable to observations 
193         WHERE( htm_i(:,:) < 0.2 .AND. htm_i(:,:) > 0. ) ; z2d = wfx_bog
194         ELSEWHERE                                       ; z2d = 0._wp
195         END WHERE
196         CALL iom_put( "vfxthin", ( wfx_opw + z2d ) * ztmp )
197      ENDIF
198
199      ztmp = rday / rhosn
200      CALL iom_put( "vfxspr"     , wfx_spr * ztmp       )        ! precip (snow)
201      CALL iom_put( "vfxsnw"     , wfx_snw * ztmp       )        ! total snw growth/melt
202      CALL iom_put( "vfxsub"     , wfx_sub * ztmp       )        ! sublimation (snow/ice)
203      CALL iom_put( "vfxsub_err" , wfx_err_sub * ztmp   )        ! "excess" of sublimation sent to ocean     
204 
205      CALL iom_put( "afxtot"     , afx_tot              )        ! concentration tendency (total)
206      CALL iom_put( "afxdyn"     , afx_dyn              )        ! concentration tendency (dynamics)
207      CALL iom_put( "afxthd"     , afx_thd              )        ! concentration tendency (thermo)
208
209      CALL iom_put ('hfxthd'     , hfx_thd(:,:)         )   
210      CALL iom_put ('hfxdyn'     , hfx_dyn(:,:)         )   
211      CALL iom_put ('hfxres'     , hfx_res(:,:)         )   
212      CALL iom_put ('hfxout'     , hfx_out(:,:)         )   
213      CALL iom_put ('hfxin'      , hfx_in(:,:)          )   
214      CALL iom_put ('hfxsnw'     , hfx_snw(:,:)         )   
215      CALL iom_put ('hfxsub'     , hfx_sub(:,:)         )   
216      CALL iom_put ('hfxerr'     , hfx_err(:,:)         )   
217      CALL iom_put ('hfxerr_rem' , hfx_err_rem(:,:)     )   
218     
219      CALL iom_put ('hfxsum'     , hfx_sum(:,:)         )   
220      CALL iom_put ('hfxbom'     , hfx_bom(:,:)         )   
221      CALL iom_put ('hfxbog'     , hfx_bog(:,:)         )   
222      CALL iom_put ('hfxdif'     , hfx_dif(:,:)         )   
223      CALL iom_put ('hfxopw'     , hfx_opw(:,:)         )   
224      CALL iom_put ('hfxtur'     , fhtur(:,:) * SUM( a_i_b(:,:,:), dim=3 ) ) ! turbulent heat flux at ice base
225      CALL iom_put ('hfxdhc'     , diag_heat(:,:)       )   ! Heat content variation in snow and ice
226      CALL iom_put ('hfxspr'     , hfx_spr(:,:)         )   ! Heat content of snow precip
227
228      !----------------------------------
229      ! Output category-dependent fields
230      !----------------------------------
231      IF ( iom_use( "iceconc_cat"  ) )  CALL iom_put( "iceconc_cat"      , a_i   * zswi2   )        ! area for categories
232      IF ( iom_use( "icethic_cat"  ) )  CALL iom_put( "icethic_cat"      , ht_i  * zswi2   )        ! thickness for categories
233      IF ( iom_use( "snowthic_cat" ) )  CALL iom_put( "snowthic_cat"     , ht_s  * zswi2   )        ! snow depth for categories
234      IF ( iom_use( "salinity_cat" ) )  CALL iom_put( "salinity_cat"     , sm_i  * zswi2   )        ! salinity for categories
235      ! ice temperature
236      IF ( iom_use( "icetemp_cat"  ) )  CALL iom_put( "icetemp_cat", ( SUM( t_i(:,:,:,:), dim=3 ) * r1_nlay_i - rt0 ) * zswi2 )
237      ! snow temperature
238      IF ( iom_use( "snwtemp_cat"  ) )  CALL iom_put( "snwtemp_cat", ( SUM( t_s(:,:,:,:), dim=3 ) * r1_nlay_s - rt0 ) * zswi2 )
239      ! ice age
240      IF ( iom_use( "iceage_cat"   ) )  CALL iom_put( "iceage_cat" , o_i * zswi2 * z1_365 )
241      ! brine volume
242      IF ( iom_use( "brinevol_cat" ) )  CALL iom_put( "brinevol_cat", bv_i * 100. * zswi2 )
243
244      !--------------------------------
245      ! Add-ons for SIMIP
246      !--------------------------------
247      zrho1 = ( rau0 - rhoic ) / rau0; zrho2 = rhosn / rau0
248
249      IF  ( iom_use( "icethic"  ) ) CALL iom_put( "icethic"     , htm_i * zswi               )          ! ice thickness
250      IF  ( iom_use( "icepres"  ) ) CALL iom_put( "icepres"     , zswi                       )          ! ice presence (1 or 0)
251      IF  ( iom_use( "snowthic" ) ) CALL iom_put( "snowthic"    , htm_s * zswi               )          ! snow thickness       
252      IF  ( iom_use( "icemass"  ) ) CALL iom_put( "icemass"     , rhoic * vt_i(:,:) * zswi   )          ! ice mass per cell area
253      IF  ( iom_use( "snomass"  ) ) CALL iom_put( "snomass"     , rhosn * vt_s(:,:) * zswi   )          ! snow mass per cell area
254      IF  ( iom_use( "icesnt"   ) ) CALL iom_put( "icesnt"      , ( tm_si - rt0 ) * zswi     )          ! snow-ice interface temperature
255      IF  ( iom_use( "icebot"   ) ) CALL iom_put( "icebot"      , ( t_bo  - rt0 ) * zswi     )          ! ice bottom temperature
256      IF  ( iom_use( "icesmass" ) ) CALL iom_put( "icesmass"    , SUM( smv_i, DIM = 3 ) * rhoic / 1000. * zswi )         ! mass of salt in sea ice per cell area
257      IF  ( iom_use( "icefb"    ) ) CALL iom_put( "icefb"       , ( zrho1 * htm_i(:,:) - zrho2 * htm_s(:,:) ) * zswi )   ! ice freeboard
258
259      IF  ( iom_use( "wfxsum"   ) ) CALL iom_put( "wfxsum"      ,   wfx_sum                  )          ! Freshwater flux from sea-ice surface
260      IF  ( iom_use( "dmithd"   ) ) CALL iom_put( "dmithd"      , - wfx_bog - wfx_bom - wfx_sum   &     ! Sea-ice mass change from thermodynamics
261              &                     - wfx_sni - wfx_opw - wfx_res )
262      IF  ( iom_use( "dmidyn"   ) ) CALL iom_put( "dmidyn"      ,   diag_dmi_dyn             )          ! Sea-ice mass change from dynamics
263      IF  ( iom_use( "dmiopw"   ) ) CALL iom_put( "dmiopw"      , - wfx_opw                  )          ! Sea-ice mass change through growth in open water
264      IF  ( iom_use( "dmibog"   ) ) CALL iom_put( "dmibog"      , - wfx_bog                  )          ! Sea-ice mass change through basal growth
265      IF  ( iom_use( "dmisni"   ) ) CALL iom_put( "dmisni"      , - wfx_sni                  )          ! Sea-ice mass change through snow-to-ice conversion
266      IF  ( iom_use( "dmisum"   ) ) CALL iom_put( "dmisum"      , - wfx_sum                  )          ! Sea-ice mass change through surface melting
267      IF  ( iom_use( "dmibom"   ) ) CALL iom_put( "dmibom"      , - wfx_bom                  )          ! Sea-ice mass change through bottom melting
268      IF  ( iom_use( "dmtsub"   ) ) CALL iom_put( "dmtsub"      , - wfx_sub                  )          ! Sea-ice mass change through evaporation and sublimation
269      IF  ( iom_use( "dmsspr"   ) ) CALL iom_put( "dmsspr"      , - wfx_spr                  )          ! snow mass change through snow fall
270      IF  ( iom_use( "dmsssi"   ) ) CALL iom_put( "dmsssi"      ,   wfx_sni*rhosn/rhoic      )          ! snow mass change through snow-to-ice conversion
271
272      IF  ( iom_use( "dmsmel"   ) ) CALL iom_put( "dmsmel"      , - wfx_snw_sum              )          ! snow mass change through melt
273      IF  ( iom_use( "dmsdyn"   ) ) CALL iom_put( "dmsdyn"      ,   diag_dms_dyn             )          ! snow mass change through dynamics
274
275      IF  ( iom_use( "hfxconbo" ) ) CALL iom_put( "hfxconbo"    ,   diag_fc_bo               )          ! bottom conduction flux
276      IF  ( iom_use( "hfxconsu" ) ) CALL iom_put( "hfxconsu"    ,   diag_fc_su               )          ! surface conduction flux
277
278      IF  ( iom_use( "wfxtot"   ) ) CALL iom_put( "wfxtot"      ,   -wfx_ice                 )          ! total freshwater flux from sea ice
279
280      IF  ( iom_use( "dmtxdyn"  ) ) CALL iom_put( "dmtxdyn"     ,   diag_dmtx_dyn            )          ! X-component of sea-ice mass transport
281      IF  ( iom_use( "dmtydyn"  ) ) CALL iom_put( "dmtydyn"     ,   diag_dmty_dyn            )          ! Y-component of sea-ice mass transport
282
283      IF  ( iom_use( "utau_oi"  ) ) CALL iom_put( "utau_oi"     ,   diag_utau_oi*zswi        )          ! X-component of ocean stress on sea ice
284      IF  ( iom_use( "vtau_oi"  ) ) CALL iom_put( "vtau_oi"     ,   diag_vtau_oi*zswi        )          ! Y-component of ocean stress on sea ice
285
286      IF  ( iom_use( "dssh_dx"  ) ) CALL iom_put( "dssh_dx"     ,   diag_dssh_dx*zswi        )          ! Sea-surface tilt term in force balance (x-component)
287      IF  ( iom_use( "dssh_dy"  ) ) CALL iom_put( "dssh_dy"     ,   diag_dssh_dy*zswi        )          ! Sea-surface tilt term in force balance (y-component)
288
289      IF  ( iom_use( "corstrx"  ) ) CALL iom_put( "corstrx"     ,   diag_corstrx*zswi        )          ! Coriolis force term in force balance (x-component)
290      IF  ( iom_use( "corstry"  ) ) CALL iom_put( "corstry"     ,   diag_corstry*zswi        )          ! Coriolis force term in force balance (y-component)
291
292      IF  ( iom_use( "intstrx"  ) ) CALL iom_put( "intstrx"     ,   diag_intstrx*zswi        )          ! Internal force term in force balance (x-component)
293      IF  ( iom_use( "intstry"  ) ) CALL iom_put( "intstry"     ,   diag_intstry*zswi        )          ! Internal force term in force balance (y-component)
294
295      IF  ( iom_use( "normstr"  ) ) CALL iom_put( "normstr"     ,   diag_sig1   *zswi        )          ! Normal stress
296      IF  ( iom_use( "sheastr"  ) ) CALL iom_put( "sheastr"     ,   diag_sig2   *zswi        )          ! Shear stress
297
298      !--------------------------------
299      ! Global ice diagnostics (SIMIP)
300      !--------------------------------
301
302      IF ( iom_use( "NH_icearea" ) .OR. iom_use( "NH_icevolu" ) .OR. iom_use( "NH_iceextt" ) )   THEN   ! NH integrated diagnostics
303 
304         WHERE( fcor > 0._wp ); zswi(:,:) = 1.0e-12
305         ELSEWHERE            ; zswi(:,:) = 0.
306         END WHERE
307
308         zdiag_area_nh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) )
309         zdiag_volu_nh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) )
310
311         WHERE( fcor > 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12
312         ELSEWHERE                              ; zswi(:,:) = 0.
313         END WHERE
314
315         zdiag_extt_nh = glob_sum( zswi(:,:) * e12t(:,:) )
316
317         IF ( iom_use( "NH_icearea" ) ) CALL iom_put( "NH_icearea" ,  zdiag_area_nh  )
318         IF ( iom_use( "NH_icevolu" ) ) CALL iom_put( "NH_icevolu" ,  zdiag_volu_nh  )
319         IF ( iom_use( "NH_iceextt" ) ) CALL iom_put( "NH_iceextt" ,  zdiag_extt_nh  )
320
321      ENDIF
322
323      IF ( iom_use( "SH_icearea" ) .OR. iom_use( "SH_icevolu" ) .OR. iom_use( "SH_iceextt" ) )   THEN   ! SH integrated diagnostics
324
325         WHERE( fcor < 0._wp ); zswi(:,:) = 1.0e-12; 
326         ELSEWHERE            ; zswi(:,:) = 0.
327         END WHERE
328
329         zdiag_area_sh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) ) 
330         zdiag_volu_sh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) )
331
332         WHERE( fcor < 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12
333         ELSEWHERE                              ; zswi(:,:) = 0.
334         END WHERE
335
336         zdiag_extt_sh = glob_sum( zswi(:,:) * e12t(:,:) )
337
338         IF ( iom_use( "SH_icearea" ) ) CALL iom_put( "SH_icearea", zdiag_area_sh )
339         IF ( iom_use( "SH_icevolu" ) ) CALL iom_put( "SH_icevolu", zdiag_volu_sh )
340         IF ( iom_use( "SH_iceextt" ) ) CALL iom_put( "SH_iceextt", zdiag_extt_sh )
341
342      ENDIF 
343
344      !--------------------------------
345      ! Fluxes through straits (SIMIP)
346      !--------------------------------
347      !
348      ! Valid only for ORCA-like grids
349      !
350      ! 4 Arctic passages are considered (Fram, CAA, Barents, Bering; see Notz et al (GMD 2016) for definitions)
351      !
352      ! Fram and Bering  straits are easy because they follow parallels
353      ! Barents and Canadian Arctic Archipelago are less easy because they do not, which is why they look so awful.
354      !
355
356      IF ( iom_use( "strait_arfl" ) .OR. iom_use( "strait_mifl" ) .OR. iom_use( "strait_msfl" ) .AND. cp_cfg == "orca" ) THEN
357
358         zdiag_area_strait(:) = 0._wp   ;   zdiag_mice_strait(:) = 0._wp   ;   zdiag_msno_strait(:) = 0._wp
359   
360         !------------------------------
361         ! === Fram & Bering Straits ===
362         !------------------------------
363         
364         SELECT CASE ( jp_cfg ) 
365         
366         CASE ( 2 )   ! --- ORCA2
367         
368            ! Fram Strait   (i_strait = 1)
369            ji0(1) = 133   ;   ji1(1) = 136
370            jj0(1) = 136 
371
372            ! Bering Strait (i_strait = 4)
373            ji0(4) = 55    ;   ji1(4) = 56
374            jj0(4) = 122
375         
376         CASE ( 1 )   ! --- eORCA1
377         
378            ! Fram Strait
379            ji0(1) = 268   ;   ji1(1) = 277
380            jj0(1) = 311
381
382            ! Bering Strait
383            ji0(4) = 113   ;   jj1(4) = 115
384            jj0(4) = 285
385           
386         END SELECT
387         
388         DO i_strait = 1, 4, 3
389
390            DO ji = mi0(ji0), mi1(ji1)
391               jj = mj0(jj0)
392   
393               zdiag_area_strait(i_strait) = zdiag_area_strait(i_strait)                                 &     ! --- ice area flux ---
394                   &                + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         &     ! northwards (positive) flow
395                   &                + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 )               ! southwards (negative) flow
396     
397               zdiag_mice_strait(i_strait) = zdiag_mice_strait(i_strait) + rhoic *                       &     ! --- ice mass flux  ---
398                   &                ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         & 
399                   &                + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 ) )         
400     
401               zdiag_msno_strait(i_strait) = zdiag_msno_strait(i_strait) + rhosn *                       &     ! --- snow mass flux ---
402                   &                ( vt_s(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         & 
403                   &                + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 ) ) 
404   
405            END DO
406
407         END DO
408
409         !---------------------
410         ! === Barents opening
411         !---------------------
412   
413         SELECT CASE ( jp_cfg ) 
414
415            CASE ( 1 )   ! 'eORCA1'
416
417               Nu = 11   ! U-Flow
418               zui(1:Nu) = (/ 282,283,284,285,286,286,287,288,289,290,292/)
419               zuj(1:Nu) = (/ 308,307,306,305,304,303,302,301,300,299,298/)
420
421               Nv = 9    ! V-Flow
422               zvi(1:Nv) = (/ 282,283,284,285,286,287,288,289,290/)
423               zvj(1:Nv) = (/ 308,307,306,305,303,302,301,300,299/)
424
425            CASE ( 2 )   ! 'ORCA2'
426
427               Nu = 5    ! U-Flow
428               zui(1:Nu) = (/ 141,142,142,143,144 /)
429               zuj(1:Nu) = (/ 134,133,132,131,130 /)
430
431               Nv = 4    ! V-Flow
432               zvi(1:Nv) = (/ 140,141,142,143 /)
433               zvj(1:Nv) = (/ 135,134,132,131 /)
434
435         END SELECT
436
437         ! Barents U-flow
438         zfarea_u = 0._wp   ;   zfmice_u = 0._wp   ;   zfmsno_u = 0._wp
439   
440         DO ii = 1, Nu
441
442            ji = mi0(zui(ii))
443            jj = mj0(zuj(ii))
444
445            zfarea_u          = zfarea_u                                           & ! --- ice area zonal flux ---
446                &             + at_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         & ! --- northward
447                &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 )           ! --- southward
448            zfmice_u          = zfmice_u + rhoic *                                 & ! --- ice mass zonal flux ---
449                &             ( vt_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &   
450                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )         
451            zfmsno_u          = zfmsno_u + rhosn *                                 & ! --- snow mass zonal flux ---
452                &             ( vt_s(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &   
453                &             + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )         
454         END DO
455   
456         ! Barents V-flow
457         zfarea_v = 0._wp   ;   zfmice_v = 0._wp   ;   zfmsno_v = 0._wp
458
459         DO ii  = 1, Nv 
460
461            ji = mi0(zvi(ii))
462            jj = mj0(zvj(ii))
463
464            zfarea_v          = zfarea_v                                           & ! --- ice area meridian flux ---
465                &             + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! --- eastward
466                &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 )           ! --- westward
467            zfmice_v          = zfmice_v + rhoic *                                 & ! --- ice mass meridian flux ---
468                &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & !
469                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         !
470            zfmsno_v          = zfmsno_v + rhosn *                                 & ! --- snow mass meridian flux ---
471                &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & !
472                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         !
473         END DO
474   
475         ! Sum Barents U-/V- contributions
476         zdiag_area_strait(3) = zfarea_u + zfarea_v 
477         zdiag_mice_strait(3) = zfmice_u + zfmice_v
478         zdiag_msno_strait(3) = zfmsno_u + zfmsno_v
479
480         !---------------------
481         ! === CAA throughflow
482         !---------------------
483   
484         SELECT CASE ( jp_cfg ) 
485
486            CASE ( 1 )   ! eORCA1
487
488               ! V-flow through Nares Strait
489               Nv = 4
490               zvi(1:Nv) = (/ 254,255,256,257 /)
491               zvj(1:Nv) = (/ 317,317,317,317 /)
492
493               ! U-flow through Queen Elisabeth Islands and McClure straits
494               Nu = 8 
495               zui(1:Nu) = (/ 231,231,231,  132,132,132,132,132  /)
496               zuj(1:Nu) = (/ 328,329,330,  318,319,320,321,322  /)
497
498               zfarea_u = 0._wp   ;   zfmice_u = 0._wp   ;   zfmsno_u = 0._wp
499
500               DO ii = 1, Nu
501
502                  ji = mi0(zui(ii))
503                  jj = mj0(zuj(ii))
504
505                  zfarea_u          = zfarea_u                                                           & ! --- ice area zonal flux ---
506                      &             + at_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         & ! --- eastward
507                      &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 )           ! --- westward
508                  zfmice_u          = zfmice_u + rhoic *                                                 & ! --- ice mass zonal flux ---
509                      &             ( vt_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &   
510                      &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )         
511                  zfmsno_u          = zfmsno_u + rhosn *                                                 & ! --- snow mass zonal flux ---
512                      &             ( vt_s(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &   
513                      &             + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )         
514
515               END DO
516
517
518            CASE ( 2 )   ! ORCA2
519
520               ! V-flow through Nares Strait
521               Nv = 2
522               zvi(1:Nv) = (/ 117,118 /)
523               zvj(1:Nv) = (/ 145,145 /)
524
525               ! U-flow through Queen Elisabeth Islands and McClure straits (not resolved in ORCA2)
526               zfarea_u = 0._wp   ;   zfmice_u = 0._wp   ;   zfmsno_u = 0._wp
527
528            END SELECT
529   
530         ! V-flow through Nares Strait
531         zfarea_v = 0._wp   ;   zfmice_v = 0._wp   ;   zfmsno_v = 0._wp
532   
533         DO ii = 1, Nv 
534
535            ji = mi0(zvi(ii))
536            jj = mj0(zvj(ii))
537
538            zfarea_v          = zfarea_v                                           & ! --- ice area meridian flux ---
539                &             + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! --- eastward
540                &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 )           ! --- westward
541            zfmice_v          = zfmice_v + rhoic *                                 & ! --- ice mass meridian flux ---
542                &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & !
543                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         !
544            zfmsno_v          = zfmsno_v + rhosn *                                 & ! --- snow mass meridian flux ---
545                &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & !
546                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         !
547
548         END DO
549
550         ! Sum U/V contributions
551         zdiag_area_strait(2) = zfarea_u + zfarea_v 
552         zdiag_mice_strait(2) = zfmice_u + zfmice_v
553         zdiag_msno_strait(2) = zfmsno_u + zfmsno_v
554   
555         ! === Ncdf output
556         IF ( iom_use("strait_arfl") ) CALL iom_put( "strait_arfl", zdiag_area_strait )
557         IF ( iom_use("strait_mifl") ) CALL iom_put( "strait_mifl", zdiag_mice_strait )
558         IF ( iom_use("strait_msfl") ) CALL iom_put( "strait_msfl", zdiag_msno_strait ) 
559
560         WRITE(numout,*) " area flx ", zdiag_area_strait(:)
561         WRITE(numout,*) " mice flx ", zdiag_mice_strait(:)
562         WRITE(numout,*) " msno flx ", zdiag_msno_strait(:)
563
564      ENDIF
565
566      !     !  Create an output files (output.lim.abort.nc) if S < 0 or u > 20 m/s
567      !     IF( kindic < 0 )   CALL lim_wri_state( 'output.abort' )
568      !     not yet implemented
569     
570      CALL wrk_dealloc( jpi, jpj, jpl, zswi2 )
571      CALL wrk_dealloc( jpi, jpj     , z2d, zswi )
572
573      IF( nn_timing == 1 )  CALL timing_stop('limwri')
574     
575   END SUBROUTINE lim_wri
576
577 
578   SUBROUTINE lim_wri_state( kt, kid, kh_i )
579      !!---------------------------------------------------------------------
580      !!                 ***  ROUTINE lim_wri_state  ***
581      !!       
582      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
583      !!      the instantaneous ice state and forcing fields for ice model
584      !!        Used to find errors in the initial state or save the last
585      !!      ocean state in case of abnormal end of a simulation
586      !!
587      !! History :
588      !!   4.0  !  2013-06  (C. Rousset)
589      !!----------------------------------------------------------------------
590      INTEGER, INTENT( in )   ::   kt               ! ocean time-step index)
591      INTEGER, INTENT( in )   ::   kid , kh_i
592      INTEGER                 ::   nz_i, jl
593      REAL(wp), DIMENSION(jpl) :: jcat
594      !!----------------------------------------------------------------------
595      DO jl = 1, jpl
596         jcat(jl) = REAL(jl)
597      ENDDO
598     
599      CALL histvert( kid, "ncatice", "Ice Categories","", jpl, jcat, nz_i, "up")
600
601      CALL histdef( kid, "sithic", "Ice thickness"           , "m"      ,   &
602      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
603      CALL histdef( kid, "siconc", "Ice concentration"       , "%"      ,   &
604      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
605      CALL histdef( kid, "sitemp", "Ice temperature"         , "C"      ,   &
606      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
607      CALL histdef( kid, "sivelu", "i-Ice speed "            , "m/s"    ,   &
608      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
609      CALL histdef( kid, "sivelv", "j-Ice speed "            , "m/s"    ,   &
610      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
611      CALL histdef( kid, "sistru", "i-Wind stress over ice " , "Pa"     ,   &
612      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
613      CALL histdef( kid, "sistrv", "j-Wind stress over ice " , "Pa"     ,   &
614      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
615      CALL histdef( kid, "sisflx", "Solar flux over ocean"     , "w/m2" ,   &
616      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
617      CALL histdef( kid, "sinflx", "Non-solar flux over ocean" , "w/m2" ,   &
618      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
619      CALL histdef( kid, "isnowpre", "Snow precipitation"      , "kg/m2/s",   &
620      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
621      CALL histdef( kid, "sisali", "Ice salinity"            , "PSU"    ,   &
622      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
623      CALL histdef( kid, "sivolu", "Ice volume"              , "m"      ,   &
624      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
625      CALL histdef( kid, "sidive", "Ice divergence"          , "10-8s-1",   &
626      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
627
628      CALL histdef( kid, "vfxbog", "Ice bottom production"   , "m/s"    ,   &
629      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
630      CALL histdef( kid, "vfxdyn", "Ice dynamic production"  , "m/s"    ,   &
631      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
632      CALL histdef( kid, "vfxopw", "Ice open water prod"     , "m/s"    ,   &
633      &       jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
634      CALL histdef( kid, "vfxsni", "Snow ice production "    , "m/s"    ,   &
635      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
636      CALL histdef( kid, "vfxres", "Ice prod from limupdate" , "m/s"    ,   &
637      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
638      CALL histdef( kid, "vfxbom", "Ice bottom melt"         , "m/s"    ,   &
639      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
640      CALL histdef( kid, "vfxsum", "Ice surface melt"        , "m/s"    ,   &
641      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
642
643      CALL histdef( kid, "sithicat", "Ice thickness"         , "m"      ,   &
644      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
645      CALL histdef( kid, "siconcat", "Ice concentration"     , "%"      ,   &
646      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
647      CALL histdef( kid, "sisalcat", "Ice salinity"           , ""      ,   &
648      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
649      CALL histdef( kid, "sitemcat", "Ice temperature"       , "C"      ,   &
650      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
651      CALL histdef( kid, "snthicat", "Snw thickness"         , "m"      ,   &
652      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
653      CALL histdef( kid, "sntemcat", "Snw temperature"       , "C"      ,   &
654      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
655
656      CALL histend( kid, snc4set )   ! end of the file definition
657
658      CALL histwrite( kid, "sithic", kt, htm_i         , jpi*jpj, (/1/) )   
659      CALL histwrite( kid, "siconc", kt, at_i          , jpi*jpj, (/1/) )
660      CALL histwrite( kid, "sitemp", kt, tm_i - rt0    , jpi*jpj, (/1/) )
661      CALL histwrite( kid, "sivelu", kt, u_ice          , jpi*jpj, (/1/) )
662      CALL histwrite( kid, "sivelv", kt, v_ice          , jpi*jpj, (/1/) )
663      CALL histwrite( kid, "sistru", kt, utau_ice       , jpi*jpj, (/1/) )
664      CALL histwrite( kid, "sistrv", kt, vtau_ice       , jpi*jpj, (/1/) )
665      CALL histwrite( kid, "sisflx", kt, qsr , jpi*jpj, (/1/) )
666      CALL histwrite( kid, "sinflx", kt, qns , jpi*jpj, (/1/) )
667      CALL histwrite( kid, "isnowpre", kt, sprecip        , jpi*jpj, (/1/) )
668      CALL histwrite( kid, "sisali", kt, smt_i          , jpi*jpj, (/1/) )
669      CALL histwrite( kid, "sivolu", kt, vt_i           , jpi*jpj, (/1/) )
670      CALL histwrite( kid, "sidive", kt, divu_i*1.0e8   , jpi*jpj, (/1/) )
671
672      CALL histwrite( kid, "vfxbog", kt, wfx_bog        , jpi*jpj, (/1/) )
673      CALL histwrite( kid, "vfxdyn", kt, wfx_dyn        , jpi*jpj, (/1/) )
674      CALL histwrite( kid, "vfxopw", kt, wfx_opw        , jpi*jpj, (/1/) )
675      CALL histwrite( kid, "vfxsni", kt, wfx_sni        , jpi*jpj, (/1/) )
676      CALL histwrite( kid, "vfxres", kt, wfx_res        , jpi*jpj, (/1/) )
677      CALL histwrite( kid, "vfxbom", kt, wfx_bom        , jpi*jpj, (/1/) )
678      CALL histwrite( kid, "vfxsum", kt, wfx_sum        , jpi*jpj, (/1/) )
679
680      CALL histwrite( kid, "sithicat", kt, ht_i        , jpi*jpj*jpl, (/1/) )   
681      CALL histwrite( kid, "siconcat", kt, a_i         , jpi*jpj*jpl, (/1/) )   
682      CALL histwrite( kid, "sisalcat", kt, sm_i        , jpi*jpj*jpl, (/1/) )   
683      CALL histwrite( kid, "sitemcat", kt, tm_i - rt0  , jpi*jpj*jpl, (/1/) )   
684      CALL histwrite( kid, "snthicat", kt, ht_s        , jpi*jpj*jpl, (/1/) )   
685      CALL histwrite( kid, "sntemcat", kt, tm_su - rt0 , jpi*jpj*jpl, (/1/) )   
686
687      ! Close the file
688      ! -----------------
689      !CALL histclo( kid )
690
691    END SUBROUTINE lim_wri_state
692
693#else
694   !!----------------------------------------------------------------------
695   !!   Default option :         Empty module          NO LIM sea-ice model
696   !!----------------------------------------------------------------------
697CONTAINS
698   SUBROUTINE lim_wri          ! Empty routine
699   END SUBROUTINE lim_wri
700#endif
701
702   !!======================================================================
703END MODULE limwri
Note: See TracBrowser for help on using the repository browser.