source: NEMO/branches/2019/dev_r11943_MERGE_2019/src/SAS/diawri.F90 @ 12182

Last change on this file since 12182 was 12182, checked in by davestorkey, 20 months ago

2019/dev_r11943_MERGE_2019: Merge in dev_ASINTER-01-05_merge.

  • Property svn:keywords set to Id
File size: 24.3 KB
Line 
1MODULE diawri
2   !!======================================================================
3   !!                     ***  MODULE  diawri  ***
4   !! Ocean diagnostics :  write ocean output files
5   !!=====================================================================
6   !! History :  OPA  ! 1991-03  (M.-A. Foujols)  Original code
7   !!            4.0  ! 1991-11  (G. Madec)
8   !!                 ! 1992-06  (M. Imbard)  correction restart file
9   !!                 ! 1992-07  (M. Imbard)  split into diawri and rstwri
10   !!                 ! 1993-03  (M. Imbard)  suppress writibm
11   !!                 ! 1998-01  (C. Levy)  NETCDF format using ioipsl INTERFACE
12   !!                 ! 1999-02  (E. Guilyardi)  name of netCDF files + variables
13   !!            8.2  ! 2000-06  (M. Imbard)  Original code (diabort.F)
14   !!   NEMO     1.0  ! 2002-06  (A.Bozec, E. Durand)  Original code (diainit.F)
15   !!             -   ! 2002-09  (G. Madec)  F90: Free form and module
16   !!             -   ! 2002-12  (G. Madec)  merge of diabort and diainit, F90
17   !!                 ! 2005-11  (V. Garnier) Surface pressure gradient organization
18   !!            3.2  ! 2008-11  (B. Lemaire) creation from old diawri
19   !!----------------------------------------------------------------------
20
21   !!----------------------------------------------------------------------
22   !!   dia_wri       : create the standart output files
23   !!   dia_wri_state : create an output NetCDF file for a single instantaeous ocean state and forcing fields
24   !!----------------------------------------------------------------------
25   USE oce             ! ocean dynamics and tracers
26   USE abl            ! abl variables in case ln_abl = .true.
27   USE dom_oce         ! ocean space and time domain
28   USE zdf_oce         ! ocean vertical physics
29   USE sbc_oce         ! Surface boundary condition: ocean fields
30   USE sbc_ice         ! Surface boundary condition: ice fields
31   USE sbcssr          ! restoring term toward SST/SSS climatology
32   USE phycst          ! physical constants
33   USE zdfmxl          ! mixed layer
34   USE dianam          ! build name of file (routine)
35   USE zdfddm          ! vertical  physics: double diffusion
36   USE diahth          ! thermocline diagnostics
37   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
38   USE in_out_manager  ! I/O manager
39   USE iom
40   USE ioipsl
41#if defined key_si3
42   USE ice
43   USE icewri
44#endif
45   USE lib_mpp         ! MPP library
46   USE timing          ! preformance summary
47
48   IMPLICIT NONE
49   PRIVATE
50
51   PUBLIC   dia_wri                 ! routines called by step.F90
52   PUBLIC   dia_wri_state
53   PUBLIC   dia_wri_alloc           ! Called by nemogcm module
54#if ! defined key_iomput   
55   PUBLIC   dia_wri_alloc_abl       ! Called by sbcabl  module (if ln_abl = .true.)
56#endif
57   INTEGER ::   nid_T, nz_T, nh_T, ndim_T, ndim_hT   ! grid_T file
58   INTEGER ::   nid_U, nz_U, nh_U, ndim_U, ndim_hU   ! grid_U file
59   INTEGER ::   nid_V, nz_V, nh_V, ndim_V, ndim_hV   ! grid_V file
60   INTEGER ::   ndim_A, ndim_hA                      ! ABL file   
61   INTEGER ::   nid_A, nz_A, nh_A                    ! grid_ABL file   
62   INTEGER ::   ndex(1)                              ! ???
63   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV
64   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hA, ndex_A ! ABL
65
66   !! * Substitutions
67#  include "vectopt_loop_substitute.h90"
68   !!----------------------------------------------------------------------
69   !! NEMO/SAS 4.0 , NEMO Consortium (2018)
70   !! $Id$
71   !! Software governed by the CeCILL license (see ./LICENSE)
72   !!----------------------------------------------------------------------
73CONTAINS
74
75# if defined key_iomput
76   !!----------------------------------------------------------------------
77   !!   'key_iomput'                                        use IOM library
78   !!----------------------------------------------------------------------
79   INTEGER FUNCTION dia_wri_alloc()
80      !
81      dia_wri_alloc = 0
82      !
83   END FUNCTION dia_wri_alloc
84
85   
86   SUBROUTINE dia_wri( kt, Kmm )
87      !!---------------------------------------------------------------------
88      !!                  ***  ROUTINE dia_wri  ***
89      !!                   
90      !! ** Purpose :   Standard output of opa: dynamics and tracer fields
91      !!      NETCDF format is used by default
92      !!      Standalone surface scheme
93      !!
94      !! ** Method  :  use iom_put
95      !!----------------------------------------------------------------------
96      !!
97      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
98      INTEGER, INTENT( in ) ::   Kmm     ! ocean time levelindex
99      !!----------------------------------------------------------------------
100      !
101      ! Output the initial state and forcings
102      IF( ninist == 1 ) THEN
103         CALL dia_wri_state( 'output.init', Kmm )
104         ninist = 0
105      ENDIF
106      !
107   END SUBROUTINE dia_wri
108
109#else
110   !!----------------------------------------------------------------------
111   !!   Default option                                  use IOIPSL  library
112   !!----------------------------------------------------------------------
113   INTEGER FUNCTION dia_wri_alloc()
114      !!----------------------------------------------------------------------
115      INTEGER :: ierr
116      !!----------------------------------------------------------------------
117      !
118      ALLOCATE( ndex_hT(jpi*jpj), ndex_hU(jpi*jpj), ndex_hV(jpi*jpj), STAT=dia_wri_alloc )
119      CALL mpp_sum( 'diawri', dia_wri_alloc )
120      !
121   END FUNCTION dia_wri_alloc
122   
123   INTEGER FUNCTION dia_wri_alloc_abl()
124      !!----------------------------------------------------------------------
125     ALLOCATE(   ndex_hA(jpi*jpj), ndex_A (jpi*jpj*jpkam1), STAT=dia_wri_alloc_abl)
126      CALL mpp_sum( 'diawri', dia_wri_alloc_abl )
127      !
128   END FUNCTION dia_wri_alloc_abl
129 
130   SUBROUTINE dia_wri( kt )
131      !!---------------------------------------------------------------------
132      !!                  ***  ROUTINE dia_wri  ***
133      !!                   
134      !! ** Purpose :   Standard output of opa: dynamics and tracer fields
135      !!      NETCDF format is used by default
136      !!
137      !! ** Method  :   At the beginning of the first time step (nit000),
138      !!      define all the NETCDF files and fields
139      !!      At each time step call histdef to compute the mean if ncessary
140      !!      Each nn_write time step, output the instantaneous or mean fields
141      !!----------------------------------------------------------------------
142      !!
143      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
144      !!
145      LOGICAL ::   ll_print = .FALSE.                        ! =T print and flush numout
146      CHARACTER (len=40) ::   clhstnam, clop, clmx           ! local names
147      INTEGER  ::   inum = 11                                ! temporary logical unit
148      INTEGER  ::   ji, jj, jk                               ! dummy loop indices
149      INTEGER  ::   ierr                                     ! error code return from allocation
150      INTEGER  ::   iimi, iima, ipk, it, itmod, ijmi, ijma   ! local integers
151      INTEGER  ::   ipka                                     ! ABL
152      REAL(wp) ::   zsto, zout, zmax, zjulian                ! local scalars
153      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d_abl   ! ABL 3D workspace
154      !!----------------------------------------------------------------------
155      !
156      ! Output the initial state and forcings
157      IF( ninist == 1 ) THEN                       
158         CALL dia_wri_state( 'output.init' )
159         ninist = 0
160      ENDIF
161      !
162      IF( nn_write == -1 )   RETURN   ! we will never do any output
163      !
164      IF( ln_timing )   CALL timing_start('dia_wri')
165      !
166      ! 0. Initialisation
167      ! -----------------
168
169      ! local variable for debugging
170      ll_print = .FALSE.
171      ll_print = ll_print .AND. lwp
172
173      ! Define frequency of output and means
174      IF( ln_mskland )   THEN   ;   clop = "only(x)"   ! put 1.e+20 on land (very expensive!!)
175      ELSE                      ;   clop = "x"         ! no use of the mask value (require less cpu time)
176      ENDIF
177#if defined key_diainstant
178      zsto = nn_write * rdt
179      clop = "inst("//TRIM(clop)//")"
180#else
181      zsto=rdt
182      clop = "ave("//TRIM(clop)//")"
183#endif
184      zout = nn_write * rdt
185      zmax = ( nitend - nit000 + 1 ) * rdt
186
187      ! Define indices of the horizontal output zoom and vertical limit storage
188      iimi = 1      ;      iima = jpi
189      ijmi = 1      ;      ijma = jpj
190      ipk = jpk
191     IF(ln_abl) ipka = jpkam1
192
193      ! define time axis
194      it = kt
195      itmod = kt - nit000 + 1
196
197
198      ! 1. Define NETCDF files and fields at beginning of first time step
199      ! -----------------------------------------------------------------
200
201      IF( kt == nit000 ) THEN
202
203         ! Define the NETCDF files (one per grid)
204
205         ! Compute julian date from starting date of the run
206         CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
207         zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
208         IF(lwp)WRITE(numout,*)
209         IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear,   &
210            &                    ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
211         IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma,   &
212                                 ' limit storage in depth = ', ipk
213
214         ! WRITE root name in date.file for use by postpro
215         IF(lwp) THEN
216            CALL dia_nam( clhstnam, nn_write,' ' )
217            CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
218            WRITE(inum,*) clhstnam
219            CLOSE(inum)
220         ENDIF
221
222         ! Define the T grid FILE ( nid_T )
223
224         CALL dia_nam( clhstnam, nn_write, 'grid_T' )
225         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename
226         CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit
227            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
228            &          nit000-1, zjulian, rdt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set )
229         CALL histvert( nid_T, "deptht", "Vertical T levels",      &  ! Vertical grid: gdept
230            &           "m", ipk, gdept_1d, nz_T, "down" )
231         !                                                            ! Index of ocean points
232         CALL wheneq( jpi*jpj    , tmask, 1, 1., ndex_hT, ndim_hT )      ! surface
233
234         ! Define the U grid FILE ( nid_U )
235
236         CALL dia_nam( clhstnam, nn_write, 'grid_U' )
237         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename
238         CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu,           &  ! Horizontal grid: glamu and gphiu
239            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
240            &          nit000-1, zjulian, rdt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set )
241         CALL histvert( nid_U, "depthu", "Vertical U levels",      &  ! Vertical grid: gdept
242            &           "m", ipk, gdept_1d, nz_U, "down" )
243         !                                                            ! Index of ocean points
244         CALL wheneq( jpi*jpj    , umask, 1, 1., ndex_hU, ndim_hU )      ! surface
245
246         ! Define the V grid FILE ( nid_V )
247
248         CALL dia_nam( clhstnam, nn_write, 'grid_V' )                   ! filename
249         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
250         CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv,           &  ! Horizontal grid: glamv and gphiv
251            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
252            &          nit000-1, zjulian, rdt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set )
253         CALL histvert( nid_V, "depthv", "Vertical V levels",      &  ! Vertical grid : gdept
254            &          "m", ipk, gdept_1d, nz_V, "down" )
255         !                                                            ! Index of ocean points
256         CALL wheneq( jpi*jpj    , vmask, 1, 1., ndex_hV, ndim_hV )      ! surface
257
258         ! No W grid FILE
259         IF( ln_abl ) THEN 
260         ! Define the ABL grid FILE ( nid_A )
261            CALL dia_nam( clhstnam, nwrite, 'grid_ABL' )
262            IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename
263            CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit
264               &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       &
265               &          nit000-1, zjulian, rdt, nh_A, nid_A, domain_id=nidom, snc4chunks=snc4set )
266            CALL histvert( nid_A, "ght_abl", "Vertical T levels",      &  ! Vertical grid: gdept
267               &           "m", ipka, ght_abl(2:jpka), nz_A, "up" )
268            !                                                            ! Index of ocean points
269         ALLOCATE( zw3d_abl(jpi,jpj,ipka) ) 
270         zw3d_abl(:,:,:) = 1._wp 
271         CALL wheneq( jpi*jpj*ipka, zw3d_abl, 1, 1., ndex_A , ndim_A  )      ! volume
272            CALL wheneq( jpi*jpj     , zw3d_abl, 1, 1., ndex_hA, ndim_hA )      ! surface
273         DEALLOCATE(zw3d_abl)
274         ENDIF
275
276         ! Declare all the output fields as NETCDF variables
277
278         !                                                                                      !!! nid_T : 3D
279         CALL histdef( nid_T, "sst_m", "Sea Surface temperature"            , "C"      ,   &  ! sst
280            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
281         CALL histdef( nid_T, "sss_m", "Sea Surface Salinity"               , "PSU"    ,   &  ! sss
282            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
283         CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux"              , "Kg/m2/s",   &  ! (emp-rnf)
284            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
285         CALL histdef( nid_T, "sosfldow", "downward salt flux"                 , "PSU/m2/s",  &  ! (sfx)
286             &         jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
287         CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux"             , "W/m2"   ,   &  ! qns + qsr
288            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
289         CALL histdef( nid_T, "soshfldo", "Shortwave Radiation"                , "W/m2"   ,   &  ! qsr
290            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
291         CALL histdef( nid_T, "soicecov", "Ice fraction"                       , "[0,1]"  ,   &  ! fr_i
292            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
293         CALL histdef( nid_T, "sowindsp", "wind speed at 10m"                  , "m/s"    ,   &  ! wndm
294            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
295!
296         IF( ln_abl ) THEN
297         !                                                                                      !!! nid_A : 3D
298         CALL histdef( nid_A, "t_abl", "Potential Temperature"     , "K"        ,       &  ! t_abl
299               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
300            CALL histdef( nid_A, "q_abl", "Humidity"                  , "kg/kg"    ,       &  ! q_abl
301               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 
302            CALL histdef( nid_A, "u_abl", "Atmospheric U-wind   "     , "m/s"        ,     &  ! u_abl
303               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )
304            CALL histdef( nid_A, "v_abl", "Atmospheric V-wind   "     , "m/s"    ,         &  ! v_abl
305               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 
306            CALL histdef( nid_A, "tke_abl", "Atmospheric TKE   "     , "m2/s2"    ,        &  ! tke_abl
307               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 
308            CALL histdef( nid_A, "avm_abl", "Atmospheric turbulent viscosity", "m2/s"   ,  &  ! avm_abl
309               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 
310            CALL histdef( nid_A, "avt_abl", "Atmospheric turbulent diffusivity", "m2/s2",  &  ! avt_abl
311               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 
312            CALL histdef( nid_A, "pblh", "Atmospheric boundary layer height "  , "m",      &  ! pblh
313               &          jpi, jpj, nh_A,  1  , 1, 1   , -99 , 32, clop, zsto, zout )                 
314#if defined key_si3
315            CALL histdef( nid_A, "oce_frac", "Fraction of open ocean"  , " ",      &  ! ato_i
316               &          jpi, jpj, nh_A,  1  , 1, 1   , -99 , 32, clop, zsto, zout )
317#endif
318          CALL histend( nid_A, snc4chunks=snc4set )
319       !
320       ENDIF
321!
322
323         CALL histend( nid_T, snc4chunks=snc4set )
324
325         !                                                                                      !!! nid_U : 3D
326         CALL histdef( nid_U, "ssu_m", "Velocity component in x-direction", "m/s"   ,         &  ! ssu
327            &          jpi, jpj, nh_U, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
328         CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis"           , "N/m2"   ,   &  ! utau
329            &          jpi, jpj, nh_U, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
330
331         CALL histend( nid_U, snc4chunks=snc4set )
332
333         !                                                                                      !!! nid_V : 3D
334         CALL histdef( nid_V, "ssv_m", "Velocity component in y-direction", "m/s",            &  ! ssv_m
335            &          jpi, jpj, nh_V, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
336         CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis"           , "N/m2"   ,   &  ! vtau
337            &          jpi, jpj, nh_V, 1  , 1, 1  , - 99, 32, clop, zsto, zout )
338
339         CALL histend( nid_V, snc4chunks=snc4set )
340
341         IF(lwp) WRITE(numout,*)
342         IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
343         IF(ll_print) CALL FLUSH(numout )
344
345      ENDIF
346
347      ! 2. Start writing data
348      ! ---------------------
349
350      ! ndex(1) est utilise ssi l'avant dernier argument est diffferent de
351      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
352      ! donne le nombre d'elements, et ndex la liste des indices a sortir
353
354      IF( lwp .AND. MOD( itmod, nn_write ) == 0 ) THEN
355         WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
356         WRITE(numout,*) '~~~~~~ '
357      ENDIF
358
359      ! Write fields on T grid
360      CALL histwrite( nid_T, "sst_m", it, sst_m, ndim_hT, ndex_hT )   ! sea surface temperature
361      CALL histwrite( nid_T, "sss_m", it, sss_m, ndim_hT, ndex_hT )   ! sea surface salinity
362      CALL histwrite( nid_T, "sowaflup", it, (emp - rnf )  , ndim_hT, ndex_hT )   ! upward water flux
363      CALL histwrite( nid_T, "sosfldow", it, sfx           , ndim_hT, ndex_hT )   ! downward salt flux
364                                                                                  ! (includes virtual salt flux beneath ice
365                                                                                  ! in linear free surface case)
366
367      CALL histwrite( nid_T, "sohefldo", it, qns + qsr     , ndim_hT, ndex_hT )   ! total heat flux
368      CALL histwrite( nid_T, "soshfldo", it, qsr           , ndim_hT, ndex_hT )   ! solar heat flux
369      CALL histwrite( nid_T, "soicecov", it, fr_i          , ndim_hT, ndex_hT )   ! ice fraction   
370      CALL histwrite( nid_T, "sowindsp", it, wndm          , ndim_hT, ndex_hT )   ! wind speed   
371!
372      IF( ln_abl ) THEN
373        ALLOCATE( zw3d_abl(jpi,jpj,jpka) )
374        IF( ln_mskland )   THEN
375          DO jk=1,jpka
376             zw3d_abl(:,:,jk) = tmask(:,:,1)
377            END DO
378       ELSE
379            zw3d_abl(:,:,:) = 1._wp     
380         ENDIF       
381       CALL histwrite( nid_A,  "pblh"   , it, pblh(:,:)                  *zw3d_abl(:,:,1     ), ndim_hA, ndex_hA )   ! pblh
382        CALL histwrite( nid_A,  "u_abl"  , it, u_abl   (:,:,2:jpka,nt_n  )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! u_abl
383        CALL histwrite( nid_A,  "v_abl"  , it, v_abl   (:,:,2:jpka,nt_n  )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! v_abl
384        CALL histwrite( nid_A,  "t_abl"  , it, tq_abl  (:,:,2:jpka,nt_n,1)*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! t_abl
385        CALL histwrite( nid_A,  "q_abl"  , it, tq_abl  (:,:,2:jpka,nt_n,2)*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! q_abl     
386        CALL histwrite( nid_A,  "tke_abl", it, tke_abl (:,:,2:jpka,nt_n  )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! tke_abl
387        CALL histwrite( nid_A,  "avm_abl", it, avm_abl (:,:,2:jpka       )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! avm_abl
388        CALL histwrite( nid_A,  "avt_abl", it, avt_abl (:,:,2:jpka       )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! avt_abl 
389#if defined key_si3
390         CALL histwrite( nid_A,  "oce_frac"   , it, ato_i(:,:)                                  , ndim_hA, ndex_hA )   ! ato_i
391#endif
392       DEALLOCATE(zw3d_abl)
393     ENDIF
394!
395
396         ! Write fields on U grid
397      CALL histwrite( nid_U, "ssu_m"   , it, ssu_m         , ndim_hU, ndex_hU )   ! i-current speed
398      CALL histwrite( nid_U, "sozotaux", it, utau          , ndim_hU, ndex_hU )   ! i-wind stress
399
400         ! Write fields on V grid
401      CALL histwrite( nid_V, "ssv_m"   , it, ssv_m         , ndim_hV, ndex_hV )   ! j-current speed
402      CALL histwrite( nid_V, "sometauy", it, vtau          , ndim_hV, ndex_hV )   ! j-wind stress
403
404      ! 3. Close all files
405      ! ---------------------------------------
406      IF( kt == nitend ) THEN
407         CALL histclo( nid_T )
408         CALL histclo( nid_U )
409         CALL histclo( nid_V )
410         IF(ln_abl) CALL histclo( nid_A )
411      ENDIF
412      !
413      IF( ln_timing )   CALL timing_stop('dia_wri')
414      !
415   END SUBROUTINE dia_wri
416#endif
417
418   SUBROUTINE dia_wri_state( cdfile_name, Kmm )
419      !!---------------------------------------------------------------------
420      !!                 ***  ROUTINE dia_wri_state  ***
421      !!       
422      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
423      !!      the instantaneous ocean state and forcing fields.
424      !!        Used to find errors in the initial state or save the last
425      !!      ocean state in case of abnormal end of a simulation
426      !!
427      !! ** Method  :   NetCDF files using ioipsl
428      !!      File 'output.init.nc'  is created if ninist = 1 (namelist)
429      !!      File 'output.abort.nc' is created in case of abnormal job end
430      !!----------------------------------------------------------------------
431      CHARACTER (len=* ), INTENT( in ) ::   cdfile_name      ! name of the file created
432      INTEGER           , INTENT( in ) ::   Kmm              ! ocean time levelindex
433      !!
434      INTEGER :: inum
435      !!----------------------------------------------------------------------
436      !
437      IF(lwp) WRITE(numout,*)
438      IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
439      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created '
440      IF(lwp) WRITE(numout,*) '                and named :', cdfile_name, '...nc'
441
442#if defined key_si3
443     CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE., kdlev = jpl )
444#else
445     CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE. )
446#endif
447
448      CALL iom_rstput( 0, 0, inum, 'votemper', ts (:,:,:,jp_tem,Kmm) )    ! now temperature
449      CALL iom_rstput( 0, 0, inum, 'vosaline', ts (:,:,:,jp_sal,Kmm) )    ! now salinity
450      CALL iom_rstput( 0, 0, inum, 'sossheig', ssh(:,:,         Kmm) )    ! sea surface height
451      CALL iom_rstput( 0, 0, inum, 'vozocrtx', uu (:,:,:,       Kmm) )    ! now i-velocity
452      CALL iom_rstput( 0, 0, inum, 'vomecrty', vv (:,:,:,       Kmm) )    ! now j-velocity
453      CALL iom_rstput( 0, 0, inum, 'vovecrtz', ww                    )    ! now k-velocity
454      CALL iom_rstput( 0, 0, inum, 'sowaflup', emp - rnf             )    ! freshwater budget
455      CALL iom_rstput( 0, 0, inum, 'sohefldo', qsr + qns             )    ! total heat flux
456      CALL iom_rstput( 0, 0, inum, 'soshfldo', qsr                   )    ! solar heat flux
457      CALL iom_rstput( 0, 0, inum, 'soicecov', fr_i                  )    ! ice fraction
458      CALL iom_rstput( 0, 0, inum, 'sozotaux', utau                  )    ! i-wind stress
459      CALL iom_rstput( 0, 0, inum, 'sometauy', vtau                  )    ! j-wind stress
460 
461#if defined key_si3
462      IF( nn_ice == 2 ) THEN   ! condition needed in case agrif + ice-model but no-ice in child grid
463         CALL ice_wri_state( inum )
464      ENDIF
465#endif
466      !
467      CALL iom_close( inum )
468      !
469   END SUBROUTINE dia_wri_state
470
471   !!======================================================================
472END MODULE diawri
Note: See TracBrowser for help on using the repository browser.