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.
diawri.F90 in NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/BENCH/MY_SRC – NEMO

source: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/BENCH/MY_SRC/diawri.F90 @ 10179

Last change on this file since 10179 was 10179, checked in by smasson, 6 years ago

dev_r10164_HPC09_ESIWACE_PREP_MERGE: action 4a: add si3 in BENCH, see #2133
File size: 27.3 KB
RevLine 
[9762]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   !!            3.7  ! 2014-01  (G. Madec) remove eddy induced velocity from no-IOM output
20   !!                 !                     change name of output variables in dia_wri_state
21   !!----------------------------------------------------------------------
22
23   !!----------------------------------------------------------------------
24   !!   dia_wri       : create the standart output files
25   !!   dia_wri_state : create an output NetCDF file for a single instantaeous ocean state and forcing fields
26   !!----------------------------------------------------------------------
27   USE oce            ! ocean dynamics and tracers
28   USE dom_oce        ! ocean space and time domain
29   USE phycst         ! physical constants
30   USE dianam         ! build name of file (routine)
31   USE diahth         ! thermocline diagnostics
32   USE dynadv   , ONLY: ln_dynadv_vec
33   USE icb_oce        ! Icebergs
34   USE icbdia         ! Iceberg budgets
35   USE ldftra         ! lateral physics: eddy diffusivity coef.
36   USE ldfdyn         ! lateral physics: eddy viscosity   coef.
37   USE sbc_oce        ! Surface boundary condition: ocean fields
38   USE sbc_ice        ! Surface boundary condition: ice fields
39   USE sbcssr         ! restoring term toward SST/SSS climatology
40   USE sbcwave        ! wave parameters
41   USE wet_dry        ! wetting and drying
42   USE zdf_oce        ! ocean vertical physics
43   USE zdfdrg         ! ocean vertical physics: top/bottom friction
44   USE zdfmxl         ! mixed layer
45   !
46   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
47   USE in_out_manager ! I/O manager
48   USE diatmb         ! Top,middle,bottom output
49   USE dia25h         ! 25h Mean output
50   USE iom            !
51   USE ioipsl         !
52
53#if defined key_si3
54   USE icewri 
55#endif
56   USE lib_mpp         ! MPP library
57   USE timing          ! preformance summary
58   USE diurnal_bulk    ! diurnal warm layer
59   USE cool_skin       ! Cool skin
60
61   IMPLICIT NONE
62   PRIVATE
63
64   PUBLIC   dia_wri                 ! routines called by step.F90
65   PUBLIC   dia_wri_state
66   PUBLIC   dia_wri_alloc           ! Called by nemogcm module
67
68   INTEGER ::   nid_T, nz_T, nh_T, ndim_T, ndim_hT   ! grid_T file
69   INTEGER ::          nb_T              , ndim_bT   ! grid_T file
70   INTEGER ::   nid_U, nz_U, nh_U, ndim_U, ndim_hU   ! grid_U file
71   INTEGER ::   nid_V, nz_V, nh_V, ndim_V, ndim_hV   ! grid_V file
72   INTEGER ::   nid_W, nz_W, nh_W                    ! grid_W file
73   INTEGER ::   ndex(1)                              ! ???
74   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV
75   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_T, ndex_U, ndex_V
76   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_bT
77
78   !! * Substitutions
79#  include "vectopt_loop_substitute.h90"
80   !!----------------------------------------------------------------------
81   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
82   !! $Id: diawri.F90 9598 2018-05-15 22:47:16Z nicolasmartin $
83   !! Software governed by the CeCILL licence     (./LICENSE)
84   !!----------------------------------------------------------------------
85CONTAINS
86
87#if defined key_iomput
88   !!----------------------------------------------------------------------
89   !!   'key_iomput'                                        use IOM library
90   !!----------------------------------------------------------------------
91   INTEGER FUNCTION dia_wri_alloc()
92      !
93      dia_wri_alloc = 0
94      !
95   END FUNCTION dia_wri_alloc
96
97   
98   SUBROUTINE dia_wri( kt )
99      !!---------------------------------------------------------------------
100      !!                  ***  ROUTINE dia_wri  ***
101      !!                   
102      !! ** Purpose :   Standard output of opa: dynamics and tracer fields
103      !!      NETCDF format is used by default
104      !!
105      !! ** Method  :  use iom_put
106      !!----------------------------------------------------------------------
107      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
108      !!
109      INTEGER ::   ji, jj, jk       ! dummy loop indices
110      INTEGER ::   ikbot            ! local integer
111      REAL(wp)::   zztmp , zztmpx   ! local scalar
112      REAL(wp)::   zztmp2, zztmpy   !   -      -
113      REAL(wp), DIMENSION(jpi,jpj)     ::   z2d   ! 2D workspace
114      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   z3d   ! 3D workspace
115      !!----------------------------------------------------------------------
116      !
117      IF( ln_timing )   CALL timing_start('dia_wri')
118      !
119      ! Output the initial state and forcings
120      IF( ninist == 1 ) THEN                       
121         CALL dia_wri_state( 'output.init', kt )
122         ninist = 0
123      ENDIF
124
125      ! Output of initial vertical scale factor
126      CALL iom_put("e3t_0", e3t_0(:,:,:) )
[10166]127      CALL iom_put("e3u_0", e3u_0(:,:,:) )
128      CALL iom_put("e3v_0", e3v_0(:,:,:) )
[9762]129      !
130      CALL iom_put( "e3t" , e3t_n(:,:,:) )
131      CALL iom_put( "e3u" , e3u_n(:,:,:) )
132      CALL iom_put( "e3v" , e3v_n(:,:,:) )
133      CALL iom_put( "e3w" , e3w_n(:,:,:) )
134      IF( iom_use("e3tdef") )   &
135         CALL iom_put( "e3tdef"  , ( ( e3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 )
136
137      IF( ll_wd ) THEN
138         CALL iom_put( "ssh" , (sshn+ssh_ref)*tmask(:,:,1) )   ! sea surface height (brought back to the reference used for wetting and drying)
139      ELSE
140         CALL iom_put( "ssh" , sshn )              ! sea surface height
141      ENDIF
142
143      IF( iom_use("wetdep") )   &                  ! wet depth
144         CALL iom_put( "wetdep" , ht_0(:,:) + sshn(:,:) )
145     
146      CALL iom_put( "toce", tsn(:,:,:,jp_tem) )    ! 3D temperature
147      CALL iom_put(  "sst", tsn(:,:,1,jp_tem) )    ! surface temperature
148      IF ( iom_use("sbt") ) THEN
149         DO jj = 1, jpj
150            DO ji = 1, jpi
151               ikbot = mbkt(ji,jj)
152               z2d(ji,jj) = tsn(ji,jj,ikbot,jp_tem)
153            END DO
154         END DO
155         CALL iom_put( "sbt", z2d )                ! bottom temperature
156      ENDIF
157     
158      CALL iom_put( "soce", tsn(:,:,:,jp_sal) )    ! 3D salinity
159      CALL iom_put(  "sss", tsn(:,:,1,jp_sal) )    ! surface salinity
160      IF ( iom_use("sbs") ) THEN
161         DO jj = 1, jpj
162            DO ji = 1, jpi
163               ikbot = mbkt(ji,jj)
164               z2d(ji,jj) = tsn(ji,jj,ikbot,jp_sal)
165            END DO
166         END DO
167         CALL iom_put( "sbs", z2d )                ! bottom salinity
168      ENDIF
169
170      IF ( iom_use("taubot") ) THEN                ! bottom stress
171         zztmp = rau0 * 0.25
172         z2d(:,:) = 0._wp
173         DO jj = 2, jpjm1
174            DO ji = fs_2, fs_jpim1   ! vector opt.
175               zztmp2 = (  ( rCdU_bot(ji+1,jj)+rCdU_bot(ji  ,jj) ) * un(ji  ,jj,mbku(ji  ,jj))  )**2   &
176                  &   + (  ( rCdU_bot(ji  ,jj)+rCdU_bot(ji-1,jj) ) * un(ji-1,jj,mbku(ji-1,jj))  )**2   &
177                  &   + (  ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj  ) ) * vn(ji,jj  ,mbkv(ji,jj  ))  )**2   &
178                  &   + (  ( rCdU_bot(ji,jj  )+rCdU_bot(ji,jj-1) ) * vn(ji,jj-1,mbkv(ji,jj-1))  )**2
179               z2d(ji,jj) = zztmp * SQRT( zztmp2 ) * tmask(ji,jj,1) 
180               !
181            END DO
182         END DO
[10170]183         CALL lbc_lnk( 'diawri', z2d, 'T', 1. )
[9762]184         CALL iom_put( "taubot", z2d )           
185      ENDIF
186         
187      CALL iom_put( "uoce", un(:,:,:) )            ! 3D i-current
188      CALL iom_put(  "ssu", un(:,:,1) )            ! surface i-current
189      IF ( iom_use("sbu") ) THEN
190         DO jj = 1, jpj
191            DO ji = 1, jpi
192               ikbot = mbku(ji,jj)
193               z2d(ji,jj) = un(ji,jj,ikbot)
194            END DO
195         END DO
196         CALL iom_put( "sbu", z2d )                ! bottom i-current
197      ENDIF
198     
199      CALL iom_put( "voce", vn(:,:,:) )            ! 3D j-current
200      CALL iom_put(  "ssv", vn(:,:,1) )            ! surface j-current
201      IF ( iom_use("sbv") ) THEN
202         DO jj = 1, jpj
203            DO ji = 1, jpi
204               ikbot = mbkv(ji,jj)
205               z2d(ji,jj) = vn(ji,jj,ikbot)
206            END DO
207         END DO
208         CALL iom_put( "sbv", z2d )                ! bottom j-current
209      ENDIF
210
211      CALL iom_put( "woce", wn )                   ! vertical velocity
212      IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN   ! vertical mass transport & its square value
213         ! Caution: in the VVL case, it only correponds to the baroclinic mass transport.
214         z2d(:,:) = rau0 * e1e2t(:,:)
215         DO jk = 1, jpk
216            z3d(:,:,jk) = wn(:,:,jk) * z2d(:,:)
217         END DO
218         CALL iom_put( "w_masstr" , z3d ) 
219         IF( iom_use('w_masstr2') )   CALL iom_put( "w_masstr2", z3d(:,:,:) * z3d(:,:,:) )
220      ENDIF
221
222      CALL iom_put( "avt" , avt )                  ! T vert. eddy diff. coef.
223      CALL iom_put( "avs" , avs )                  ! S vert. eddy diff. coef.
224      CALL iom_put( "avm" , avm )                  ! T vert. eddy visc. coef.
225
226      IF( iom_use('logavt') )   CALL iom_put( "logavt", LOG( MAX( 1.e-20_wp, avt(:,:,:) ) ) )
227      IF( iom_use('logavs') )   CALL iom_put( "logavs", LOG( MAX( 1.e-20_wp, avs(:,:,:) ) ) )
228
229      IF ( iom_use("sstgrad") .OR. iom_use("sstgrad2") ) THEN
230         DO jj = 2, jpjm1                                    ! sst gradient
231            DO ji = fs_2, fs_jpim1   ! vector opt.
232               zztmp  = tsn(ji,jj,1,jp_tem)
233               zztmpx = ( tsn(ji+1,jj,1,jp_tem) - zztmp ) * r1_e1u(ji,jj) + ( zztmp - tsn(ji-1,jj  ,1,jp_tem) ) * r1_e1u(ji-1,jj)
234               zztmpy = ( tsn(ji,jj+1,1,jp_tem) - zztmp ) * r1_e2v(ji,jj) + ( zztmp - tsn(ji  ,jj-1,1,jp_tem) ) * r1_e2v(ji,jj-1)
235               z2d(ji,jj) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy )   &
236                  &              * umask(ji,jj,1) * umask(ji-1,jj,1) * vmask(ji,jj,1) * umask(ji,jj-1,1)
237            END DO
238         END DO
[10170]239         CALL lbc_lnk( 'diawri', z2d, 'T', 1. )
[9762]240         CALL iom_put( "sstgrad2",  z2d )          ! square of module of sst gradient
241         z2d(:,:) = SQRT( z2d(:,:) )
242         CALL iom_put( "sstgrad" ,  z2d )          ! module of sst gradient
243      ENDIF
244         
245      ! heat and salt contents
246      IF( iom_use("heatc") ) THEN
247         z2d(:,:)  = 0._wp 
248         DO jk = 1, jpkm1
249            DO jj = 1, jpj
250               DO ji = 1, jpi
251                  z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk)
252               END DO
253            END DO
254         END DO
255         CALL iom_put( "heatc", rau0_rcp * z2d )   ! vertically integrated heat content (J/m2)
256      ENDIF
257
258      IF( iom_use("saltc") ) THEN
259         z2d(:,:)  = 0._wp 
260         DO jk = 1, jpkm1
261            DO jj = 1, jpj
262               DO ji = 1, jpi
263                  z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_sal) * tmask(ji,jj,jk)
264               END DO
265            END DO
266         END DO
267         CALL iom_put( "saltc", rau0 * z2d )          ! vertically integrated salt content (PSU*kg/m2)
268      ENDIF
269      !
270      IF ( iom_use("eken") ) THEN
271         z3d(:,:,jpk) = 0._wp 
272         DO jk = 1, jpkm1
273            DO jj = 2, jpjm1
274               DO ji = fs_2, fs_jpim1   ! vector opt.
275                  zztmp  = 0.25_wp * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk)
276                  z3d(ji,jj,jk) = zztmp * (  un(ji-1,jj,jk)**2 * e2u(ji-1,jj) * e3u_n(ji-1,jj,jk)   &
277                     &                     + un(ji  ,jj,jk)**2 * e2u(ji  ,jj) * e3u_n(ji  ,jj,jk)   &
278                     &                     + vn(ji,jj-1,jk)**2 * e1v(ji,jj-1) * e3v_n(ji,jj-1,jk)   &
279                     &                     + vn(ji,jj  ,jk)**2 * e1v(ji,jj  ) * e3v_n(ji,jj  ,jk)   )
280               END DO
281            END DO
282         END DO
[10170]283         CALL lbc_lnk( 'diawri', z3d, 'T', 1. )
[9762]284         CALL iom_put( "eken", z3d )                 ! kinetic energy
285      ENDIF
286      !
287      CALL iom_put( "hdiv", hdivn )                  ! Horizontal divergence
288      !
289      IF( iom_use("u_masstr") .OR. iom_use("u_masstr_vint") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN
290         z3d(:,:,jpk) = 0.e0
291         z2d(:,:) = 0.e0
292         DO jk = 1, jpkm1
293            z3d(:,:,jk) = rau0 * un(:,:,jk) * e2u(:,:) * e3u_n(:,:,jk) * umask(:,:,jk)
294            z2d(:,:) = z2d(:,:) + z3d(:,:,jk)
295         END DO
296         CALL iom_put( "u_masstr"     , z3d )         ! mass transport in i-direction
297         CALL iom_put( "u_masstr_vint", z2d )         ! mass transport in i-direction vertical sum
298      ENDIF
299     
300      IF( iom_use("u_heattr") ) THEN
301         z2d(:,:) = 0._wp 
302         DO jk = 1, jpkm1
303            DO jj = 2, jpjm1
304               DO ji = fs_2, fs_jpim1   ! vector opt.
305                  z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_tem) + tsn(ji+1,jj,jk,jp_tem) )
306               END DO
307            END DO
308         END DO
[10170]309         CALL lbc_lnk( 'diawri', z2d, 'U', -1. )
[9762]310         CALL iom_put( "u_heattr", 0.5*rcp * z2d )    ! heat transport in i-direction
311      ENDIF
312
313      IF( iom_use("u_salttr") ) THEN
314         z2d(:,:) = 0.e0 
315         DO jk = 1, jpkm1
316            DO jj = 2, jpjm1
317               DO ji = fs_2, fs_jpim1   ! vector opt.
318                  z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_sal) + tsn(ji+1,jj,jk,jp_sal) )
319               END DO
320            END DO
321         END DO
[10170]322         CALL lbc_lnk( 'diawri', z2d, 'U', -1. )
[9762]323         CALL iom_put( "u_salttr", 0.5 * z2d )        ! heat transport in i-direction
324      ENDIF
325
326     
327      IF( iom_use("v_masstr") .OR. iom_use("v_heattr") .OR. iom_use("v_salttr") ) THEN
328         z3d(:,:,jpk) = 0.e0
329         DO jk = 1, jpkm1
330            z3d(:,:,jk) = rau0 * vn(:,:,jk) * e1v(:,:) * e3v_n(:,:,jk) * vmask(:,:,jk)
331         END DO
332         CALL iom_put( "v_masstr", z3d )              ! mass transport in j-direction
333      ENDIF
334     
335      IF( iom_use("v_heattr") ) THEN
336         z2d(:,:) = 0.e0 
337         DO jk = 1, jpkm1
338            DO jj = 2, jpjm1
339               DO ji = fs_2, fs_jpim1   ! vector opt.
340                  z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_tem) + tsn(ji,jj+1,jk,jp_tem) )
341               END DO
342            END DO
343         END DO
[10170]344         CALL lbc_lnk( 'diawri', z2d, 'V', -1. )
[9762]345         CALL iom_put( "v_heattr", 0.5*rcp * z2d )    !  heat transport in j-direction
346      ENDIF
347
348      IF( iom_use("v_salttr") ) THEN
349         z2d(:,:) = 0._wp 
350         DO jk = 1, jpkm1
351            DO jj = 2, jpjm1
352               DO ji = fs_2, fs_jpim1   ! vector opt.
353                  z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_sal) + tsn(ji,jj+1,jk,jp_sal) )
354               END DO
355            END DO
356         END DO
[10170]357         CALL lbc_lnk( 'diawri', z2d, 'V', -1. )
[9762]358         CALL iom_put( "v_salttr", 0.5 * z2d )        !  heat transport in j-direction
359      ENDIF
360
361      IF( iom_use("tosmint") ) THEN
362         z2d(:,:) = 0._wp
363         DO jk = 1, jpkm1
364            DO jj = 2, jpjm1
365               DO ji = fs_2, fs_jpim1   ! vector opt.
366                  z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) *  tsn(ji,jj,jk,jp_tem)
367               END DO
368            END DO
369         END DO
[10170]370         CALL lbc_lnk( 'diawri', z2d, 'T', -1. )
[9762]371         CALL iom_put( "tosmint", rau0 * z2d )        ! Vertical integral of temperature
372      ENDIF
373      IF( iom_use("somint") ) THEN
374         z2d(:,:)=0._wp
375         DO jk = 1, jpkm1
376            DO jj = 2, jpjm1
377               DO ji = fs_2, fs_jpim1   ! vector opt.
378                  z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_sal)
379               END DO
380            END DO
381         END DO
[10170]382         CALL lbc_lnk( 'diawri', z2d, 'T', -1. )
[9762]383         CALL iom_put( "somint", rau0 * z2d )         ! Vertical integral of salinity
384      ENDIF
385
386      CALL iom_put( "bn2", rn2 )                      ! Brunt-Vaisala buoyancy frequency (N^2)
387      !
388
389      IF (ln_diatmb)   CALL dia_tmb                   ! tmb values
390         
391      IF (ln_dia25h)   CALL dia_25h( kt )             ! 25h averaging
392
393      IF( ln_timing )   CALL timing_stop('dia_wri')
394      !
395   END SUBROUTINE dia_wri
396
397#else
398   !!----------------------------------------------------------------------
399   !!   Default option                                  use IOIPSL  library
400   !!----------------------------------------------------------------------
401
402   INTEGER FUNCTION dia_wri_alloc()
403      !
404      dia_wri_alloc = 0
405      !
406   END FUNCTION dia_wri_alloc
407
408   SUBROUTINE dia_wri( kt )
409      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
[10179]410
411      IF( ninist == 1 ) THEN     !==  Output the initial state and forcings  ==!
412         CALL dia_wri_state( 'output.init', kt )
413         ninist = 0
414      ENDIF
415
[9762]416   END SUBROUTINE dia_wri
417
418#endif
419
420   SUBROUTINE dia_wri_state( cdfile_name, kt )
421      !!---------------------------------------------------------------------
422      !!                 ***  ROUTINE dia_wri_state  ***
423      !!       
424      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
425      !!      the instantaneous ocean state and forcing fields.
426      !!        Used to find errors in the initial state or save the last
427      !!      ocean state in case of abnormal end of a simulation
428      !!
429      !! ** Method  :   NetCDF files using ioipsl
430      !!      File 'output.init.nc'  is created if ninist = 1 (namelist)
431      !!      File 'output.abort.nc' is created in case of abnormal job end
432      !!----------------------------------------------------------------------
433      CHARACTER (len=* ), INTENT( in ) ::   cdfile_name      ! name of the file created
434      INTEGER           , INTENT( in ) ::   kt               ! ocean time-step index
435      !!
436      CHARACTER (len=32) :: clname
437      CHARACTER (len=40) :: clop
438      INTEGER  ::   id_i , nz_i, nh_i       
439      INTEGER, DIMENSION(1) ::   idex             ! local workspace
440      REAL(wp) ::   zsto, zout, zmax, zjulian
441      !!----------------------------------------------------------------------
442      !
443      ! 0. Initialisation
444      ! -----------------
445
446      ! Define name, frequency of output and means
447      clname = cdfile_name
448      IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
449      zsto = rdt
450      clop = "inst(x)"           ! no use of the mask value (require less cpu time)
451      zout = rdt
452      zmax = ( nitend - nit000 + 1 ) * rdt
453
454      IF(lwp) WRITE(numout,*)
455      IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
456      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created '
457      IF(lwp) WRITE(numout,*) '                and named :', clname, '.nc'
458
459
460      ! 1. Define NETCDF files and fields at beginning of first time step
461      ! -----------------------------------------------------------------
462
463      ! Compute julian date from starting date of the run
464      CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )         ! time axis
465      zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
466      CALL histbeg( clname, jpi, glamt, jpj, gphit,   &
467          1, jpi, 1, jpj, nit000-1, zjulian, rdt, nh_i, id_i, domain_id=nidom, snc4chunks=snc4set ) ! Horizontal grid : glamt and gphit
468      CALL histvert( id_i, "deptht", "Vertical T levels",   &    ! Vertical grid : gdept
469          "m", jpk, gdept_1d, nz_i, "down")
470
471      ! Declare all the output fields as NetCDF variables
472
473      CALL histdef( id_i, "vosaline", "Salinity"              , "PSU"    ,   &   ! salinity
474         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
475      CALL histdef( id_i, "votemper", "Temperature"           , "C"      ,   &   ! temperature
476         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
477      CALL histdef( id_i, "sossheig", "Sea Surface Height"    , "m"      ,   &  ! ssh
478         &          jpi, jpj, nh_i, 1  , 1, 1  , nz_i, 32, clop, zsto, zout )
479      CALL histdef( id_i, "vozocrtx", "Zonal Current"         , "m/s"    ,   &   ! zonal current
480         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
481      CALL histdef( id_i, "vomecrty", "Meridional Current"    , "m/s"    ,   &   ! meridonal current
482         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
483      CALL histdef( id_i, "vovecrtz", "Vertical Velocity"     , "m/s"    ,   &   ! vertical current
484         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
485         !
486      IF( ALLOCATED(ahtu) ) THEN
487         CALL histdef( id_i, "ahtu"    , "u-eddy diffusivity"    , "m2/s"    ,   &   ! zonal current
488            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
489         CALL histdef( id_i, "ahtv"    , "v-eddy diffusivity"    , "m2/s"    ,   &   ! meridonal current
490            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
491      ENDIF
492      IF( ALLOCATED(ahmt) ) THEN
493         CALL histdef( id_i, "ahmt"    , "t-eddy viscosity"      , "m2/s"    ,   &   ! zonal current
494            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
495         CALL histdef( id_i, "ahmf"    , "f-eddy viscosity"      , "m2/s"    ,   &   ! meridonal current
496            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
497      ENDIF
498         !
499      CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S",   &   ! net freshwater
500         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
501      CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2"   ,   &   ! net heat flux
502         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
503      CALL histdef( id_i, "soshfldo", "Shortwave Radiation"   , "W/m2"   ,   &   ! solar flux
504         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
505      CALL histdef( id_i, "soicecov", "Ice fraction"          , "[0,1]"  ,   &   ! fr_i
506         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
507      CALL histdef( id_i, "sozotaux", "Zonal Wind Stress"     , "N/m2"   ,   &   ! i-wind stress
508         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
509      CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2"   ,   &   ! j-wind stress
510         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
511      IF( .NOT.ln_linssh ) THEN
512         CALL histdef( id_i, "vovvldep", "T point depth"         , "m"      , &   ! t-point depth
513            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
514         CALL histdef( id_i, "vovvle3t", "T point thickness"     , "m"      , &   ! t-point depth
515            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
516      ENDIF
517      !
518      IF( ln_wave .AND. ln_sdw ) THEN
519         CALL histdef( id_i, "sdzocrtx", "Stokes Drift Zonal"    , "m/s"    , &   ! StokesDrift zonal current
520            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
521         CALL histdef( id_i, "sdmecrty", "Stokes Drift Merid"    , "m/s"    , &   ! StokesDrift meridonal current
522            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
523         CALL histdef( id_i, "sdvecrtz", "Stokes Drift Vert"     , "m/s"    , &   ! StokesDrift vertical current
524            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
525      ENDIF
526
527#if defined key_si3
528      IF( nn_ice == 2 ) THEN   ! condition needed in case agrif + ice-model but no-ice in child grid
529         CALL ice_wri_state( kt, id_i, nh_i )
530      ENDIF
531#else
532      CALL histend( id_i, snc4chunks=snc4set )
533#endif
534
535      ! 2. Start writing data
536      ! ---------------------
537      ! idex(1) est utilise ssi l'avant dernier argument est diffferent de
538      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
539      ! donne le nombre d'elements, et idex la liste des indices a sortir
540      idex(1) = 1   ! init to avoid compil warning
541
542      ! Write all fields on T grid
543      CALL histwrite( id_i, "votemper", kt, tsn(:,:,:,jp_tem), jpi*jpj*jpk, idex )    ! now temperature
544      CALL histwrite( id_i, "vosaline", kt, tsn(:,:,:,jp_sal), jpi*jpj*jpk, idex )    ! now salinity
545      CALL histwrite( id_i, "sossheig", kt, sshn             , jpi*jpj    , idex )    ! sea surface height
546      CALL histwrite( id_i, "vozocrtx", kt, un               , jpi*jpj*jpk, idex )    ! now i-velocity
547      CALL histwrite( id_i, "vomecrty", kt, vn               , jpi*jpj*jpk, idex )    ! now j-velocity
548      CALL histwrite( id_i, "vovecrtz", kt, wn               , jpi*jpj*jpk, idex )    ! now k-velocity
549      !
550      IF( ALLOCATED(ahtu) ) THEN
551         CALL histwrite( id_i, "ahtu"    , kt, ahtu             , jpi*jpj*jpk, idex )    ! aht at u-point
552         CALL histwrite( id_i, "ahtv"    , kt, ahtv             , jpi*jpj*jpk, idex )    !  -  at v-point
553      ENDIF
554      IF( ALLOCATED(ahmt) ) THEN
555         CALL histwrite( id_i, "ahmt"    , kt, ahmt             , jpi*jpj*jpk, idex )    ! ahm at t-point
556         CALL histwrite( id_i, "ahmf"    , kt, ahmf             , jpi*jpj*jpk, idex )    !  -  at f-point
557      ENDIF
558      !
559      CALL histwrite( id_i, "sowaflup", kt, emp - rnf        , jpi*jpj    , idex )    ! freshwater budget
560      CALL histwrite( id_i, "sohefldo", kt, qsr + qns        , jpi*jpj    , idex )    ! total heat flux
561      CALL histwrite( id_i, "soshfldo", kt, qsr              , jpi*jpj    , idex )    ! solar heat flux
562      CALL histwrite( id_i, "soicecov", kt, fr_i             , jpi*jpj    , idex )    ! ice fraction
563      CALL histwrite( id_i, "sozotaux", kt, utau             , jpi*jpj    , idex )    ! i-wind stress
564      CALL histwrite( id_i, "sometauy", kt, vtau             , jpi*jpj    , idex )    ! j-wind stress
565
566      IF(  .NOT.ln_linssh  ) THEN             
567         CALL histwrite( id_i, "vovvldep", kt, gdept_n(:,:,:), jpi*jpj*jpk, idex )!  T-cell depth
568         CALL histwrite( id_i, "vovvle3t", kt, e3t_n (:,:,:) , jpi*jpj*jpk, idex )!  T-cell thickness 
569      END IF
570 
571      IF( ln_wave .AND. ln_sdw ) THEN
572         CALL histwrite( id_i, "sdzocrtx", kt, usd           , jpi*jpj*jpk, idex)     ! now StokesDrift i-velocity
573         CALL histwrite( id_i, "sdmecrty", kt, vsd           , jpi*jpj*jpk, idex)     ! now StokesDrift j-velocity
574         CALL histwrite( id_i, "sdvecrtz", kt, wsd           , jpi*jpj*jpk, idex)     ! now StokesDrift k-velocity
575      ENDIF
576
577      ! 3. Close the file
578      ! -----------------
579      CALL histclo( id_i )
580#if ! defined key_iomput
581      IF( ninist /= 1  ) THEN
582         CALL histclo( nid_T )
583         CALL histclo( nid_U )
584         CALL histclo( nid_V )
585         CALL histclo( nid_W )
586      ENDIF
587#endif
588      !
589   END SUBROUTINE dia_wri_state
590
591   !!======================================================================
592END MODULE diawri
Note: See TracBrowser for help on using the repository browser.