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 @ 8258

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

fix temporal averages for SIMIP
  • Property svn:keywords set to Id
File size: 32.6 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) ::  z2da, z2db, ztmp, zrho1, zrho2, zmiss_val
57      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zswi2
58      REAL(wp), POINTER, DIMENSION(:,:)   ::  z2d, zswi    ! 2D workspace
59      REAL(wp), POINTER, DIMENSION(:,:)   ::  zfb          ! ice freeboard
60      REAL(wp), POINTER, DIMENSION(:,:)   ::  zamask, zamask15 ! 15% concentration mask
61      REAL(wp), POINTER, DIMENSION(:,:)   ::  zmiss        ! missing value array
62
63      ! Global ice diagnostics (SIMIP)
64      REAL(wp) ::  zdiag_area_nh, &   ! area, extent, volume
65         &         zdiag_extt_nh, &
66         &         zdiag_area_sh, & 
67         &         zdiag_extt_sh, & 
68         &         zdiag_volu_nh, & 
69         &         zdiag_volu_sh 
70
71      !!-------------------------------------------------------------------
72
73      IF( nn_timing == 1 )  CALL timing_start('limwri')
74
75      CALL wrk_alloc( jpi, jpj, jpl, zswi2 )
76      CALL wrk_alloc( jpi, jpj     , z2d, zswi )
77      CALL wrk_alloc( jpi, jpj     , zfb, zamask, zamask15 )
78      CALL wrk_alloc( jpi, jpj     , zmiss                 )
79
80      !-----------------------------
81      ! Missing value array
82      !-----------------------------
83      zmiss_val  = -1.0e20
84      zmiss(:,:) = zmiss_val * ( 1. - zswi(:,:) )
85
86      !-----------------------------
87      ! Mean category values
88      !-----------------------------
89
90      ! brine volume
91      CALL lim_var_bv 
92
93      ! tresholds for outputs
94      DO jj = 1, jpj
95         DO ji = 1, jpi
96            zswi(ji,jj)  = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice, 0 if no ice
97         END DO
98      END DO
99      DO jj = 1, jpj
100         DO ji = 1, jpi
101            zamask(ji,jj)  = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - 0.05 ) )   ! 1 if 5% ice, 0 if less - required to mask thickness and snow depth
102            zamask15(ji,jj)  = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - 0.15 ) ) ! 1 if 15% ice, 0 if less
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') ) CALL iom_put( "uice_ipa"     , u_ice      )       ! ice velocity u component
130      IF ( iom_use('vice_ipa') ) CALL iom_put( "vice_ipa"     , v_ice      )       ! ice velocity v component
131
132      IF ( ( iom_use( "icevel" ) .OR. ( iom_use( "icevel_mv" ) ) ) THEN
133         DO jj = 2 , jpjm1
134            DO ji = 2 , jpim1
135               z2da  = ( u_ice(ji,jj) * umask(ji,jj,1) + u_ice(ji-1,jj) * umask(ji-1,jj,1) ) * 0.5_wp
136               z2db  = ( v_ice(ji,jj) * vmask(ji,jj,1) + v_ice(ji,jj-1) * vmask(ji,jj-1,1) ) * 0.5_wp
137               z2d(ji,jj) = SQRT( z2da * z2da + z2db * z2db )
138           END DO
139         END DO
140         CALL lbc_lnk( z2d, 'T', 1. )
141         IF ( iom_use( "icevel" )    CALL iom_put( "icevel"       , z2d        )                                            ! ice velocity module
142         IF ( iom_use( "icevel_mv" ) CALL iom_put( "icevel_mv"    , z2d(:,:) * zswi(:,:) + zmiss(:,:)                       ! ice velocity module (missing value)
143      ENDIF
144      !
145      IF ( iom_use( "miceage" ) )       CALL iom_put( "miceage"     , om_i * zswi * zamask15 )  ! mean ice age
146      IF ( iom_use( "micet" ) )         CALL iom_put( "micet"       , ( tm_i  - rt0 ) * zswi )  ! ice mean    temperature
147      IF ( iom_use( "icest" ) )         CALL iom_put( "icest"       , ( tm_su - rt0 ) * zswi )  ! ice surface temperature
148      IF ( iom_use( "icecolf" ) )       CALL iom_put( "icecolf"     , hicol                  )  ! frazil ice collection thickness
149      !
150      CALL iom_put( "isst"        , sst_m               )        ! sea surface temperature
151      CALL iom_put( "isss"        , sss_m               )        ! sea surface salinity
152      CALL iom_put( "iceconc"     , at_i  * zswi        )        ! ice concentration
153      CALL iom_put( "icevolu"     , vt_i  * zswi        )        ! ice volume = mean ice thickness over the cell
154      CALL iom_put( "icehc"       , et_i  * zswi        )        ! ice total heat content
155      CALL iom_put( "isnowhc"     , et_s  * zswi        )        ! snow total heat content
156      CALL iom_put( "ibrinv"      , bvm_i * zswi * 100. )        ! brine volume
157      CALL iom_put( "utau_ice"    , utau_ice*zswi       )        ! wind stress over ice along i-axis at I-point
158      CALL iom_put( "vtau_ice"    , vtau_ice*zswi       )        ! wind stress over ice along j-axis at I-point
159      CALL iom_put( "snowpre"     , sprecip * 86400.    )        ! snow precipitation
160      CALL iom_put( "micesalt"    , smt_i   * zswi      )        ! mean ice salinity
161
162      CALL iom_put( "icestr"      , strength * zswi )            ! ice strength
163      CALL iom_put( "idive"       , divu_i              )        ! divergence
164      CALL iom_put( "ishear"      , shear_i             )        ! shear
165      CALL iom_put( "snowvol"     , vt_s   * zswi       )        ! snow volume
166     
167      CALL iom_put( "icetrp"      , diag_trp_vi * rday  )        ! ice volume transport
168      CALL iom_put( "snwtrp"      , diag_trp_vs * rday  )        ! snw volume transport
169      CALL iom_put( "saltrp"      , diag_trp_smv * rday * rhoic ) ! salt content transport
170      CALL iom_put( "deitrp"      , diag_trp_ei         )        ! advected ice enthalpy (W/m2)
171      CALL iom_put( "destrp"      , diag_trp_es         )        ! advected snw enthalpy (W/m2)
172
173      CALL iom_put( "sfxbog"      , sfx_bog * rday      )        ! salt flux from bottom growth
174      CALL iom_put( "sfxbom"      , sfx_bom * rday      )        ! salt flux from bottom melting
175      CALL iom_put( "sfxsum"      , sfx_sum * rday      )        ! salt flux from surface melting
176      CALL iom_put( "sfxsni"      , sfx_sni * rday      )        ! salt flux from snow ice formation
177      CALL iom_put( "sfxopw"      , sfx_opw * rday      )        ! salt flux from open water formation
178      CALL iom_put( "sfxdyn"      , sfx_dyn * rday      )        ! salt flux from ridging rafting
179      CALL iom_put( "sfxres"      , sfx_res * rday      )        ! salt flux from limupdate (resultant)
180      CALL iom_put( "sfxbri"      , sfx_bri * rday      )        ! salt flux from brines
181      CALL iom_put( "sfxsub"      , sfx_sub * rday      )        ! salt flux from sublimation
182      CALL iom_put( "sfx"         , sfx     * rday      )        ! total salt flux
183
184      ztmp = rday / rhoic
185      CALL iom_put( "vfxres"     , wfx_res * ztmp       )        ! daily prod./melting due to limupdate
186      CALL iom_put( "vfxopw"     , wfx_opw * ztmp       )        ! daily lateral thermodynamic ice production
187      CALL iom_put( "vfxsni"     , wfx_sni * ztmp       )        ! daily snowice ice production
188      CALL iom_put( "vfxbog"     , wfx_bog * ztmp       )        ! daily bottom thermodynamic ice production
189      CALL iom_put( "vfxdyn"     , wfx_dyn * ztmp       )        ! daily dynamic ice production (rid/raft)
190      CALL iom_put( "vfxsum"     , wfx_sum * ztmp       )        ! surface melt
191      CALL iom_put( "vfxbom"     , wfx_bom * ztmp       )        ! bottom melt
192      CALL iom_put( "vfxice"     , wfx_ice * ztmp       )        ! total ice growth/melt
193
194      IF ( iom_use( "vfxthin" ) ) THEN   ! ice production for open water + thin ice (<20cm) => comparable to observations 
195         WHERE( htm_i(:,:) < 0.2 .AND. htm_i(:,:) > 0. ) ; z2d = wfx_bog
196         ELSEWHERE                                       ; z2d = 0._wp
197         END WHERE
198         CALL iom_put( "vfxthin", ( wfx_opw + z2d ) * ztmp )
199      ENDIF
200
201      ztmp = rday / rhosn
202      CALL iom_put( "vfxspr"     , wfx_spr * ztmp       )        ! precip (snow)
203      CALL iom_put( "vfxsnw"     , wfx_snw * ztmp       )        ! total snw growth/melt
204      CALL iom_put( "vfxsub"     , wfx_sub * ztmp       )        ! sublimation (snow/ice)
205      CALL iom_put( "vfxsub_err" , wfx_err_sub * ztmp   )        ! "excess" of sublimation sent to ocean     
206 
207      CALL iom_put( "afxtot"     , afx_tot              )        ! concentration tendency (total)
208      CALL iom_put( "afxdyn"     , afx_dyn              )        ! concentration tendency (dynamics)
209      CALL iom_put( "afxthd"     , afx_thd              )        ! concentration tendency (thermo)
210
211      CALL iom_put ('hfxthd'     , hfx_thd(:,:)         )   ! 
212      CALL iom_put ('hfxdyn'     , hfx_dyn(:,:)         )   ! 
213      CALL iom_put ('hfxres'     , hfx_res(:,:)         )   ! 
214      CALL iom_put ('hfxout'     , hfx_out(:,:)         )   ! 
215      CALL iom_put ('hfxin'      , hfx_in(:,:)          )   ! 
216      CALL iom_put ('hfxsnw'     , hfx_snw(:,:)         )   ! 
217      CALL iom_put ('hfxsub'     , hfx_sub(:,:)         )   ! 
218      CALL iom_put ('hfxerr'     , hfx_err(:,:)         )   ! 
219      CALL iom_put ('hfxerr_rem' , hfx_err_rem(:,:)     )   ! 
220     
221      CALL iom_put ('hfxsum'     , hfx_sum(:,:)         )   ! 
222      CALL iom_put ('hfxbom'     , hfx_bom(:,:)         )   ! 
223      CALL iom_put ('hfxbog'     , hfx_bog(:,:)         )   ! 
224      CALL iom_put ('hfxdif'     , hfx_dif(:,:)         )   ! 
225      CALL iom_put ('hfxopw'     , hfx_opw(:,:)         )   ! 
226      CALL iom_put ('hfxtur'     , fhtur(:,:) * SUM( a_i_b(:,:,:), dim=3 ) ) ! turbulent heat flux at ice base
227      CALL iom_put ('hfxdhc'     , diag_heat(:,:)       )   ! Heat content variation in snow and ice
228      CALL iom_put ('hfxspr'     , hfx_spr(:,:)         )   ! Heat content of snow precip
229
230      !----------------------------------
231      ! Output category-dependent fields
232      !----------------------------------
233      IF ( iom_use( "iceconc_cat"  ) )  CALL iom_put( "iceconc_cat"      , a_i   * zswi2   )        ! area for categories
234      IF ( iom_use( "icethic_cat"  ) )  CALL iom_put( "icethic_cat"      , ht_i  * zswi2   )        ! thickness for categories
235      IF ( iom_use( "snowthic_cat" ) )  CALL iom_put( "snowthic_cat"     , ht_s  * zswi2   )        ! snow depth for categories
236      IF ( iom_use( "salinity_cat" ) )  CALL iom_put( "salinity_cat"     , sm_i  * zswi2   )        ! salinity for categories
237      ! ice temperature
238      IF ( iom_use( "icetemp_cat"  ) )  CALL iom_put( "icetemp_cat", ( SUM( t_i(:,:,:,:), dim=3 ) * r1_nlay_i - rt0 ) * zswi2 )
239      ! snow temperature
240      IF ( iom_use( "snwtemp_cat"  ) )  CALL iom_put( "snwtemp_cat", ( SUM( t_s(:,:,:,:), dim=3 ) * r1_nlay_s - rt0 ) * zswi2 )
241      ! ice age
242      IF ( iom_use( "iceage_cat"   ) )  CALL iom_put( "iceage_cat" , o_i * zswi2 ) 
243      ! brine volume
244      IF ( iom_use( "brinevol_cat" ) )  CALL iom_put( "brinevol_cat", bv_i * 100. * zswi2 )
245
246      !--------------------------------
247      ! Add-ons for SIMIP
248      !--------------------------------
249      zrho1 = ( rau0 - rhoic ) / rau0; zrho2 = rhosn / rau0
250
251      IF  ( iom_use( "icepres"  ) ) CALL iom_put( "icepres"     , zswi(:,:)                     )                         ! Ice presence (1 or 0)
252      IF  ( iom_use( "icemass"  ) ) CALL iom_put( "icemass"     , rhoic * vt_i(:,:) * zswi(:,:) )                         ! Ice mass per cell area
253      IF  ( iom_use( "icethic"  ) ) CALL iom_put( "icethic"     , htm_i(:,:) * zamask(:,:)  * zswi(:,:) + zmiss(:,:) )    ! Ice thickness
254      IF  ( iom_use( "snomass"  ) ) CALL iom_put( "snomass"     , rhosn * vt_s(:,:)         * zswi(:,:) + zmiss(:,:) )    ! Snow mass per cell area
255      IF  ( iom_use( "snowthic" ) ) CALL iom_put( "snowthic"    , htm_s(:,:) * zamask(:,:)  * zswi(:,:) + zmiss(:,:) )    ! Snow thickness       
256      IF  ( iom_use( "icesnt"   ) ) CALL iom_put( "icesnt"      , ( tm_si(:,:) - rt0 )      * zswi(:,:) + zmiss(:,:) )    ! Snow-ice interface temperature
257      IF  ( iom_use( "icebot"   ) ) CALL iom_put( "icebot"      , ( t_bo(:,:)  - rt0 )      * zswi(:,:) + zmiss(:,:) )    ! Ice bottom temperature
258      IF  ( iom_use( "iceage_mv") ) CALL iom_put( "iceage_mv"   , om_i(:,:) * zamask15(:,:) * zswi(:,:) + zmiss(:,:) )    ! Ice age
259      IF  ( iom_use( "icesmass" ) ) CALL iom_put( "icesmass"    , SUM( smv_i, DIM = 3 ) * rhoic * 1.0e-3 * zswi(:,:) )    ! Mass of salt in sea ice per cell area
260      IF  ( iom_use( "icesal_mv") ) CALL iom_put( "icesal_mv"   , smt_i(:,:)                * zswi(:,:) + zmiss(:,:) )    ! Ice salinity
261
262      IF  ( iom_use( "icefb"    ) ) THEN
263         zfb(:,:) = ( zrho1 * htm_i(:,:) - zrho2 * htm_s(:,:) )                                         
264         WHERE( zfb < 0._wp ) ;   zfb = 0._wp ;   END WHERE
265                                    CALL iom_put( "icefb"       , zfb(:,:)                  * zswi(:,:) + zmiss(:,:) )    ! Ice freeboard
266      ENDIF
267
268      IF  ( iom_use( "snhc_mv"  ) ) CALL iom_put( "snhc_mv"     , et_s(:,:)                 * zswi(:,:) + zmiss(:,:) )    ! Snow total heat content
269
270      IF  ( iom_use( "dmithd"   ) ) CALL iom_put( "dmithd"      , - wfx_bog - wfx_bom - wfx_sum   &                       ! Sea-ice mass change from thermodynamics
271              &                     - wfx_sni - wfx_opw - wfx_res )
272      IF  ( iom_use( "dmidyn"   ) ) CALL iom_put( "dmidyn"      ,   diag_dmi_dyn             )                            ! Sea-ice mass change from dynamics
273      IF  ( iom_use( "dmiopw"   ) ) CALL iom_put( "dmiopw"      , - wfx_opw                  )                            ! Sea-ice mass change through growth in open water
274      IF  ( iom_use( "dmibog"   ) ) CALL iom_put( "dmibog"      , - wfx_bog                  )                            ! Sea-ice mass change through basal growth
275      IF  ( iom_use( "dmisni"   ) ) CALL iom_put( "dmisni"      , - wfx_sni                  )                            ! Sea-ice mass change through snow-to-ice conversion
276      IF  ( iom_use( "dmisum"   ) ) CALL iom_put( "dmisum"      , - wfx_sum                  )                            ! Sea-ice mass change through surface melting
277      IF  ( iom_use( "dmibom"   ) ) CALL iom_put( "dmibom"      , - wfx_bom                  )                            ! Sea-ice mass change through bottom melting
278
279      IF  ( iom_use( "dmtsub"   ) ) CALL iom_put( "dmtsub"      , - wfx_sub                  )                            ! Sea-ice mass change through evaporation and sublimation
280      IF  ( iom_use( "dmssub"   ) ) CALL iom_put( "dmssub"      , - wfx_snw_sub              )                            ! Snow mass change through sublimation
281      IF  ( iom_use( "dmisub"   ) ) CALL iom_put( "dmisub"      , - wfx_ice_sub              )                            ! Sea-ice mass change through sublimation
282
283      IF  ( iom_use( "dmsspr"   ) ) CALL iom_put( "dmsspr"      , - wfx_spr                  )                            ! Snow mass change through snow fall
284      IF  ( iom_use( "dmsssi"   ) ) CALL iom_put( "dmsssi"      ,   wfx_sni*rhosn/rhoic      )                            ! Snow mass change through snow-to-ice conversion
285
286      IF  ( iom_use( "dmsmel"   ) ) CALL iom_put( "dmsmel"      , - wfx_snw_sum              )                            ! Snow mass change through melt
287      IF  ( iom_use( "dmsdyn"   ) ) CALL iom_put( "dmsdyn"      ,   diag_dms_dyn             )                            ! Snow mass change through dynamics
288
289      IF  ( iom_use( "hfxconsu" ) ) CALL iom_put( "hfxsenso"    ,   -fhtur(:,:)              * zswi(:,:) + zmiss(:,:) )   ! Sensible oceanic heat flux
290      IF  ( iom_use( "hfxconbo" ) ) CALL iom_put( "hfxconbo"    ,   diag_fc_bo               * zswi(:,:) + zmiss(:,:) )   ! Bottom conduction flux
291      IF  ( iom_use( "hfxconsu" ) ) CALL iom_put( "hfxconsu"    ,   diag_fc_su               * zswi(:,:) + zmiss(:,:) )   ! Surface conduction flux
292
293      IF  ( iom_use( "wfxtot"   ) ) CALL iom_put( "wfxtot"      ,   wfx_ice(:,:)             * zswi(:,:) + zmiss(:,:) )   ! Total freshwater flux from sea ice
294      IF  ( iom_use( "wfxsum"   ) ) CALL iom_put( "wfxsum"      ,   wfx_sum(:,:)             * zswi(:,:) + zmiss(:,:) )   ! Freshwater flux from sea-ice surface
295      IF  ( iom_use( "sfx_mv"   ) ) CALL iom_put( "sfx_mv"      ,   sfx(:,:)                 * zswi(:,:) + zmiss(:,:) )   ! Total salt flux
296
297      IF  ( iom_use( "uice_mv"  ) ) CALL iom_put( "uice_mv"     ,   u_ice(:,:)               * zswi(:,:) + zmiss(:,:) )   ! ice velocity u component
298      IF  ( iom_use( "vice_mv"  ) ) CALL iom_put( "vice_mv"     ,   v_ice(:,:)               * zswi(:,:) + zmiss(:,:) )   ! ice velocity v component
299     
300      IF  ( iom_use( "xmtrpice" ) ) CALL iom_put( "xmtrpice"     ,  diag_xmtrp_ice(:,:)      )                            ! X-component of sea-ice mass transport (kg/s)
301      IF  ( iom_use( "ymtrpice" ) ) CALL iom_put( "ymtrpice"     ,  diag_ymtrp_ice(:,:)      )                            ! Y-component of sea-ice mass transport
302
303      IF  ( iom_use( "xmtrpsnw" ) ) CALL iom_put( "xmtrpsnw"     ,  diag_xmtrp_snw(:,:)      )                            ! X-component of snow mass transport (kg/s)
304      IF  ( iom_use( "ymtrpsnw" ) ) CALL iom_put( "ymtrpsnw"     ,  diag_ymtrp_snw(:,:)      )                            ! Y-component of snow mass transport
305
306      IF  ( iom_use( "xatrp"    ) ) CALL iom_put( "xatrp"        ,  diag_xatrp(:,:)              )                        ! X-component of ice area transport
307      IF  ( iom_use( "yatrp"    ) ) CALL iom_put( "yatrp"        ,  diag_yatrp(:,:)              )                        ! Y-component of ice area transport
308
309      IF  ( iom_use("utau_ai_mv") ) CALL iom_put( "utau_ai_mv"   ,  utau_ice(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Wind stress term in force balance (x)
310      IF  ( iom_use("vtau_ai_mv") ) CALL iom_put( "vtau_ai_mv"   ,  vtau_ice(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Wind stress term in force balance (y)
311
312      IF  ( iom_use( "utau_oi"  ) ) CALL iom_put( "utau_oi"     ,   diag_utau_oi(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Ocean stress term in force balance (x)
313      IF  ( iom_use( "vtau_oi"  ) ) CALL iom_put( "vtau_oi"     ,   diag_vtau_oi(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Ocean stress term in force balance (y)
314
315      IF  ( iom_use( "icestr_mv") ) CALL iom_put( "icestr_mv"   ,   strength(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Ice strength
316
317      IF  ( iom_use( "dssh_dx"  ) ) CALL iom_put( "dssh_dx"     ,   diag_dssh_dx(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Sea-surface tilt term in force balance (x)
318      IF  ( iom_use( "dssh_dy"  ) ) CALL iom_put( "dssh_dy"     ,   diag_dssh_dy(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Sea-surface tilt term in force balance (y)
319
320      IF  ( iom_use( "corstrx"  ) ) CALL iom_put( "corstrx"     ,   diag_corstrx(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Coriolis force term in force balance (x)
321      IF  ( iom_use( "corstry"  ) ) CALL iom_put( "corstry"     ,   diag_corstry(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Coriolis force term in force balance (y)
322
323      IF  ( iom_use( "intstrx"  ) ) CALL iom_put( "intstrx"     ,   diag_intstrx(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Internal force term in force balance (x)
324      IF  ( iom_use( "intstry"  ) ) CALL iom_put( "intstry"     ,   diag_intstry(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Internal force term in force balance (y)
325
326      IF  ( iom_use( "normstr"  ) ) CALL iom_put( "normstr"     ,   diag_sig1(:,:) * zswi(:,:)   )                        ! Normal stress
327      IF  ( iom_use( "sheastr"  ) ) CALL iom_put( "sheastr"     ,   diag_sig2(:,:) * zswi(:,:)   )                        ! Shear stress
328
329      !--------------------------------
330      ! Global ice diagnostics (SIMIP)
331      !--------------------------------
332
333      IF ( iom_use( "NH_icearea" ) .OR. iom_use( "NH_icevolu" ) .OR. iom_use( "NH_iceextt" ) )   THEN   ! NH integrated diagnostics
334 
335         WHERE( fcor > 0._wp ); zswi(:,:) = 1.0e-12
336         ELSEWHERE            ; zswi(:,:) = 0.
337         END WHERE
338
339         zdiag_area_nh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) )
340         zdiag_volu_nh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) )
341
342         WHERE( fcor > 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12
343         ELSEWHERE                              ; zswi(:,:) = 0.
344         END WHERE
345
346         zdiag_extt_nh = glob_sum( zswi(:,:) * e12t(:,:) )
347
348         IF ( iom_use( "NH_icearea" ) ) CALL iom_put( "NH_icearea" ,  zdiag_area_nh  )
349         IF ( iom_use( "NH_icevolu" ) ) CALL iom_put( "NH_icevolu" ,  zdiag_volu_nh  )
350         IF ( iom_use( "NH_iceextt" ) ) CALL iom_put( "NH_iceextt" ,  zdiag_extt_nh  )
351
352      ENDIF
353
354      IF ( iom_use( "SH_icearea" ) .OR. iom_use( "SH_icevolu" ) .OR. iom_use( "SH_iceextt" ) )   THEN   ! SH integrated diagnostics
355
356         WHERE( fcor < 0._wp ); zswi(:,:) = 1.0e-12; 
357         ELSEWHERE            ; zswi(:,:) = 0.
358         END WHERE
359
360         zdiag_area_sh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) ) 
361         zdiag_volu_sh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) )
362
363         WHERE( fcor < 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12
364         ELSEWHERE                              ; zswi(:,:) = 0.
365         END WHERE
366
367         zdiag_extt_sh = glob_sum( zswi(:,:) * e12t(:,:) )
368
369         IF ( iom_use( "SH_icearea" ) ) CALL iom_put( "SH_icearea", zdiag_area_sh )
370         IF ( iom_use( "SH_icevolu" ) ) CALL iom_put( "SH_icevolu", zdiag_volu_sh )
371         IF ( iom_use( "SH_iceextt" ) ) CALL iom_put( "SH_iceextt", zdiag_extt_sh )
372
373      ENDIF 
374
375      !     !  Create an output files (output.lim.abort.nc) if S < 0 or u > 20 m/s
376      !     IF( kindic < 0 )   CALL lim_wri_state( 'output.abort' )
377      !     not yet implemented
378     
379      CALL wrk_dealloc( jpi, jpj, jpl, zswi2 )
380      CALL wrk_dealloc( jpi, jpj     , z2d, zswi )
381      CALL wrk_dealloc( jpi, jpj     , zfb, zamask, zamask15 )
382      CALL wrk_dealloc( jpi, jpj     , zmiss                 )
383
384      IF( nn_timing == 1 )  CALL timing_stop('limwri')
385     
386   END SUBROUTINE lim_wri
387
388 
389   SUBROUTINE lim_wri_state( kt, kid, kh_i )
390      !!---------------------------------------------------------------------
391      !!                 ***  ROUTINE lim_wri_state  ***
392      !!       
393      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
394      !!      the instantaneous ice state and forcing fields for ice model
395      !!        Used to find errors in the initial state or save the last
396      !!      ocean state in case of abnormal end of a simulation
397      !!
398      !! History :
399      !!   4.0  !  2013-06  (C. Rousset)
400      !!----------------------------------------------------------------------
401      INTEGER, INTENT( in )   ::   kt               ! ocean time-step index)
402      INTEGER, INTENT( in )   ::   kid , kh_i
403      INTEGER                 ::   nz_i, jl
404      REAL(wp), DIMENSION(jpl) :: jcat
405      !!----------------------------------------------------------------------
406      DO jl = 1, jpl
407         jcat(jl) = REAL(jl)
408      ENDDO
409     
410      CALL histvert( kid, "ncatice", "Ice Categories","", jpl, jcat, nz_i, "up")
411
412      CALL histdef( kid, "sithic", "Ice thickness"           , "m"      ,   &
413      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
414      CALL histdef( kid, "siconc", "Ice concentration"       , "%"      ,   &
415      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
416      CALL histdef( kid, "sitemp", "Ice temperature"         , "C"      ,   &
417      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
418      CALL histdef( kid, "sivelu", "i-Ice speed "            , "m/s"    ,   &
419      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
420      CALL histdef( kid, "sivelv", "j-Ice speed "            , "m/s"    ,   &
421      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
422      CALL histdef( kid, "sistru", "i-Wind stress over ice " , "Pa"     ,   &
423      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
424      CALL histdef( kid, "sistrv", "j-Wind stress over ice " , "Pa"     ,   &
425      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
426      CALL histdef( kid, "sisflx", "Solar flux over ocean"     , "w/m2" ,   &
427      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
428      CALL histdef( kid, "sinflx", "Non-solar flux over ocean" , "w/m2" ,   &
429      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
430      CALL histdef( kid, "isnowpre", "Snow precipitation"      , "kg/m2/s",   &
431      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
432      CALL histdef( kid, "sisali", "Ice salinity"            , "PSU"    ,   &
433      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
434      CALL histdef( kid, "sivolu", "Ice volume"              , "m"      ,   &
435      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
436      CALL histdef( kid, "sidive", "Ice divergence"          , "10-8s-1",   &
437      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
438
439      CALL histdef( kid, "vfxbog", "Ice bottom production"   , "m/s"    ,   &
440      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
441      CALL histdef( kid, "vfxdyn", "Ice dynamic production"  , "m/s"    ,   &
442      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
443      CALL histdef( kid, "vfxopw", "Ice open water prod"     , "m/s"    ,   &
444      &       jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
445      CALL histdef( kid, "vfxsni", "Snow ice production "    , "m/s"    ,   &
446      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
447      CALL histdef( kid, "vfxres", "Ice prod from limupdate" , "m/s"    ,   &
448      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
449      CALL histdef( kid, "vfxbom", "Ice bottom melt"         , "m/s"    ,   &
450      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
451      CALL histdef( kid, "vfxsum", "Ice surface melt"        , "m/s"    ,   &
452      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
453
454      CALL histdef( kid, "sithicat", "Ice thickness"         , "m"      ,   &
455      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
456      CALL histdef( kid, "siconcat", "Ice concentration"     , "%"      ,   &
457      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
458      CALL histdef( kid, "sisalcat", "Ice salinity"           , ""      ,   &
459      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
460      CALL histdef( kid, "sitemcat", "Ice temperature"       , "C"      ,   &
461      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
462      CALL histdef( kid, "snthicat", "Snw thickness"         , "m"      ,   &
463      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
464      CALL histdef( kid, "sntemcat", "Snw temperature"       , "C"      ,   &
465      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
466
467      CALL histend( kid, snc4set )   ! end of the file definition
468
469      CALL histwrite( kid, "sithic", kt, htm_i         , jpi*jpj, (/1/) )   
470      CALL histwrite( kid, "siconc", kt, at_i          , jpi*jpj, (/1/) )
471      CALL histwrite( kid, "sitemp", kt, tm_i - rt0    , jpi*jpj, (/1/) )
472      CALL histwrite( kid, "sivelu", kt, u_ice          , jpi*jpj, (/1/) )
473      CALL histwrite( kid, "sivelv", kt, v_ice          , jpi*jpj, (/1/) )
474      CALL histwrite( kid, "sistru", kt, utau_ice       , jpi*jpj, (/1/) )
475      CALL histwrite( kid, "sistrv", kt, vtau_ice       , jpi*jpj, (/1/) )
476      CALL histwrite( kid, "sisflx", kt, qsr , jpi*jpj, (/1/) )
477      CALL histwrite( kid, "sinflx", kt, qns , jpi*jpj, (/1/) )
478      CALL histwrite( kid, "isnowpre", kt, sprecip        , jpi*jpj, (/1/) )
479      CALL histwrite( kid, "sisali", kt, smt_i          , jpi*jpj, (/1/) )
480      CALL histwrite( kid, "sivolu", kt, vt_i           , jpi*jpj, (/1/) )
481      CALL histwrite( kid, "sidive", kt, divu_i*1.0e8   , jpi*jpj, (/1/) )
482
483      CALL histwrite( kid, "vfxbog", kt, wfx_bog        , jpi*jpj, (/1/) )
484      CALL histwrite( kid, "vfxdyn", kt, wfx_dyn        , jpi*jpj, (/1/) )
485      CALL histwrite( kid, "vfxopw", kt, wfx_opw        , jpi*jpj, (/1/) )
486      CALL histwrite( kid, "vfxsni", kt, wfx_sni        , jpi*jpj, (/1/) )
487      CALL histwrite( kid, "vfxres", kt, wfx_res        , jpi*jpj, (/1/) )
488      CALL histwrite( kid, "vfxbom", kt, wfx_bom        , jpi*jpj, (/1/) )
489      CALL histwrite( kid, "vfxsum", kt, wfx_sum        , jpi*jpj, (/1/) )
490
491      CALL histwrite( kid, "sithicat", kt, ht_i        , jpi*jpj*jpl, (/1/) )   
492      CALL histwrite( kid, "siconcat", kt, a_i         , jpi*jpj*jpl, (/1/) )   
493      CALL histwrite( kid, "sisalcat", kt, sm_i        , jpi*jpj*jpl, (/1/) )   
494      CALL histwrite( kid, "sitemcat", kt, tm_i - rt0  , jpi*jpj*jpl, (/1/) )   
495      CALL histwrite( kid, "snthicat", kt, ht_s        , jpi*jpj*jpl, (/1/) )   
496      CALL histwrite( kid, "sntemcat", kt, tm_su - rt0 , jpi*jpj*jpl, (/1/) )   
497
498      ! Close the file
499      ! -----------------
500      !CALL histclo( kid )
501
502    END SUBROUTINE lim_wri_state
503
504#else
505   !!----------------------------------------------------------------------
506   !!   Default option :         Empty module          NO LIM sea-ice model
507   !!----------------------------------------------------------------------
508CONTAINS
509   SUBROUTINE lim_wri          ! Empty routine
510   END SUBROUTINE lim_wri
511#endif
512
513   !!======================================================================
514END MODULE limwri
Note: See TracBrowser for help on using the repository browser.