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/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90 @ 8272

Last change on this file since 8272 was 8272, checked in by clem, 7 years ago

update following r8270 & r8271
  • Property svn:keywords set to Id
File size: 34.9 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 limvar
21   USE in_out_manager
22   USE lbclnk
23   USE lib_mpp         ! MPP library
24   USE wrk_nemo        ! work arrays
25   USE iom
26   USE timing          ! Timing
27   USE lib_fortran     ! Fortran utilities
28
29   IMPLICIT NONE
30   PRIVATE
31
32   PUBLIC lim_wri        ! routine called by lim_step.F90
33   PUBLIC lim_wri_state  ! called by dia_wri_state
34
35   !!----------------------------------------------------------------------
36   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
37   !! $Id$
38   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
39   !!----------------------------------------------------------------------
40CONTAINS
41
42   SUBROUTINE lim_wri( kindic )
43      !!-------------------------------------------------------------------
44      !!  This routine computes the average of some variables and write it
45      !!  on the ouput files.
46      !!  ATTENTION cette routine n'est valable que si le pas de temps est
47      !!  egale a une fraction entiere de 1 jours.
48      !!  Diff 1-D 3-D : suppress common also included in etat
49      !!                 suppress cmoymo 11-18
50      !!  modif : 03/06/98
51      !!-------------------------------------------------------------------
52      INTEGER, INTENT(in) ::   kindic   ! if kindic < 0 there has been an error somewhere
53      !
54      INTEGER  ::  ji, jj, jk, jl  ! dummy loop indices
55      REAL(wp) ::  z2da, z2db, ztmp, zrho1, zrho2, zmiss_val
56      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zswi2, zmiss2
57      REAL(wp), POINTER, DIMENSION(:,:)   ::  z2d, zswi, zmiss !  2D workspace
58      REAL(wp), POINTER, DIMENSION(:,:)   ::  zfb              ! ice freeboard
59      REAL(wp), POINTER, DIMENSION(:,:)   ::  zamask, zamask15 ! 15% concentration mask
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      !!-------------------------------------------------------------------
70
71      IF( nn_timing == 1 )  CALL timing_start('limwri')
72
73      CALL wrk_alloc( jpi, jpj, jpl, zswi2, zmiss2 )
74      CALL wrk_alloc( jpi, jpj     , z2d, zswi, zmiss )
75      CALL wrk_alloc( jpi, jpj     , zfb, zamask, zamask15 )
76
77      !----------------------------------------
78      ! Brine volume, switches, missing values
79      !----------------------------------------
80
81      ! brine volume
82      CALL lim_var_bv 
83
84      ! tresholds for outputs
85      DO jj = 1, jpj
86         DO ji = 1, jpi
87            zswi(ji,jj)      = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice, 0 if no ice
88            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
89            zamask15(ji,jj)  = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - 0.15   ) ) ! 1 if 15% ice, 0 if less
90         END DO
91      END DO
92      DO jl = 1, jpl
93         DO jj = 1, jpj
94            DO ji = 1, jpi
95               zswi2(ji,jj,jl)  = MAX( 0._wp , SIGN( 1._wp , a_i(ji,jj,jl) - epsi06 ) )
96            END DO
97         END DO
98      END DO
99
100      zmiss_val     = -1.0e20
101      zmiss(:,:)    = zmiss_val * ( 1. - zswi(:,:) )
102      zmiss2(:,:,:) = zmiss_val * ( 1. - zswi2(:,:,:) )
103
104      !----------------------------------------
105      ! Standard outputs
106      !----------------------------------------
107      ! fluxes
108      ! pfrld is the lead fraction at the previous time step (actually between TRP and THD)
109      IF( iom_use('qsr_oce') )   CALL iom_put( "qsr_oce" , qsr_oce(:,:) * pfrld(:,:) )                                   !     solar flux at ocean surface
110      IF( iom_use('qns_oce') )   CALL iom_put( "qns_oce" , qns_oce(:,:) * pfrld(:,:) + qemp_oce(:,:) )                   ! non-solar flux at ocean surface
111      IF( iom_use('qsr_ice') )   CALL iom_put( "qsr_ice" , SUM( qsr_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) )                 !     solar flux at ice surface
112      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
113      IF( iom_use('qtr_ice') )   CALL iom_put( "qtr_ice" , SUM( ftr_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) )                 !     solar flux transmitted thru ice
114      IF( iom_use('qt_oce' ) )   CALL iom_put( "qt_oce"  , ( qsr_oce(:,:) + qns_oce(:,:) ) * pfrld(:,:) + qemp_oce(:,:) ) 
115      IF( iom_use('qt_ice' ) )   CALL iom_put( "qt_ice"  , SUM( ( qns_ice(:,:,:) + qsr_ice(:,:,:) )   &
116         &                                                      * a_i_b(:,:,:), dim=3 ) + qemp_ice(:,:) )
117      IF( iom_use('qemp_oce') )  CALL iom_put( "qemp_oce" , qemp_oce(:,:) ) 
118      IF( iom_use('qemp_ice') )  CALL iom_put( "qemp_ice" , qemp_ice(:,:) ) 
119      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)
120      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)
121
122      ! velocity
123      IF ( iom_use('uice_ipa') ) CALL iom_put( "uice_ipa"     , u_ice      )       ! ice velocity u component
124      IF ( iom_use('vice_ipa') ) CALL iom_put( "vice_ipa"     , v_ice      )       ! ice velocity v component
125
126      IF ( ( iom_use( "icevel" ) ) .OR. ( iom_use( "icevel_mv" ) ) ) THEN
127         DO jj = 2 , jpjm1
128            DO ji = 2 , jpim1
129               z2da  = ( u_ice(ji,jj) * umask(ji,jj,1) + u_ice(ji-1,jj) * umask(ji-1,jj,1) ) * 0.5_wp
130               z2db  = ( v_ice(ji,jj) * vmask(ji,jj,1) + v_ice(ji,jj-1) * vmask(ji,jj-1,1) ) * 0.5_wp
131               z2d(ji,jj) = SQRT( z2da * z2da + z2db * z2db )
132           END DO
133         END DO
134         CALL lbc_lnk( z2d, 'T', 1. )
135         IF ( iom_use( "icevel" )  )   CALL iom_put( "icevel"       , z2d        )                                            ! ice velocity module
136         IF ( iom_use( "icevel_mv" ) ) CALL iom_put( "icevel_mv"    , z2d(:,:) * zswi(:,:) + zmiss(:,:) )                     ! ice velocity module (missing value)
137      ENDIF
138
139      IF ( iom_use( "tau_icebfr" ) )    CALL iom_put( "tau_icebfr"  , tau_icebfr             )  ! ice friction with ocean bottom (landfast ice) 
140      !
141      IF ( iom_use( "miceage" ) )       CALL iom_put( "miceage"     , om_i * zswi * zamask15 )  ! mean ice age
142      IF ( iom_use( "micet" ) )         CALL iom_put( "micet"       , ( tm_i  - rt0 ) * zswi )  ! ice mean    temperature
143      IF ( iom_use( "icest" ) )         CALL iom_put( "icest"       , ( tm_su - rt0 ) * zswi )  ! ice surface temperature
144      IF ( iom_use( "icecolf" ) )       CALL iom_put( "icecolf"     , hicol                  )  ! frazil ice collection thickness
145      !
146      CALL iom_put( "isst"        , sst_m               )        ! sea surface temperature
147      CALL iom_put( "isss"        , sss_m               )        ! sea surface salinity
148      CALL iom_put( "iceconc"     , at_i  * zswi        )        ! ice concentration
149      CALL iom_put( "icevolu"     , vt_i  * zswi        )        ! ice volume = mean ice thickness over the cell
150      CALL iom_put( "icehc"       , et_i  * zswi        )        ! ice total heat content
151      CALL iom_put( "isnowhc"     , et_s  * zswi        )        ! snow total heat content
152      CALL iom_put( "ibrinv"      , bvm_i * zswi * 100. )        ! brine volume
153      CALL iom_put( "snowpre"     , sprecip * 86400.    )        ! snow precipitation
154      CALL iom_put( "micesalt"    , smt_i   * zswi      )        ! mean ice salinity
155      CALL iom_put( "snowvol"     , vt_s    * zswi      )        ! snow volume
156
157      CALL iom_put( "idive"       , divu_i(:,:)  * zswi(:,:) + zmiss(:,:) )   ! divergence
158      CALL iom_put( "ishear"      , shear_i(:,:) * zswi(:,:) + zmiss(:,:) )   ! shear
159     
160      CALL iom_put( "icetrp"      , diag_trp_vi * rday  )        ! ice volume transport
161      CALL iom_put( "snwtrp"      , diag_trp_vs * rday  )        ! snw volume transport
162      CALL iom_put( "saltrp"      , diag_trp_smv * rday * rhoic ) ! salt content transport
163      CALL iom_put( "deitrp"      , diag_trp_ei         )        ! advected ice enthalpy (W/m2)
164      CALL iom_put( "destrp"      , diag_trp_es         )        ! advected snw enthalpy (W/m2)
165
166      CALL iom_put( "sfxbog"      , sfx_bog * rday      )        ! salt flux from bottom growth
167      CALL iom_put( "sfxbom"      , sfx_bom * rday      )        ! salt flux from bottom melting
168      CALL iom_put( "sfxsum"      , sfx_sum * rday      )        ! salt flux from surface melting
169      CALL iom_put( "sfxlam"      , sfx_lam * rday      )        ! salt flux from lateral melting
170      CALL iom_put( "sfxsni"      , sfx_sni * rday      )        ! salt flux from snow ice formation
171      CALL iom_put( "sfxopw"      , sfx_opw * rday      )        ! salt flux from open water formation
172      CALL iom_put( "sfxdyn"      , sfx_dyn * rday      )        ! salt flux from ridging rafting
173      CALL iom_put( "sfxres"      , sfx_res * rday      )        ! salt flux from limupdate (resultant)
174      CALL iom_put( "sfxbri"      , sfx_bri * rday      )        ! salt flux from brines
175      CALL iom_put( "sfxsub"      , sfx_sub * rday      )        ! salt flux from sublimation
176      CALL iom_put( "sfx"         , sfx     * rday      )        ! total salt flux
177
178      ztmp = rday / rhoic
179      CALL iom_put( "vfxres"     , wfx_res * ztmp       )        ! daily prod./melting due to limupdate
180      CALL iom_put( "vfxopw"     , wfx_opw * ztmp       )        ! daily lateral thermodynamic ice production
181      CALL iom_put( "vfxsni"     , wfx_sni * ztmp       )        ! daily snowice ice production
182      CALL iom_put( "vfxbog"     , wfx_bog * ztmp       )        ! daily bottom thermodynamic ice production
183      CALL iom_put( "vfxdyn"     , wfx_dyn * ztmp       )        ! daily dynamic ice production (rid/raft)
184      CALL iom_put( "vfxsum"     , wfx_sum * ztmp       )        ! surface melt
185      CALL iom_put( "vfxbom"     , wfx_bom * ztmp       )        ! bottom melt
186      CALL iom_put( "vfxlam"     , wfx_lam * ztmp       )        ! lateral melt
187      CALL iom_put( "vfxice"     , wfx_ice * ztmp       )        ! total ice growth/melt
188
189      IF ( ln_pnd ) &
190         CALL iom_put( "vfxpnd"  , wfx_pnd * ztmp       )        ! melt pond water flux
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(:,:) * at_i_b(:,:) ) ! 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      ! MV MP 2016
229      IF ( ln_pnd ) THEN
230         CALL iom_put( "iceamp"  , at_ip  * zswi        )   ! melt pond total fraction
231         CALL iom_put( "icevmp"  , vt_ip  * zswi        )   ! melt pond total volume per unit area
232      ENDIF
233      ! END MV MP 2016
234
235      !----------------------------------
236      ! Output category-dependent fields
237      !----------------------------------
238      IF ( iom_use( "iceconc_cat"  ) )  CALL iom_put( "iceconc_cat"      , a_i   * zswi2   )        ! area for categories
239      IF ( iom_use( "icethic_cat"  ) )  CALL iom_put( "icethic_cat"      , ht_i  * zswi2   )        ! thickness for categories
240      IF ( iom_use( "snowthic_cat" ) )  CALL iom_put( "snowthic_cat"     , ht_s  * zswi2   )        ! snow depth for categories
241      IF ( iom_use( "salinity_cat" ) )  CALL iom_put( "salinity_cat"     , sm_i  * zswi2   )        ! salinity for categories
242      ! ice temperature
243      IF ( iom_use( "icetemp_cat"  ) )  CALL iom_put( "icetemp_cat", ( SUM( t_i(:,:,:,:), dim=3 ) * r1_nlay_i - rt0 ) * zswi2 )
244      ! snow temperature
245      IF ( iom_use( "snwtemp_cat"  ) )  CALL iom_put( "snwtemp_cat", ( SUM( t_s(:,:,:,:), dim=3 ) * r1_nlay_s - rt0 ) * zswi2 )
246      ! ice age
247      IF ( iom_use( "iceage_cat"   ) )  CALL iom_put( "iceage_cat" , o_i * zswi2 ) 
248      ! brine volume
249      IF ( iom_use( "brinevol_cat" ) )  CALL iom_put( "brinevol_cat", bv_i * 100. * zswi2 )
250
251      ! MV MP 2016
252      IF ( ln_pnd ) THEN
253         IF ( iom_use( "iceamp_cat"  ) )  CALL iom_put( "iceamp_cat"     , a_ip       * zswi2   )       ! melt pond frac for categories
254         IF ( iom_use( "icevmp_cat"  ) )  CALL iom_put( "icevmp_cat"     , v_ip       * zswi2   )       ! melt pond frac for categories
255         IF ( iom_use( "icehmp_cat"  ) )  CALL iom_put( "icehmp_cat"     , h_ip       * zswi2   )       ! melt pond frac for categories
256         IF ( iom_use( "iceafp_cat"  ) )  CALL iom_put( "iceafp_cat"     , a_ip_frac  * zswi2   )       ! melt pond frac for categories
257      ENDIF
258      ! END MV MP 2016
259
260      !--------------------------------
261      ! Add-ons for SIMIP
262      !--------------------------------
263      zrho1 = ( rau0 - rhoic ) / rau0; zrho2 = rhosn / rau0
264
265      IF  ( iom_use( "icepres"  ) ) CALL iom_put( "icepres"     , zswi(:,:)                     )                                ! Ice presence (1 or 0)
266      IF  ( iom_use( "icemass"  ) ) CALL iom_put( "icemass"     , rhoic * vt_i(:,:) * zswi(:,:) )                                ! Ice mass per cell area
267      IF  ( iom_use( "icethic"  ) ) CALL iom_put( "icethic"     , htm_i(:,:) * zamask(:,:)  + ( 1. - zamask(:,:) ) * zmiss_val )     ! Ice thickness
268      IF  ( iom_use( "snomass"  ) ) CALL iom_put( "snomass"     , rhosn * vt_s(:,:)         * zswi(:,:) + zmiss(:,:) )           ! Snow mass per cell area
269      IF  ( iom_use( "snothic"  ) ) CALL iom_put( "snothic"     , htm_s(:,:) * zamask(:,:)  + ( 1. - zamask(:,:) ) * zmiss_val )     ! Snow thickness       
270
271      IF  ( iom_use( "iceconc_cat_mv"  ) )  CALL iom_put( "iceconc_cat_mv" , a_i(:,:,:)  * zswi2(:,:,:) + zmiss2(:,:,:) )        ! Area for categories
272      IF  ( iom_use( "icethic_cat_mv"  ) )  CALL iom_put( "icethic_cat_mv" , ht_i(:,:,:) * zswi2(:,:,:) + zmiss2(:,:,:) )        ! Thickness for categories
273      IF  ( iom_use( "snowthic_cat_mv" ) )  CALL iom_put( "snowthic_cat_mv", ht_s(:,:,:) * zswi2(:,:,:) + zmiss2(:,:,:) )        ! Snow depth for categories
274
275      IF  ( iom_use( "icestK"   ) ) CALL iom_put( "icestK"      , tm_su(:,:)                * zswi(:,:) + zmiss(:,:) )           ! Ice surface temperature
276      IF  ( iom_use( "icesntK"  ) ) CALL iom_put( "icesntK"     , tm_si(:,:)                * zswi(:,:) + zmiss(:,:) )           ! Snow-ice interface temperature
277      IF  ( iom_use( "icebotK"  ) ) CALL iom_put( "icebotK"     , t_bo(:,:)                 * zswi(:,:) + zmiss(:,:) )           ! Ice bottom temperature
278      IF  ( iom_use( "iceage"   ) ) CALL iom_put( "iceage"      , om_i(:,:) * zamask15(:,:) + ( 1. - zamask15(:,:) ) * zmiss_val )   ! Ice age
279      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
280      IF  ( iom_use( "icesal"   ) ) CALL iom_put( "icesal"      , smt_i(:,:)                * zswi(:,:) + zmiss(:,:) )           ! Ice salinity
281
282      IF  ( iom_use( "icefb"    ) ) THEN
283         zfb(:,:) = ( zrho1 * htm_i(:,:) - zrho2 * htm_s(:,:) )                                         
284         WHERE( zfb < 0._wp ) ;   zfb = 0._wp ;   END WHERE
285                                    CALL iom_put( "icefb"       , zfb(:,:)                  * zswi(:,:) + zmiss(:,:) )           ! Ice freeboard
286      ENDIF
287
288      IF  ( iom_use( "isnhcneg" ) ) CALL iom_put( "isnhcneg"    , - et_s(:,:)               * zswi(:,:) + zmiss(:,:) )           ! Snow total heat content
289
290      IF  ( iom_use( "dmithd"   ) ) CALL iom_put( "dmithd"      , - wfx_bog - wfx_bom - wfx_sum   &                       ! Sea-ice mass change from thermodynamics
291              &                     - wfx_sni - wfx_opw - wfx_res )
292      IF  ( iom_use( "dmidyn"   ) ) CALL iom_put( "dmidyn"      ,   diag_dmi_dyn             )                            ! Sea-ice mass change from dynamics
293      IF  ( iom_use( "dmiopw"   ) ) CALL iom_put( "dmiopw"      , - wfx_opw                  )                            ! Sea-ice mass change through growth in open water
294      IF  ( iom_use( "dmibog"   ) ) CALL iom_put( "dmibog"      , - wfx_bog                  )                            ! Sea-ice mass change through basal growth
295      IF  ( iom_use( "dmisni"   ) ) CALL iom_put( "dmisni"      , - wfx_sni                  )                            ! Sea-ice mass change through snow-to-ice conversion
296      IF  ( iom_use( "dmisum"   ) ) CALL iom_put( "dmisum"      , - wfx_sum                  )                            ! Sea-ice mass change through surface melting
297      IF  ( iom_use( "dmibom"   ) ) CALL iom_put( "dmibom"      , - wfx_bom                  )                            ! Sea-ice mass change through bottom melting
298
299      IF  ( iom_use( "dmtsub"   ) ) CALL iom_put( "dmtsub"      , - wfx_sub                  )                            ! Sea-ice mass change through evaporation and sublimation
300      IF  ( iom_use( "dmssub"   ) ) CALL iom_put( "dmssub"      , - wfx_snw_sub              )                            ! Snow mass change through sublimation
301      IF  ( iom_use( "dmisub"   ) ) CALL iom_put( "dmisub"      , - wfx_ice_sub              )                            ! Sea-ice mass change through sublimation
302
303      IF  ( iom_use( "dmsspr"   ) ) CALL iom_put( "dmsspr"      , - wfx_spr                  )                            ! Snow mass change through snow fall
304      IF  ( iom_use( "dmsssi"   ) ) CALL iom_put( "dmsssi"      ,   wfx_sni*rhosn/rhoic      )                            ! Snow mass change through snow-to-ice conversion
305
306      IF  ( iom_use( "dmsmel"   ) ) CALL iom_put( "dmsmel"      , - wfx_snw_sum              )                            ! Snow mass change through melt
307      IF  ( iom_use( "dmsdyn"   ) ) CALL iom_put( "dmsdyn"      ,   diag_dms_dyn             )                            ! Snow mass change through dynamics
308
309      IF  ( iom_use( "hfxsenso" ) ) CALL iom_put( "hfxsenso"    ,   -fhtur(:,:)              * zswi(:,:) + zmiss(:,:) )   ! Sensible oceanic heat flux
310      IF  ( iom_use( "hfxconbo" ) ) CALL iom_put( "hfxconbo"    ,   diag_fc_bo               * zswi(:,:) + zmiss(:,:) )   ! Bottom conduction flux
311      IF  ( iom_use( "hfxconsu" ) ) CALL iom_put( "hfxconsu"    ,   diag_fc_su               * zswi(:,:) + zmiss(:,:) )   ! Surface conduction flux
312
313      IF  ( iom_use( "wfxtot"   ) ) CALL iom_put( "wfxtot"      ,   wfx_ice(:,:)             * zswi(:,:) + zmiss(:,:) )   ! Total freshwater flux from sea ice
314      IF  ( iom_use( "wfxsum"   ) ) CALL iom_put( "wfxsum"      ,   wfx_sum(:,:)             * zswi(:,:) + zmiss(:,:) )   ! Freshwater flux from sea-ice surface
315      IF  ( iom_use( "sfx_mv"   ) ) CALL iom_put( "sfx_mv"      ,   sfx(:,:) * 0.001         * zswi(:,:) + zmiss(:,:) )   ! Total salt flux
316
317      IF  ( iom_use( "uice_mv"  ) ) CALL iom_put( "uice_mv"     ,   u_ice(:,:)               * zswi(:,:) + zmiss(:,:) )   ! ice velocity u component
318      IF  ( iom_use( "vice_mv"  ) ) CALL iom_put( "vice_mv"     ,   v_ice(:,:)               * zswi(:,:) + zmiss(:,:) )   ! ice velocity v component
319     
320      IF  ( iom_use( "xmtrpice" ) ) CALL iom_put( "xmtrpice"     ,  diag_xmtrp_ice(:,:)      )                            ! X-component of sea-ice mass transport (kg/s)
321      IF  ( iom_use( "ymtrpice" ) ) CALL iom_put( "ymtrpice"     ,  diag_ymtrp_ice(:,:)      )                            ! Y-component of sea-ice mass transport
322
323      IF  ( iom_use( "xmtrpsnw" ) ) CALL iom_put( "xmtrpsnw"     ,  diag_xmtrp_snw(:,:)      )                            ! X-component of snow mass transport (kg/s)
324      IF  ( iom_use( "ymtrpsnw" ) ) CALL iom_put( "ymtrpsnw"     ,  diag_ymtrp_snw(:,:)      )                            ! Y-component of snow mass transport
325
326      IF  ( iom_use( "xatrp"    ) ) CALL iom_put( "xatrp"        ,  diag_xatrp(:,:)              )                        ! X-component of ice area transport
327      IF  ( iom_use( "yatrp"    ) ) CALL iom_put( "yatrp"        ,  diag_yatrp(:,:)              )                        ! Y-component of ice area transport
328
329      IF  ( iom_use( "utau_ice" ) ) CALL iom_put( "utau_ice"     ,  utau_ice(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Wind stress term in force balance (x)
330      IF  ( iom_use( "vtau_ice" ) ) CALL iom_put( "vtau_ice"     ,  vtau_ice(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Wind stress term in force balance (y)
331
332      IF  ( iom_use( "utau_oi"  ) ) CALL iom_put( "utau_oi"     ,   diag_utau_oi(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Ocean stress term in force balance (x)
333      IF  ( iom_use( "vtau_oi"  ) ) CALL iom_put( "vtau_oi"     ,   diag_vtau_oi(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Ocean stress term in force balance (y)
334
335      IF  ( iom_use( "icestr"   ) ) CALL iom_put( "icestr"      ,   strength(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Ice strength
336
337      IF  ( iom_use( "dssh_dx"  ) ) CALL iom_put( "dssh_dx"     ,   diag_dssh_dx(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Sea-surface tilt term in force balance (x)
338      IF  ( iom_use( "dssh_dy"  ) ) CALL iom_put( "dssh_dy"     ,   diag_dssh_dy(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Sea-surface tilt term in force balance (y)
339
340      IF  ( iom_use( "corstrx"  ) ) CALL iom_put( "corstrx"     ,   diag_corstrx(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Coriolis force term in force balance (x)
341      IF  ( iom_use( "corstry"  ) ) CALL iom_put( "corstry"     ,   diag_corstry(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Coriolis force term in force balance (y)
342
343      IF  ( iom_use( "intstrx"  ) ) CALL iom_put( "intstrx"     ,   diag_intstrx(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Internal force term in force balance (x)
344      IF  ( iom_use( "intstry"  ) ) CALL iom_put( "intstry"     ,   diag_intstry(:,:)        * zswi(:,:) + zmiss(:,:) )   ! Internal force term in force balance (y)
345
346      IF  ( iom_use( "normstr"  ) ) CALL iom_put( "normstr"     ,   diag_sig1(:,:)           * zswi(:,:) + zmiss(:,:) )   ! Normal stress
347      IF  ( iom_use( "sheastr"  ) ) CALL iom_put( "sheastr"     ,   diag_sig2(:,:)           * zswi(:,:) + zmiss(:,:) )   ! Shear stress
348
349      !--------------------------------
350      ! Global ice diagnostics (SIMIP)
351      !--------------------------------
352
353      IF ( iom_use( "NH_icearea" ) .OR. iom_use( "NH_icevolu" ) .OR. iom_use( "NH_iceextt" ) )   THEN   ! NH integrated diagnostics
354 
355         WHERE( ff_t > 0._wp ); zswi(:,:) = 1.0e-12
356         ELSEWHERE            ; zswi(:,:) = 0.
357         END WHERE
358
359         zdiag_area_nh = glob_sum( at_i(:,:) * zswi(:,:) * e1e2t(:,:) )
360         zdiag_volu_nh = glob_sum( vt_i(:,:) * zswi(:,:) * e1e2t(:,:) )
361
362         WHERE( ff_t > 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12
363         ELSEWHERE                              ; zswi(:,:) = 0.
364         END WHERE
365
366         zdiag_extt_nh = glob_sum( zswi(:,:) * e1e2t(:,:) )
367
368         IF ( iom_use( "NH_icearea" ) ) CALL iom_put( "NH_icearea" ,  zdiag_area_nh  )
369         IF ( iom_use( "NH_icevolu" ) ) CALL iom_put( "NH_icevolu" ,  zdiag_volu_nh  )
370         IF ( iom_use( "NH_iceextt" ) ) CALL iom_put( "NH_iceextt" ,  zdiag_extt_nh  )
371
372      ENDIF
373
374      IF ( iom_use( "SH_icearea" ) .OR. iom_use( "SH_icevolu" ) .OR. iom_use( "SH_iceextt" ) )   THEN   ! SH integrated diagnostics
375
376         WHERE( ff_t < 0._wp ); zswi(:,:) = 1.0e-12; 
377         ELSEWHERE            ; zswi(:,:) = 0.
378         END WHERE
379
380         zdiag_area_sh = glob_sum( at_i(:,:) * zswi(:,:) * e1e2t(:,:) ) 
381         zdiag_volu_sh = glob_sum( vt_i(:,:) * zswi(:,:) * e1e2t(:,:) )
382
383         WHERE( ff_t < 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12
384         ELSEWHERE                              ; zswi(:,:) = 0.
385         END WHERE
386
387         zdiag_extt_sh = glob_sum( zswi(:,:) * e1e2t(:,:) )
388
389         IF ( iom_use( "SH_icearea" ) ) CALL iom_put( "SH_icearea", zdiag_area_sh )
390         IF ( iom_use( "SH_icevolu" ) ) CALL iom_put( "SH_icevolu", zdiag_volu_sh )
391         IF ( iom_use( "SH_iceextt" ) ) CALL iom_put( "SH_iceextt", zdiag_extt_sh )
392
393      ENDIF 
394
395      !     !  Create an output files (output.lim.abort.nc) if S < 0 or u > 20 m/s
396      !     IF( kindic < 0 )   CALL lim_wri_state( 'output.abort' )
397      !     not yet implemented
398     
399      CALL wrk_dealloc( jpi, jpj, jpl, zswi2, zmiss2 )
400      CALL wrk_dealloc( jpi, jpj     , z2d, zswi, zmiss )
401      CALL wrk_dealloc( jpi, jpj     , zfb, zamask, zamask15 )
402
403      IF( nn_timing == 1 )  CALL timing_stop('limwri')
404     
405   END SUBROUTINE lim_wri
406
407 
408   SUBROUTINE lim_wri_state( kt, kid, kh_i )
409      !!---------------------------------------------------------------------
410      !!                 ***  ROUTINE lim_wri_state  ***
411      !!       
412      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
413      !!      the instantaneous ice state and forcing fields for ice model
414      !!        Used to find errors in the initial state or save the last
415      !!      ocean state in case of abnormal end of a simulation
416      !!
417      !! History :
418      !!   4.0  !  2013-06  (C. Rousset)
419      !!----------------------------------------------------------------------
420      INTEGER, INTENT( in )   ::   kt               ! ocean time-step index)
421      INTEGER, INTENT( in )   ::   kid , kh_i
422      INTEGER                 ::   nz_i, jl
423      REAL(wp), DIMENSION(jpl) :: jcat
424      !!----------------------------------------------------------------------
425      DO jl = 1, jpl
426         jcat(jl) = REAL(jl)
427      ENDDO
428     
429      CALL histvert( kid, "ncatice", "Ice Categories","", jpl, jcat, nz_i, "up")
430
431      CALL histdef( kid, "sithic", "Ice thickness"           , "m"      ,   &
432      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
433      CALL histdef( kid, "siconc", "Ice concentration"       , "%"      ,   &
434      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
435      CALL histdef( kid, "sitemp", "Ice temperature"         , "C"      ,   &
436      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
437      CALL histdef( kid, "sivelu", "i-Ice speed "            , "m/s"    ,   &
438      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
439      CALL histdef( kid, "sivelv", "j-Ice speed "            , "m/s"    ,   &
440      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
441      CALL histdef( kid, "sistru", "i-Wind stress over ice " , "Pa"     ,   &
442      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
443      CALL histdef( kid, "sistrv", "j-Wind stress over ice " , "Pa"     ,   &
444      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
445      CALL histdef( kid, "sisflx", "Solar flux over ocean"     , "w/m2" ,   &
446      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
447      CALL histdef( kid, "sinflx", "Non-solar flux over ocean" , "w/m2" ,   &
448      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
449      CALL histdef( kid, "isnowpre", "Snow precipitation"      , "kg/m2/s",   &
450      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
451      CALL histdef( kid, "sisali", "Ice salinity"            , "PSU"    ,   &
452      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
453      CALL histdef( kid, "sivolu", "Ice volume"              , "m"      ,   &
454      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
455      CALL histdef( kid, "sidive", "Ice divergence"          , "10-8s-1",   &
456      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
457
458      ! MV MP 2016
459      IF ( ln_pnd ) THEN
460         CALL histdef( kid, "si_amp", "Melt pond fraction"      , "%"      ,   &
461      &         jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
462         CALL histdef( kid, "si_vmp", "Melt pond volume"        ,  "m"     ,   &
463      &         jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
464      ENDIF
465      ! END MV MP 2016
466
467      CALL histdef( kid, "vfxbog", "Ice bottom production"   , "m/s"    ,   &
468      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
469      CALL histdef( kid, "vfxdyn", "Ice dynamic production"  , "m/s"    ,   &
470      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
471      CALL histdef( kid, "vfxopw", "Ice open water prod"     , "m/s"    ,   &
472      &       jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
473      CALL histdef( kid, "vfxsni", "Snow ice production "    , "m/s"    ,   &
474      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
475      CALL histdef( kid, "vfxres", "Ice prod from limupdate" , "m/s"    ,   &
476      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
477      CALL histdef( kid, "vfxbom", "Ice bottom melt"         , "m/s"    ,   &
478      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
479      CALL histdef( kid, "vfxsum", "Ice surface melt"        , "m/s"    ,   &
480      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
481
482      CALL histdef( kid, "sithicat", "Ice thickness"         , "m"      ,   &
483      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
484      CALL histdef( kid, "siconcat", "Ice concentration"     , "%"      ,   &
485      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
486      CALL histdef( kid, "sisalcat", "Ice salinity"           , ""      ,   &
487      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
488      CALL histdef( kid, "sitemcat", "Ice temperature"       , "C"      ,   &
489      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
490      CALL histdef( kid, "snthicat", "Snw thickness"         , "m"      ,   &
491      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
492      CALL histdef( kid, "sntemcat", "Snw temperature"       , "C"      ,   &
493      &      jpi, jpj, kh_i, jpl, 1, jpl, nz_i, 32, "inst(x)", rdt, rdt )
494
495      CALL histend( kid, snc4set )   ! end of the file definition
496
497      CALL histwrite( kid, "sithic", kt, htm_i         , jpi*jpj, (/1/) )   
498      CALL histwrite( kid, "siconc", kt, at_i          , jpi*jpj, (/1/) )
499      CALL histwrite( kid, "sitemp", kt, tm_i - rt0    , jpi*jpj, (/1/) )
500      CALL histwrite( kid, "sivelu", kt, u_ice          , jpi*jpj, (/1/) )
501      CALL histwrite( kid, "sivelv", kt, v_ice          , jpi*jpj, (/1/) )
502      CALL histwrite( kid, "sistru", kt, utau_ice       , jpi*jpj, (/1/) )
503      CALL histwrite( kid, "sistrv", kt, vtau_ice       , jpi*jpj, (/1/) )
504      CALL histwrite( kid, "sisflx", kt, qsr , jpi*jpj, (/1/) )
505      CALL histwrite( kid, "sinflx", kt, qns , jpi*jpj, (/1/) )
506      CALL histwrite( kid, "isnowpre", kt, sprecip        , jpi*jpj, (/1/) )
507      CALL histwrite( kid, "sisali", kt, smt_i          , jpi*jpj, (/1/) )
508      CALL histwrite( kid, "sivolu", kt, vt_i           , jpi*jpj, (/1/) )
509      CALL histwrite( kid, "sidive", kt, divu_i*1.0e8   , jpi*jpj, (/1/) )
510
511      ! MV MP 2016
512      IF ( ln_pnd ) THEN
513         CALL histwrite( kid, "si_amp", kt, at_ip         , jpi*jpj, (/1/) )
514         CALL histwrite( kid, "si_vmp", kt, vt_ip         , jpi*jpj, (/1/) )
515      ENDIF
516      ! END MV MP 2016
517
518      CALL histwrite( kid, "vfxbog", kt, wfx_bog        , jpi*jpj, (/1/) )
519      CALL histwrite( kid, "vfxdyn", kt, wfx_dyn        , jpi*jpj, (/1/) )
520      CALL histwrite( kid, "vfxopw", kt, wfx_opw        , jpi*jpj, (/1/) )
521      CALL histwrite( kid, "vfxsni", kt, wfx_sni        , jpi*jpj, (/1/) )
522      CALL histwrite( kid, "vfxres", kt, wfx_res        , jpi*jpj, (/1/) )
523      CALL histwrite( kid, "vfxbom", kt, wfx_bom        , jpi*jpj, (/1/) )
524      CALL histwrite( kid, "vfxsum", kt, wfx_sum        , jpi*jpj, (/1/) )
525      IF ( ln_pnd ) &
526         CALL histwrite( kid, "vfxpnd", kt, wfx_pnd     , jpi*jpj, (/1/) )
527
528      CALL histwrite( kid, "sithicat", kt, ht_i        , jpi*jpj*jpl, (/1/) )   
529      CALL histwrite( kid, "siconcat", kt, a_i         , jpi*jpj*jpl, (/1/) )   
530      CALL histwrite( kid, "sisalcat", kt, sm_i        , jpi*jpj*jpl, (/1/) )   
531      CALL histwrite( kid, "sitemcat", kt, tm_i - rt0  , jpi*jpj*jpl, (/1/) )   
532      CALL histwrite( kid, "snthicat", kt, ht_s        , jpi*jpj*jpl, (/1/) )   
533      CALL histwrite( kid, "sntemcat", kt, tm_su - rt0 , jpi*jpj*jpl, (/1/) )   
534
535      ! Close the file
536      ! -----------------
537      !CALL histclo( kid )
538
539    END SUBROUTINE lim_wri_state
540
541#else
542   !!----------------------------------------------------------------------
543   !!   Default option :         Empty module          NO LIM sea-ice model
544   !!----------------------------------------------------------------------
545CONTAINS
546   SUBROUTINE lim_wri          ! Empty routine
547   END SUBROUTINE lim_wri
548#endif
549
550   !!======================================================================
551END MODULE limwri
Note: See TracBrowser for help on using the repository browser.