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

source: branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90 @ 4525

Last change on this file since 4525 was 4525, checked in by flavoni, 9 years ago

change limwri.F90, use iom_put for all variables
  • Property svn:keywords set to Id
File size: 18.0 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 dom_ice
20   USE 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 par_ice
27   USE iom
28   USE timing          ! Timing
29   USE lib_fortran     ! Fortran utilities
30
31   IMPLICIT NONE
32   PRIVATE
33
34   PUBLIC lim_wri        ! routine called by lim_step.F90
35   PUBLIC lim_wri_state  ! called by dia_wri_state
36
37   REAL(wp)  ::   epsi06 = 1.e-6_wp
38   REAL(wp)  ::   zzero  = 0._wp
39   REAL(wp)  ::   zone   = 1._wp     
40   !!----------------------------------------------------------------------
41   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
42   !! $Id$
43   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
44   !!----------------------------------------------------------------------
45CONTAINS
46
47#if defined key_dimgout
48# include "limwri_dimg.h90"
49#else
50
51   SUBROUTINE lim_wri( kindic )
52      !!-------------------------------------------------------------------
53      !!  This routine computes the average of some variables and write it
54      !!  on the ouput files.
55      !!  ATTENTION cette routine n'est valable que si le pas de temps est
56      !!  egale a une fraction entiere de 1 jours.
57      !!  Diff 1-D 3-D : suppress common also included in etat
58      !!                 suppress cmoymo 11-18
59      !!  modif : 03/06/98
60      !!-------------------------------------------------------------------
61      INTEGER, INTENT(in) ::   kindic   ! if kindic < 0 there has been an error somewhere
62      !
63      INTEGER ::  ji, jj, jk, jl  ! dummy loop indices
64      REAL(wp) ::  zinda, zindb, z1_365
65      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zoi, zei
66      REAL(wp), POINTER, DIMENSION(:,:)   :: z2d, z2da, z2db, zind    ! 2D workspace
67      !!-------------------------------------------------------------------
68
69      IF( nn_timing == 1 )  CALL timing_start('limwri')
70
71      CALL wrk_alloc( jpi, jpj, jpl, zoi, zei )
72      CALL wrk_alloc( jpi, jpj     , z2d, z2da, z2db, zind )
73
74
75      !-----------------------------
76      ! Mean category values
77      !-----------------------------
78
79      CALL lim_var_icetm      ! mean sea ice temperature
80
81      CALL lim_var_bv         ! brine volume
82
83      DO jj = 1, jpj          ! presence indicator of ice
84         DO ji = 1, jpi
85            zind(ji,jj)  = MAX( zzero , SIGN( zone , at_i(ji,jj) - epsi06 ) )
86         END DO
87      END DO
88      !
89      !
90      CALL iom_put( "iceconc"      , at_i               )        ! ice concentration
91      !                                             
92      DO jj = 1, jpj                                             ! mean ice thickness
93         DO ji = 1, jpi
94            z2d(ji,jj)  = vt_i(ji,jj) / MAX( at_i(ji,jj), epsi06 ) * zind(ji,jj)
95         END DO
96      END DO
97      CALL iom_put( "icethic_cea"  , z2d                 )       ! ice thickness (i.e. icethi(:,:))
98      CALL iom_put( "icevolu"      , vt_i                )       ! ice volume = mean ice thickness over the cell
99      DO jj = 1, jpj                                           
100         DO ji = 1, jpi
101            z2d(ji,jj)  = vt_s(ji,jj) / MAX( at_i(ji,jj), epsi06 ) * zind(ji,jj)
102         END DO
103      END DO
104      CALL iom_put( "snowthic_cea" , z2d                 )       ! snow thickness = mean snow thickness over the cell
105      CALL iom_put( "ioceflxb"     , fbif * at_i         )       ! oceanic flux at the ice base
106      CALL iom_put( "isst"         , sst_m               )       ! sea surface temperature
107      CALL iom_put( "isss"         , sss_m               )       ! sea surface salinity
108      CALL iom_put( "qt_oce"       , qns + qsr           )       ! total flux at ocean surface
109      !
110      DO jj = 2 , jpjm1
111         DO ji = 2 , jpim1
112            z2da(ji,jj)  = (  u_ice(ji,jj) * tmu(ji,jj) + u_ice(ji-1,jj) * tmu(ji-1,jj) ) * 0.5_wp
113            z2db(ji,jj)  = (  v_ice(ji,jj) * tmv(ji,jj) + v_ice(ji,jj-1) * tmv(ji,jj-1) ) * 0.5_wp
114        END DO
115      END DO
116      CALL lbc_lnk( z2da, 'T', -1. )
117      CALL lbc_lnk( z2db, 'T', -1. )
118      DO jj = 1, jpj                                 
119         DO ji = 1, jpi
120            z2d(ji,jj)  = SQRT( z2da(ji,jj) * z2da(ji,jj) + z2db(ji,jj) * z2db(ji,jj) ) 
121         END DO
122      END DO
123      CALL iom_put( "uice_ipa"     , z2da                )       ! ice velocity u component
124      CALL iom_put( "vice_ipa"     , z2db                )       ! ice velocity v component
125      CALL iom_put( "icevel"       , z2d                 )       ! ice velocity module
126      CALL iom_put( "icebopr"      , diag_bot_gr * rday  )       ! daily bottom thermodynamic ice production
127      CALL iom_put( "icedypr"      , diag_dyn_gr * rday  )       ! daily dynamic ice production (rid/raft)
128!!SF BE CAREFUL : qsr_oce qnd qns_oce are after penetration over ice
129      CALL iom_put( "qsr_oce"      , qsr                 )       ! solar flux at ocean surface
130      CALL iom_put( "qns_oce"      , qns                 )       ! non-solar flux at ocean surface
131!!SF end be careful
132      CALL iom_put( "hfbri"        , fhbri               )       ! heat flux due to brine release
133      CALL iom_put( "utau_ice"     , utau_ice            )       ! wind stress over ice along i-axis at I-point
134      CALL iom_put( "vtau_ice"     , vtau_ice            )       ! wind stress over ice along j-axis at I-point
135!!SF commented because this computation is not ok
136 !SF because qsr is not qsr_ocean but it contains already qsr_ice
137      !SF
138      !SF DO jj = 1 , jpj
139      !SF    DO ji = 1 , jpi
140      !SF       z2d(ji,jj)  = ( 1._wp - at_i(ji,jj) ) * qsr(ji,jj)
141      !SF   END DO
142      !SF END DO
143      !SF CALL iom_put( "qsr_io"       , z2d                 )        ! solar flux at ice/ocean surface
144      !SF DO jj = 1 , jpj
145      !SF    DO ji = 1 , jpi
146      !SF       z2d(ji,jj)  = ( 1._wp - at_i(ji,jj) ) * qns(ji,jj)
147      !SF   END DO
148      !SF END DO
149      !SF CALL iom_put( "qns_io"       , z2d                 )        ! non-solar flux at ice/ocean surface
150      CALL iom_put( "snowpre"      , sprecip             )        ! snow precipitation
151      CALL iom_put( "micesalt"     , smt_i               )        ! mean ice salinity
152      !
153      z2d(:,:) = 0.e0
154      DO jl = 1, jpl
155         DO jj = 1, jpj
156            DO ji = 1, jpi
157               z2d(ji,jj) = z2d(ji,jj) + zind(ji,jj) * oa_i(ji,jj,jl)
158            END DO
159         END DO
160      END DO
161      z1_365 = 1._wp / 365._wp
162      CALL iom_put( "miceage"     , z2d * z1_365         )        ! mean ice age
163      CALL iom_put( "icelapr"     , diag_lat_gr * rday   )        ! daily lateral thermodynamic ice production
164      CALL iom_put( "icesipr"     , diag_sni_gr * rday   )        ! daily snowice ice production
165      DO jj = 1, jpj
166         DO ji = 1, jpi
167            z2d(ji,jj) = ( tm_i(ji,jj) - rtt ) * zind(ji,jj)
168         END DO
169      END DO
170
171      CALL iom_put( "micet"       , z2d                  )        ! mean ice temperature
172      CALL iom_put( "icehc"       , et_i                 )        ! ice total heat content
173      CALL iom_put( "isnowhc"     , et_s                 )        ! snow total heat content
174      !
175      z2d(:,:) = 0.e0
176      DO jl = 1, jpl
177         DO jj = 1, jpj
178            DO ji = 1, jpi
179               z2d(ji,jj) = z2d(ji,jj) + zind(ji,jj) * ( t_su(ji,jj,jl) - rtt ) * a_i(ji,jj,jl) / MAX( at_i(ji,jj) , epsi06 )
180            END DO
181         END DO
182      END DO
183      CALL iom_put( "icest"       , z2d                 )        ! ice surface temperature
184      CALL iom_put( "sfxthd"      , sfx_thd * rday      )        ! equivalent FW salt flux
185      CALL iom_put( "ibrinv"      , bv_i * 100._wp      )        ! brine volume
186      DO jj = 1, jpj
187         DO ji = 1, jpi
188            zindb  = MAX( zzero , SIGN( zone , at_i(ji,jj) ) )
189            z2d(ji,jj) = hicol(ji,jj) * zindb
190         END DO
191      END DO
192      CALL iom_put( "icecolf"     , z2d                 )        ! frazil ice collection thickness
193      CALL iom_put( "icestr"      , strength * 0.001    )        ! ice strength
194      CALL iom_put( "isume"       , diag_sur_me * rday  )        ! surface melt
195      CALL iom_put( "ibome"       , diag_bot_me * rday  )        ! bottom melt
196      CALL iom_put( "idive"       , divu_i * 1.0e8      )        ! divergence
197      CALL iom_put( "ishear"      , shear_i * 1.0e8     )        ! shear
198      CALL iom_put( "isume"       , diag_sur_me * rday  )        ! surface melt
199      CALL iom_put( "icerepr"     , diag_res_pr * rday  )        ! daily prod./melting due to limupdate
200      CALL iom_put( "snowvol"     , vt_s                )        ! snow volume
201      CALL iom_put( "sfxmec"      , sfx_mec * rday      )        ! salt flux from ridging rafting
202      CALL iom_put( "sfxres"      , sfx_res * rday      )        ! salt flux from limupdate (resultant)
203      CALL iom_put( "icetrp"      , diag_trp_vi * rday  )        ! ice volume transport
204
205
206      !--------------------------------
207      ! Output values for each category
208      !--------------------------------
209
210         DO jl = 1, jpl 
211            CALL lbc_lnk( a_i(:,:,jl)  , 'T' ,  1. )
212            CALL lbc_lnk( sm_i(:,:,jl) , 'T' ,  1. )
213            CALL lbc_lnk( oa_i(:,:,jl) , 'T' ,  1. )
214            CALL lbc_lnk( ht_i(:,:,jl) , 'T' ,  1. )
215            CALL lbc_lnk( ht_s(:,:,jl) , 'T' ,  1. )
216         END DO
217
218         DO jl = 1, jpl 
219            DO jj = 1, jpj
220               DO ji = 1, jpi
221                  zinda = MAX( zzero , SIGN( zone , a_i(ji,jj,jl) - epsi06 ) )
222                  zoi(ji,jj,jl) = oa_i(ji,jj,jl)  / MAX( a_i(ji,jj,jl) , epsi06 ) * zinda
223               END DO
224            END DO
225         END DO
226 
227         CALL iom_put( "iceage_cat"     , zoi         )        ! ice age for categories
228
229         ! Compute brine volume
230         zei(:,:,:) = 0._wp
231         DO jl = 1, jpl 
232            DO jk = 1, nlay_i
233               DO jj = 1, jpj
234                  DO ji = 1, jpi
235                     zinda = MAX( zzero , SIGN( zone , a_i(ji,jj,jl) - epsi06 ) )
236                     zei(ji,jj,jl) = zei(ji,jj,jl) + 100.0* &
237                        ( - tmut * s_i(ji,jj,jk,jl) / MIN( ( t_i(ji,jj,jk,jl) - rtt ), - epsi06 ) ) * &
238                        zinda / nlay_i
239                  END DO
240               END DO
241            END DO
242         END DO
243
244         DO jl = 1, jpl 
245            CALL lbc_lnk( zei(:,:,jl) , 'T' ,  1. )
246         END DO
247
248         CALL iom_put( "iceconc_cat"      , a_i         )        ! area for categories
249         CALL iom_put( "icethic_cat"      , ht_i        )        ! thickness for categories
250         CALL iom_put( "snowthic_cat"     , ht_s        )        ! snow depth for categories
251         CALL iom_put( "salinity_cat"     , sm_i        )        ! salinity for categories
252         CALL iom_put( "brinevol_cat"     , zei         )        ! brine volume for categories
253
254      !     !  Create an output files (output.lim.abort.nc) if S < 0 or u > 20
255      !     m/s
256      !     IF( kindic < 0 )   CALL lim_wri_state( 'output.abort' )
257      !     not yet implemented
258
259      CALL wrk_dealloc( jpi, jpj, jpl, zoi, zei )
260      CALL wrk_dealloc( jpi, jpj     , z2d, zind, z2da, z2db )
261
262      IF( nn_timing == 1 )  CALL timing_stop('limwri')
263     
264   END SUBROUTINE lim_wri
265#endif
266
267 
268   SUBROUTINE lim_wri_state( kt, kid, kh_i )
269      !!---------------------------------------------------------------------
270      !!                 ***  ROUTINE lim_wri_state  ***
271      !!       
272      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
273      !!      the instantaneous ice state and forcing fields for ice model
274      !!        Used to find errors in the initial state or save the last
275      !!      ocean state in case of abnormal end of a simulation
276      !!
277      !! History :
278      !!   4.1  !  2013-06  (C. Rousset)
279      !!----------------------------------------------------------------------
280      INTEGER, INTENT( in ) ::   kt               ! ocean time-step index)
281      INTEGER, INTENT( in ) ::   kid , kh_i       
282      !!----------------------------------------------------------------------
283   
284
285   
286      CALL histdef( kid, "iicethic", "Ice thickness"           , "m"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
287      CALL histdef( kid, "iiceconc", "Ice concentration"       , "%"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
288      CALL histdef( kid, "iicetemp", "Ice temperature"         , "C"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
289      CALL histdef( kid, "iicevelu", "i-Ice speed (I-point)"   , "m/s"    , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
290      CALL histdef( kid, "iicevelv", "j-Ice speed (I-point)"   , "m/s"    , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
291      CALL histdef( kid, "iicestru", "i-Wind stress over ice (I-pt)", "Pa", jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
292      CALL histdef( kid, "iicestrv", "j-Wind stress over ice (I-pt)", "Pa", jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
293      CALL histdef( kid, "iicesflx", "Solar flux over ocean"     , "w/m2"   , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
294      CALL histdef( kid, "iicenflx", "Non-solar flux over ocean" , "w/m2"   , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
295      CALL histdef( kid, "isnowpre", "Snow precipitation"      , "kg/m2/s", jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
296      CALL histdef( kid, "iicesali", "Ice salinity"            , "PSU"    , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
297      CALL histdef( kid, "iicevolu", "Ice volume"              , "m"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
298      CALL histdef( kid, "iicedive", "Ice divergence"          , "10-8s-1", jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
299      CALL histdef( kid, "iicebopr", "Ice bottom production"   , "m/s"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
300      CALL histdef( kid, "iicedypr", "Ice dynamic production"  , "m/s"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
301      CALL histdef( kid, "iicelapr", "Ice open water prod"     , "m/s"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
302      CALL histdef( kid, "iicesipr", "Snow ice production "    , "m/s"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
303      CALL histdef( kid, "iicerepr", "Ice prod from limupdate" , "m/s"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
304      CALL histdef( kid, "iicebome", "Ice bottom melt"         , "m/s"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
305      CALL histdef( kid, "iicesume", "Ice surface melt"        , "m/s"      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
306      CALL histdef( kid, "iisfxthd", "Salt flux from thermo"   , ""      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
307      CALL histdef( kid, "iisfxmec", "Salt flux from dynmics"  , ""      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
308      CALL histdef( kid, "iisfxres", "Salt flux from limupdate", ""      , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
309
310      CALL histend( kid, snc4set )   ! end of the file definition
311
312      CALL histwrite( kid, "iicethic", kt, icethi        , jpi*jpj, (/1/) )   
313      CALL histwrite( kid, "iiceconc", kt, at_i          , jpi*jpj, (/1/) )
314      CALL histwrite( kid, "iicetemp", kt, tm_i - rtt    , jpi*jpj, (/1/) )
315      CALL histwrite( kid, "iicevelu", kt, u_ice          , jpi*jpj, (/1/) )
316      CALL histwrite( kid, "iicevelv", kt, v_ice          , jpi*jpj, (/1/) )
317      CALL histwrite( kid, "iicestru", kt, utau_ice       , jpi*jpj, (/1/) )
318      CALL histwrite( kid, "iicestrv", kt, vtau_ice       , jpi*jpj, (/1/) )
319      CALL histwrite( kid, "iicesflx", kt, qsr , jpi*jpj, (/1/) )
320      CALL histwrite( kid, "iicenflx", kt, qns , jpi*jpj, (/1/) )
321      CALL histwrite( kid, "isnowpre", kt, sprecip        , jpi*jpj, (/1/) )
322      CALL histwrite( kid, "iicesali", kt, smt_i          , jpi*jpj, (/1/) )
323      CALL histwrite( kid, "iicevolu", kt, vt_i           , jpi*jpj, (/1/) )
324      CALL histwrite( kid, "iicedive", kt, divu_i*1.0e8   , jpi*jpj, (/1/) )
325
326      CALL histwrite( kid, "iicebopr", kt, diag_bot_gr        , jpi*jpj, (/1/) )
327      CALL histwrite( kid, "iicedypr", kt, diag_dyn_gr        , jpi*jpj, (/1/) )
328      CALL histwrite( kid, "iicelapr", kt, diag_lat_gr        , jpi*jpj, (/1/) )
329      CALL histwrite( kid, "iicesipr", kt, diag_sni_gr        , jpi*jpj, (/1/) )
330      CALL histwrite( kid, "iicerepr", kt, diag_res_pr        , jpi*jpj, (/1/) )
331      CALL histwrite( kid, "iicebome", kt, diag_bot_me        , jpi*jpj, (/1/) )
332      CALL histwrite( kid, "iicesume", kt, diag_sur_me        , jpi*jpj, (/1/) )
333      CALL histwrite( kid, "iisfxthd", kt, sfx_thd        , jpi*jpj, (/1/) )
334      CALL histwrite( kid, "iisfxmec", kt, sfx_mec        , jpi*jpj, (/1/) )
335      CALL histwrite( kid, "iisfxres", kt, sfx_res        , jpi*jpj, (/1/) )
336
337      ! Close the file
338      ! -----------------
339      CALL histclo( kid )
340
341    END SUBROUTINE lim_wri_state
342
343#else
344   !!----------------------------------------------------------------------
345   !!   Default option :         Empty module          NO LIM sea-ice model
346   !!----------------------------------------------------------------------
347CONTAINS
348   SUBROUTINE lim_wri          ! Empty routine
349   END SUBROUTINE lim_wri
350#endif
351
352   !!======================================================================
353END MODULE limwri
Note: See TracBrowser for help on using the repository browser.