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 branches/UKMO/r5518_rm_um_cpl/NEMOGCM/NEMO/SAS_SRC – NEMO

Context Navigation

source: branches/UKMO/r5518_rm_um_cpl/NEMOGCM/NEMO/SAS_SRC/diawri.F90 @ 7141

Last change on this file since 7141 was 7141, checked in by jcastill, 7 years ago
Remove svn keywords
File size: 21.9 KB

Line
1	MODULE 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 dom_oce ! ocean space and time domain
27	USE zdf_oce ! ocean vertical physics
28	USE sbc_oce ! Surface boundary condition: ocean fields
29	USE sbc_ice ! Surface boundary condition: ice fields
30	USE sbcssr ! restoring term toward SST/SSS climatology
31	USE phycst ! physical constants
32	USE zdfmxl ! mixed layer
33	USE dianam ! build name of file (routine)
34	USE zdfddm ! vertical physics: double diffusion
35	USE diahth ! thermocline diagnostics
36	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
37	USE in_out_manager ! I/O manager
38	USE diadimg ! dimg direct access file format output
39	USE diaar5, ONLY : lk_diaar5
40	USE iom
41	USE ioipsl
42	#if defined key_lim2
43	USE limwri_2
44	#endif
45	USE lib_mpp ! MPP library
46	USE timing ! preformance summary
47	USE wrk_nemo ! working array
48
49	IMPLICIT NONE
50	PRIVATE
51
52	PUBLIC dia_wri ! routines called by step.F90
53	PUBLIC dia_wri_state
54	PUBLIC dia_wri_alloc ! Called by nemogcm module
55
56	INTEGER :: nid_T, nz_T, nh_T, ndim_T, ndim_hT ! grid_T file
57	INTEGER :: nid_U, nz_U, nh_U, ndim_U, ndim_hU ! grid_U file
58	INTEGER :: nid_V, nz_V, nh_V, ndim_V, ndim_hV ! grid_V file
59	INTEGER :: ndex(1) ! ???
60	INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV
61
62	!! * Substitutions
63	# include "zdfddm_substitute.h90"
64	# include "domzgr_substitute.h90"
65	# include "vectopt_loop_substitute.h90"
66	!!----------------------------------------------------------------------
67	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
68	!! $Id$
69	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
70	!!----------------------------------------------------------------------
71	CONTAINS
72
73	INTEGER FUNCTION dia_wri_alloc()
74	!!----------------------------------------------------------------------
75	INTEGER :: ierr
76	!!----------------------------------------------------------------------
77	!
78	ALLOCATE( ndex_hT(jpijpj), ndex_hU(jpijpj), ndex_hV(jpi*jpj), STAT=dia_wri_alloc )
79	IF( lk_mpp ) CALL mpp_sum( dia_wri_alloc )
80	!
81	END FUNCTION dia_wri_alloc
82
83	#if defined key_dimgout
84	!!----------------------------------------------------------------------
85	!! 'key_dimgout' DIMG output file
86	!!----------------------------------------------------------------------
87	# include "diawri_dimg.h90"
88
89	#else
90	!!----------------------------------------------------------------------
91	!! Default option NetCDF output file
92	!!----------------------------------------------------------------------
93	# if defined key_iomput
94	!!----------------------------------------------------------------------
95	!! 'key_iomput' use IOM library
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	!! Standalone surface scheme
105	!!
106	!! ** Method : use iom_put
107	!!----------------------------------------------------------------------
108	!!
109	INTEGER, INTENT( in ) :: kt ! ocean time-step index
110	!!----------------------------------------------------------------------
111	!
112	!! no relevant 2D arrays to write in iomput case
113	!
114	END SUBROUTINE dia_wri
115
116	#else
117	!!----------------------------------------------------------------------
118	!! Default option use IOIPSL library
119	!!----------------------------------------------------------------------
120
121	SUBROUTINE dia_wri( kt )
122	!!---------------------------------------------------------------------
123	!! * ROUTINE dia_wri *
124	!!
125	!! ** Purpose : Standard output of opa: dynamics and tracer fields
126	!! NETCDF format is used by default
127	!!
128	!! ** Method : At the beginning of the first time step (nit000),
129	!! define all the NETCDF files and fields
130	!! At each time step call histdef to compute the mean if ncessary
131	!! Each nwrite time step, output the instantaneous or mean fields
132	!!----------------------------------------------------------------------
133	!!
134	INTEGER, INTENT( in ) :: kt ! ocean time-step index
135	!!
136	LOGICAL :: ll_print = .FALSE. ! =T print and flush numout
137	CHARACTER (len=40) :: clhstnam, clop, clmx ! local names
138	INTEGER :: inum = 11 ! temporary logical unit
139	INTEGER :: ji, jj, jk ! dummy loop indices
140	INTEGER :: ierr ! error code return from allocation
141	INTEGER :: iimi, iima, ipk, it, itmod, ijmi, ijma ! local integers
142	REAL(wp) :: zsto, zout, zmax, zjulian, zdt ! local scalars
143	!!----------------------------------------------------------------------
144	!
145	IF( nn_timing == 1 ) CALL timing_start('dia_wri')
146	!
147	! Output the initial state and forcings
148	IF( ninist == 1 ) THEN
149	CALL dia_wri_state( 'output.init', kt )
150	ninist = 0
151	ENDIF
152	!
153	! 0. Initialisation
154	! -----------------
155
156	! local variable for debugging
157	ll_print = .FALSE.
158	ll_print = ll_print .AND. lwp
159
160	! Define frequency of output and means
161	zdt = rdt
162	IF( nacc == 1 ) zdt = rdtmin
163	IF( ln_mskland ) THEN ; clop = "only(x)" ! put 1.e+20 on land (very expensive!!)
164	ELSE ; clop = "x" ! no use of the mask value (require less cpu time)
165	ENDIF
166	#if defined key_diainstant
167	zsto = nwrite * zdt
168	clop = "inst("//TRIM(clop)//")"
169	#else
170	zsto=zdt
171	clop = "ave("//TRIM(clop)//")"
172	#endif
173	zout = nwrite * zdt
174	zmax = ( nitend - nit000 + 1 ) * zdt
175
176	! Define indices of the horizontal output zoom and vertical limit storage
177	iimi = 1 ; iima = jpi
178	ijmi = 1 ; ijma = jpj
179	ipk = jpk
180
181	! define time axis
182	it = kt
183	itmod = kt - nit000 + 1
184
185
186	! 1. Define NETCDF files and fields at beginning of first time step
187	! -----------------------------------------------------------------
188
189	IF( kt == nit000 ) THEN
190
191	! Define the NETCDF files (one per grid)
192
193	! Compute julian date from starting date of the run
194	CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
195	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
196	IF(lwp)WRITE(numout,*)
197	IF(lwp)WRITE(numout,*) 'Date 0 used :', nit000, ' YEAR ', nyear, &
198	& ' MONTH ', nmonth, ' DAY ', nday, 'Julian day : ', zjulian
199	IF(lwp)WRITE(numout,*) ' indexes of zoom = ', iimi, iima, ijmi, ijma, &
200	' limit storage in depth = ', ipk
201
202	! WRITE root name in date.file for use by postpro
203	IF(lwp) THEN
204	CALL dia_nam( clhstnam, nwrite,' ' )
205	CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
206	WRITE(inum,*) clhstnam
207	CLOSE(inum)
208	ENDIF
209
210	! Define the T grid FILE ( nid_T )
211
212	CALL dia_nam( clhstnam, nwrite, 'grid_T' )
213	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
214	CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit
215	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
216	& nit000-1, zjulian, zdt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set )
217	CALL histvert( nid_T, "deptht", "Vertical T levels", & ! Vertical grid: gdept
218	& "m", ipk, gdept_1d, nz_T, "down" )
219	! ! Index of ocean points
220	CALL wheneq( jpi*jpj , tmask, 1, 1., ndex_hT, ndim_hT ) ! surface
221
222	! Define the U grid FILE ( nid_U )
223
224	CALL dia_nam( clhstnam, nwrite, 'grid_U' )
225	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename
226	CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu, & ! Horizontal grid: glamu and gphiu
227	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
228	& nit000-1, zjulian, zdt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set )
229	CALL histvert( nid_U, "depthu", "Vertical U levels", & ! Vertical grid: gdept
230	& "m", ipk, gdept_1d, nz_U, "down" )
231	! ! Index of ocean points
232	CALL wheneq( jpi*jpj , umask, 1, 1., ndex_hU, ndim_hU ) ! surface
233
234	! Define the V grid FILE ( nid_V )
235
236	CALL dia_nam( clhstnam, nwrite, 'grid_V' ) ! filename
237	IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam
238	CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv, & ! Horizontal grid: glamv and gphiv
239	& iimi, iima-iimi+1, ijmi, ijma-ijmi+1, &
240	& nit000-1, zjulian, zdt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set )
241	CALL histvert( nid_V, "depthv", "Vertical V levels", & ! Vertical grid : gdept
242	& "m", ipk, gdept_1d, nz_V, "down" )
243	! ! Index of ocean points
244	CALL wheneq( jpi*jpj , vmask, 1, 1., ndex_hV, ndim_hV ) ! surface
245
246	! No W grid FILE
247
248	! Declare all the output fields as NETCDF variables
249
250	! !!! nid_T : 3D
251	CALL histdef( nid_T, "sst_m", "Sea Surface temperature" , "C" , & ! sst
252	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
253	CALL histdef( nid_T, "sss_m", "Sea Surface Salinity" , "PSU" , & ! sss
254	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
255	CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux" , "Kg/m2/s", & ! (emp-rnf)
256	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
257	CALL histdef( nid_T, "sosfldow", "downward salt flux" , "PSU/m2/s", & ! (sfx)
258	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
259	CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! qns + qsr
260	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
261	CALL histdef( nid_T, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! qsr
262	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
263	CALL histdef( nid_T, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
264	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
265	CALL histdef( nid_T, "sowindsp", "wind speed at 10m" , "m/s" , & ! wndm
266	& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
267
268	CALL histend( nid_T, snc4chunks=snc4set )
269
270	! !!! nid_U : 3D
271	CALL histdef( nid_U, "ssu_m", "Velocity component in x-direction", "m/s" , & ! ssu
272	& jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
273	CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau
274	& jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
275
276	CALL histend( nid_U, snc4chunks=snc4set )
277
278	! !!! nid_V : 3D
279	CALL histdef( nid_V, "ssv_m", "Velocity component in y-direction", "m/s", & ! ssv_m
280	& jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
281	CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau
282	& jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
283
284	CALL histend( nid_V, snc4chunks=snc4set )
285
286	IF(lwp) WRITE(numout,*)
287	IF(lwp) WRITE(numout,*) 'End of NetCDF Initialization'
288	IF(ll_print) CALL FLUSH(numout )
289
290	ENDIF
291
292	! 2. Start writing data
293	! ---------------------
294
295	! ndex(1) est utilise ssi l'avant dernier argument est diffferent de
296	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
297	! donne le nombre d'elements, et ndex la liste des indices a sortir
298
299	IF( lwp .AND. MOD( itmod, nwrite ) == 0 ) THEN
300	WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step'
301	WRITE(numout,*) '~~~~~~ '
302	ENDIF
303
304	! Write fields on T grid
305	CALL histwrite( nid_T, "sst_m", it, sst_m, ndim_hT, ndex_hT ) ! sea surface temperature
306	CALL histwrite( nid_T, "sss_m", it, sss_m, ndim_hT, ndex_hT ) ! sea surface salinity
307	CALL histwrite( nid_T, "sowaflup", it, (emp - rnf ) , ndim_hT, ndex_hT ) ! upward water flux
308	CALL histwrite( nid_T, "sosfldow", it, sfx , ndim_hT, ndex_hT ) ! downward salt flux
309	! (includes virtual salt flux beneath ice
310	! in linear free surface case)
311
312	CALL histwrite( nid_T, "sohefldo", it, qns + qsr , ndim_hT, ndex_hT ) ! total heat flux
313	CALL histwrite( nid_T, "soshfldo", it, qsr , ndim_hT, ndex_hT ) ! solar heat flux
314	CALL histwrite( nid_T, "soicecov", it, fr_i , ndim_hT, ndex_hT ) ! ice fraction
315	CALL histwrite( nid_T, "sowindsp", it, wndm , ndim_hT, ndex_hT ) ! wind speed
316
317	! Write fields on U grid
318	CALL histwrite( nid_U, "ssu_m" , it, ssu_m , ndim_hU, ndex_hU ) ! i-current speed
319	CALL histwrite( nid_U, "sozotaux", it, utau , ndim_hU, ndex_hU ) ! i-wind stress
320
321	! Write fields on V grid
322	CALL histwrite( nid_V, "ssv_m" , it, ssv_m , ndim_hV, ndex_hV ) ! j-current speed
323	CALL histwrite( nid_V, "sometauy", it, vtau , ndim_hV, ndex_hV ) ! j-wind stress
324
325	! 3. Close all files
326	! ---------------------------------------
327	IF( kt == nitend ) THEN
328	CALL histclo( nid_T )
329	CALL histclo( nid_U )
330	CALL histclo( nid_V )
331	ENDIF
332	!
333	IF( nn_timing == 1 ) CALL timing_stop('dia_wri')
334	!
335	END SUBROUTINE dia_wri
336	# endif
337
338	#endif
339
340	SUBROUTINE dia_wri_state( cdfile_name, kt )
341	!!---------------------------------------------------------------------
342	!! * ROUTINE dia_wri_state *
343	!!
344	!! ** Purpose : create a NetCDF file named cdfile_name which contains
345	!! the instantaneous ocean state and forcing fields.
346	!! Used to find errors in the initial state or save the last
347	!! ocean state in case of abnormal end of a simulation
348	!!
349	!! ** Method : NetCDF files using ioipsl
350	!! File 'output.init.nc' is created if ninist = 1 (namelist)
351	!! File 'output.abort.nc' is created in case of abnormal job end
352	!!----------------------------------------------------------------------
353	CHARACTER (len=* ), INTENT( in ) :: cdfile_name ! name of the file created
354	INTEGER , INTENT( in ) :: kt ! ocean time-step index
355	!!
356	CHARACTER (len=32) :: clname
357	CHARACTER (len=40) :: clop
358	INTEGER :: id_i , nz_i, nh_i
359	INTEGER, DIMENSION(1) :: idex ! local workspace
360	REAL(wp) :: zsto, zout, zmax, zjulian, zdt
361	!!----------------------------------------------------------------------
362	!
363	IF( nn_timing == 1 ) CALL timing_start('dia_wri_state')
364
365	! 0. Initialisation
366	! -----------------
367
368	! Define name, frequency of output and means
369	clname = cdfile_name
370	IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
371	zdt = rdt
372	zsto = rdt
373	clop = "inst(x)" ! no use of the mask value (require less cpu time)
374	zout = rdt
375	zmax = ( nitend - nit000 + 1 ) * zdt
376
377	IF(lwp) WRITE(numout,*)
378	IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
379	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ and forcing fields file created '
380	IF(lwp) WRITE(numout,*) ' and named :', clname, '.nc'
381
382
383	! 1. Define NETCDF files and fields at beginning of first time step
384	! -----------------------------------------------------------------
385
386	! Compute julian date from starting date of the run
387	CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian ) ! time axis
388	zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment
389	CALL histbeg( clname, jpi, glamt, jpj, gphit, &
390	1, jpi, 1, jpj, nit000-1, zjulian, zdt, nh_i, id_i, domain_id=nidom, snc4chunks=snc4set ) ! Horizontal grid : glamt and gphit
391	CALL histvert( id_i, "deptht", "Vertical T levels", & ! Vertical grid : gdept
392	"m", jpk, gdept_1d, nz_i, "down")
393
394	! Declare all the output fields as NetCDF variables
395
396	CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S", & ! net freshwater
397	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
398	CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2" , & ! net heat flux
399	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
400	CALL histdef( id_i, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! solar flux
401	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
402	CALL histdef( id_i, "soicecov", "Ice fraction" , "[0,1]" , & ! fr_i
403	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
404	CALL histdef( id_i, "sozotaux", "Zonal Wind Stress" , "N/m2" , & ! i-wind stress
405	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
406	CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2" , & ! j-wind stress
407	& jpi, jpj, nh_i, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
408
409	#if defined key_lim2
410	CALL lim_wri_state_2( kt, id_i, nh_i )
411	#else
412	CALL histend( id_i, snc4chunks=snc4set )
413	#endif
414
415	! 2. Start writing data
416	! ---------------------
417	! idex(1) est utilise ssi l'avant dernier argument est diffferent de
418	! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
419	! donne le nombre d'elements, et idex la liste des indices a sortir
420	idex(1) = 1 ! init to avoid compil warning
421
422	! Write all fields on T grid
423	CALL histwrite( id_i, "sowaflup", kt, emp , jpi*jpj , idex ) ! freshwater budget
424	CALL histwrite( id_i, "sohefldo", kt, qsr + qns , jpi*jpj , idex ) ! total heat flux
425	CALL histwrite( id_i, "soshfldo", kt, qsr , jpi*jpj , idex ) ! solar heat flux
426	CALL histwrite( id_i, "soicecov", kt, fr_i , jpi*jpj , idex ) ! ice fraction
427	CALL histwrite( id_i, "sozotaux", kt, utau , jpi*jpj , idex ) ! i-wind stress
428	CALL histwrite( id_i, "sometauy", kt, vtau , jpi*jpj , idex ) ! j-wind stress
429
430	! 3. Close the file
431	! -----------------
432	CALL histclo( id_i )
433	#if ! defined key_iomput && ! defined key_dimgout
434	IF( ninist /= 1 ) THEN
435	CALL histclo( nid_T )
436	CALL histclo( nid_U )
437	CALL histclo( nid_V )
438	ENDIF
439	#endif
440
441	IF( nn_timing == 1 ) CALL timing_stop('dia_wri_state')
442	!
443
444	END SUBROUTINE dia_wri_state
445	!!======================================================================
446	END MODULE diawri

Note: See TracBrowser for help on using the repository browser.

Download in other formats: