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.

trcrst.F90 in branches/UKMO/dev_r5518_GO6_package_FOAMv14_phytobal3d/NEMOGCM/NEMO/TOP_SRC – NEMO

Context Navigation

source: branches/UKMO/dev_r5518_GO6_package_FOAMv14_phytobal3d/NEMOGCM/NEMO/TOP_SRC/trcrst.F90 @ 13097

Last change on this file since 13097 was 13097, checked in by dford, 4 years ago
Allow nitrogen balancing for both 2D and 3D chlorophyll increments.
File size: 44.1 KB

Line
1	MODULE trcrst
2	!!======================================================================
3	!! * MODULE trcrst *
4	!! TOP : Manage the passive tracer restart
5	!!======================================================================
6	!! History : - ! 1991-03 () original code
7	!! 1.0 ! 2005-03 (O. Aumont, A. El Moussaoui) F90
8	!! - ! 2005-10 (C. Ethe) print control
9	!! 2.0 ! 2005-10 (C. Ethe, G. Madec) revised architecture
10	!!----------------------------------------------------------------------
11	#if defined key_top
12	!!----------------------------------------------------------------------
13	!! 'key_top' TOP models
14	!!----------------------------------------------------------------------
15	!!----------------------------------------------------------------------
16	!! trc_rst : Restart for passive tracer
17	!!----------------------------------------------------------------------
18	!!----------------------------------------------------------------------
19	!! 'key_top' TOP models
20	!!----------------------------------------------------------------------
21	!! trc_rst_opn : open restart file
22	!! trc_rst_read : read restart file
23	!! trc_rst_wri : write restart file
24	!!----------------------------------------------------------------------
25	USE oce_trc
26	USE trc
27	USE trcnam_trp
28	USE iom
29	USE ioipsl, ONLY : ju2ymds ! for calendar
30	USE daymod
31	!! AXY (05/11/13): need these for MEDUSA to input/output benthic reservoirs
32	USE par_medusa
33	USE sms_medusa
34	USE trcsms_medusa
35	!!
36	#if defined key_idtra
37	USE trcsms_idtra
38	#endif
39	!!
40	#if defined key_cfc
41	USE trcsms_cfc
42	#endif
43	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
44	USE sbc_oce, ONLY: lk_oasis
45	USE oce, ONLY: CO2Flux_out_cpl, DMS_out_cpl, chloro_out_cpl !! Coupling variable
46	USE trcstat
47	USE obs_const, ONLY: obfillflt ! Observation operator fill value
48
49	IMPLICIT NONE
50	PRIVATE
51
52	PUBLIC trc_rst_opn ! called by ???
53	PUBLIC trc_rst_read ! called by ???
54	PUBLIC trc_rst_wri ! called by ???
55	PUBLIC trc_rst_cal
56	PUBLIC trc_rst_stat
57	#if defined key_medusa && defined key_roam
58	PUBLIC trc_rst_conserve
59	#endif
60
61	!! * Substitutions
62	# include "top_substitute.h90"
63
64	CONTAINS
65
66	SUBROUTINE trc_rst_opn( kt )
67	!!----------------------------------------------------------------------
68	!! * trc_rst_opn *
69	!!
70	!! ** purpose : output of sea-trc variable in a netcdf file
71	!!----------------------------------------------------------------------
72	INTEGER, INTENT(in) :: kt ! number of iteration
73	INTEGER :: iyear, imonth, iday
74	REAL (wp) :: zsec
75	REAL (wp) :: zfjulday
76	!
77	CHARACTER(LEN=20) :: clkt ! ocean time-step define as a character
78	CHARACTER(LEN=50) :: clname ! trc output restart file name
79	CHARACTER(LEN=256) :: clpath ! full path to ocean output restart file
80	!!----------------------------------------------------------------------
81	!
82	IF( lk_offline ) THEN
83	IF( kt == nittrc000 ) THEN
84	lrst_trc = .FALSE.
85	IF( ln_rst_list ) THEN
86	nrst_lst = 1
87	nitrst = nstocklist( nrst_lst )
88	ELSE
89	nitrst = nitend
90	ENDIF
91	ENDIF
92
93	IF( .NOT. ln_rst_list .AND. MOD( kt - 1, nstock ) == 0 ) THEN
94	! we use kt - 1 and not kt - nittrc000 to keep the same periodicity from the beginning of the experiment
95	nitrst = kt + nstock - 1 ! define the next value of nitrst for restart writing
96	IF( nitrst > nitend ) nitrst = nitend ! make sure we write a restart at the end of the run
97	ENDIF
98	ELSE
99	IF( kt == nittrc000 ) lrst_trc = .FALSE.
100	ENDIF
101
102	! to get better performances with NetCDF format:
103	! we open and define the tracer restart file one tracer time step before writing the data (-> at nitrst - 2*nn_dttrc + 1)
104	! except if we write tracer restart files every tracer time step or if a tracer restart file was writen at nitend - 2*nn_dttrc + 1
105	IF( kt == nitrst - 2*nn_dttrc .OR. nstock == nn_dttrc .OR. ( kt == nitend - nn_dttrc .AND. .NOT. lrst_trc ) ) THEN
106	IF ( ln_rstdate ) THEN
107	!! JPALM -- 22-12-2015 -- modif to get the good date on restart trc file name
108	!! -- the condition to open the rst file is not the same than for the dynamic rst.
109	!! -- here it - for an obscure reason - is open 2 time-step before the restart writing process
110	!! instead of 1.
111	!! -- i am not sure if someone forgot +1 in the if loop condition as
112	!! it is writen in all comments nitrst -2*nn_dttrc + 1 and the condition is
113	!! nitrst - 2*nn_dttrc
114	!! -- nevertheless we didn't wanted to broke something already working
115	!! and just adapted the part we added.
116	!! -- So instead of calling ju2ymds( fjulday + (rdttra(1))
117	!! we call ju2ymds( fjulday + (2*rdttra(1))
118	!!--------------------------------------------------------------------
119	zfjulday = fjulday + (2*rdttra(1)) / rday
120	IF( ABS(zfjulday - REAL(NINT(zfjulday),wp)) < 0.1 / rday ) zfjulday = REAL(NINT(zfjulday),wp) ! avoid truncation error
121	CALL ju2ymds( zfjulday + (2*rdttra(1)) / rday, iyear, imonth, iday, zsec )
122	WRITE(clkt, '(i4.4,2i2.2)') iyear, imonth, iday
123	ELSE
124	! beware of the format used to write kt (default is i8.8, that should be large enough)
125	IF( nitrst > 1.0e9 ) THEN ; WRITE(clkt,* ) nitrst
126	ELSE ; WRITE(clkt,'(i8.8)') nitrst
127	ENDIF
128	ENDIF
129	! create the file
130	IF(lwp) WRITE(numout,*)
131	clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_"//TRIM(cn_trcrst_out)
132	clpath = TRIM(cn_trcrst_outdir)
133	IF( clpath(LEN_TRIM(clpath):) /= '/' ) clpath = TRIM(clpath) // '/'
134	IF(lwp) WRITE(numout,*) &
135	' open trc restart.output NetCDF file: ',TRIM(clpath)//clname
136	CALL iom_open( TRIM(clpath)//TRIM(clname), numrtw, ldwrt = .TRUE., kiolib = jprstlib )
137	lrst_trc = .TRUE.
138	ENDIF
139	!
140	END SUBROUTINE trc_rst_opn
141
142	SUBROUTINE trc_rst_read
143	!!----------------------------------------------------------------------
144	!! * trc_rst_opn *
145	!!
146	!! ** purpose : read passive tracer fields in restart files
147	!!----------------------------------------------------------------------
148	INTEGER :: jn, jl
149	!! AXY (05/11/13): temporary variables
150	REAL(wp) :: fq0,fq1,fq2
151
152	!!----------------------------------------------------------------------
153	!
154	IF(lwp) WRITE(numout,*)
155	IF(lwp) WRITE(numout,*) 'trc_rst_read : read data in the TOP restart file'
156	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
157
158	! READ prognostic variables and computes diagnostic variable
159	DO jn = 1, jptra
160	CALL iom_get( numrtr, jpdom_autoglo, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
161	trn(:,:,:,jn) = trn(:,:,:,jn) * tmask(:,:,:)
162	END DO
163
164	DO jn = 1, jptra
165	CALL iom_get( numrtr, jpdom_autoglo, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
166	trb(:,:,:,jn) = trb(:,:,:,jn) * tmask(:,:,:)
167	END DO
168	!
169	!! AXY (09/06/14): the ARCHER version of MEDUSA-2 does not include an equivalent
170	!! call to MEDUSA-2 at this point; this suggests that the FCM
171	!! version of NEMO date significantly earlier than the current
172	!! version
173
174	#if defined key_medusa
175	!! AXY (13/01/12): check if the restart contains sediment fields;
176	!! this is only relevant for simulations that include
177	!! biogeochemistry and are restarted from earlier runs
178	!! in which there was no sediment component
179	!!
180	IF( iom_varid( numrtr, 'B_SED_N', ldstop = .FALSE. ) > 0 ) THEN
181	!! YES; in which case read them
182	!!
183	IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields present - reading in ...'
184	CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_N', zb_sed_n(:,:) )
185	CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_N', zn_sed_n(:,:) )
186	CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_FE', zb_sed_fe(:,:) )
187	CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_FE', zn_sed_fe(:,:) )
188	CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_SI', zb_sed_si(:,:) )
189	CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_SI', zn_sed_si(:,:) )
190	CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_C', zb_sed_c(:,:) )
191	CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_C', zn_sed_c(:,:) )
192	CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_CA', zb_sed_ca(:,:) )
193	CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_CA', zn_sed_ca(:,:) )
194	ELSE
195	!! NO; in which case set them to zero
196	!!
197	IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields absent - setting to zero ...'
198	zb_sed_n(:,:) = 0.0 !! organic N
199	zn_sed_n(:,:) = 0.0
200	zb_sed_fe(:,:) = 0.0 !! organic Fe
201	zn_sed_fe(:,:) = 0.0
202	zb_sed_si(:,:) = 0.0 !! inorganic Si
203	zn_sed_si(:,:) = 0.0
204	zb_sed_c(:,:) = 0.0 !! organic C
205	zn_sed_c(:,:) = 0.0
206	zb_sed_ca(:,:) = 0.0 !! inorganic C
207	zn_sed_ca(:,:) = 0.0
208	ENDIF
209	!!
210	!! calculate stats on these fields
211	IF(lwp) WRITE(numout,*) ' MEDUSA sediment field stats (min, max, sum) ...'
212	call trc_rst_dia_stat(zn_sed_n(:,:), 'Sediment N')
213	call trc_rst_dia_stat(zn_sed_fe(:,:), 'Sediment Fe')
214	call trc_rst_dia_stat(zn_sed_si(:,:), 'Sediment Si')
215	call trc_rst_dia_stat(zn_sed_c(:,:), 'Sediment C')
216	call trc_rst_dia_stat(zn_sed_ca(:,:), 'Sediment Ca')
217	!!
218	!! AXY (07/07/15): read in temporally averaged fields for DMS
219	!! calculations
220	!!
221	IF( iom_varid( numrtr, 'B_DMS_CHN', ldstop = .FALSE. ) > 0 ) THEN
222	!! YES; in which case read them
223	!!
224	IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS present - reading in ...'
225	CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_CHN', zb_dms_chn(:,:) )
226	CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_CHN', zn_dms_chn(:,:) )
227	CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_CHD', zb_dms_chd(:,:) )
228	CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_CHD', zn_dms_chd(:,:) )
229	CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_MLD', zb_dms_mld(:,:) )
230	CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_MLD', zn_dms_mld(:,:) )
231	CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_QSR', zb_dms_qsr(:,:) )
232	CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_QSR', zn_dms_qsr(:,:) )
233	CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_DIN', zb_dms_din(:,:) )
234	CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_DIN', zn_dms_din(:,:) )
235	ELSE
236	!! NO; in which case set them to zero
237	!!
238	IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS absent - setting to zero ...'
239	zb_dms_chn(:,:) = 0.0 !! CHN
240	zn_dms_chn(:,:) = 0.0
241	zb_dms_chd(:,:) = 0.0 !! CHD
242	zn_dms_chd(:,:) = 0.0
243	zb_dms_mld(:,:) = 0.0 !! MLD
244	zn_dms_mld(:,:) = 0.0
245	zb_dms_qsr(:,:) = 0.0 !! QSR
246	zn_dms_qsr(:,:) = 0.0
247	zb_dms_din(:,:) = 0.0 !! DIN
248	zn_dms_din(:,:) = 0.0
249	ENDIF
250	!!
251	!! JPALM 14-06-2016 -- add CO2 flux and DMS surf through the restart
252	!! -- needed for the coupling with atm
253	IF( iom_varid( numrtr, 'N_DMS_srf', ldstop = .FALSE. ) > 0 ) THEN
254	IF(lwp) WRITE(numout,*) 'DMS surf concentration - reading in ...'
255	CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_srf', zb_dms_srf(:,:) )
256	CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_srf', zn_dms_srf(:,:) )
257	ELSE
258	IF(lwp) WRITE(numout,*) 'DMS surf concentration - setting to zero ...'
259	zb_dms_srf(:,:) = 0.0 !! DMS
260	zn_dms_srf(:,:) = 0.0
261	ENDIF
262	IF (lk_oasis) THEN
263	DMS_out_cpl(:,:) = zn_dms_srf(:,:) !! Coupling variable
264	END IF
265	!!
266	IF( iom_varid( numrtr, 'B_CO2_flx', ldstop = .FALSE. ) > 0 ) THEN
267	IF(lwp) WRITE(numout,*) 'CO2 air-sea flux - reading in ...'
268	CALL iom_get( numrtr, jpdom_autoglo, 'B_CO2_flx', zb_co2_flx(:,:) )
269	CALL iom_get( numrtr, jpdom_autoglo, 'N_CO2_flx', zn_co2_flx(:,:) )
270	ELSE
271	IF(lwp) WRITE(numout,*) 'CO2 air-sea flux - setting to zero ...'
272	zb_co2_flx(:,:) = 0.0 !! CO2 flx
273	zn_co2_flx(:,:) = 0.0
274	ENDIF
275	IF (lk_oasis) THEN
276	CO2Flux_out_cpl(:,:) = zn_co2_flx(:,:) !! Coupling variable
277	END IF
278	!!
279	!! JPALM 02-06-2017 -- in complement to DMS surf
280	!! -- the atm model needs surf Chl
281	!! as proxy of org matter from the ocean
282	!! -- needed for the coupling with atm
283	!! 07-12-2017 -- To make things cleaner, we want to store an
284	!! unscaled Chl field in the restart and only
285	!! scale it when reading it in.
286
287	IF( iom_varid( numrtr, 'N_CHL_srf', ldstop = .FALSE. ) > 0 ) THEN
288	IF(lwp) WRITE(numout,*) 'Chl cpl concentration - reading in ... - scale by ', scl_chl
289	CALL iom_get( numrtr, jpdom_autoglo, 'N_CHL_srf', zn_chl_srf(:,:) )
290	ELSE
291	IF(lwp) WRITE(numout,*) 'set Chl coupled concentration - scaled by ', scl_chl
292	zn_chl_srf(:,:) = MAX( 0.0, (trn(:,:,1,jpchn) + trn(:,:,1,jpchd)) * 1.E-6 )
293	ENDIF
294	IF (lk_oasis) THEN
295	chloro_out_cpl(:,:) = zn_chl_srf(:,:) * scl_chl !! Coupling variable
296	END IF
297	!!
298	!! calculate stats on these fields
299	IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS stats (min, max, sum) ...'
300	call trc_rst_dia_stat(zn_dms_chn(:,:), 'DMS, CHN')
301	call trc_rst_dia_stat(zn_dms_chd(:,:), 'DMS, CHD')
302	call trc_rst_dia_stat(zn_dms_mld(:,:), 'DMS, MLD')
303	call trc_rst_dia_stat(zn_dms_qsr(:,:), 'DMS, QSR')
304	call trc_rst_dia_stat(zn_dms_din(:,:), 'DMS, DIN')
305	call trc_rst_dia_stat(zn_dms_srf(:,:), 'DMS surf')
306	call trc_rst_dia_stat(zn_co2_flx(:,:), 'CO2 flux')
307	IF (lk_oasis) THEN
308	call trc_rst_dia_stat(chloro_out_cpl(:,:), 'CHL cpl')
309	END IF
310	!!
311	!! JPALM 14-06-2016 -- add Carbonate chenistry variables through the restart
312	!! -- needed for monthly call of carb-chem routine and better reproducibility
313	# if defined key_roam
314	IF( iom_varid( numrtr, 'pH_3D', ldstop = .FALSE. ) > 0 ) THEN
315	IF(lwp) WRITE(numout,*) 'Carbonate chem variable - reading in ...'
316	CALL iom_get( numrtr, jpdom_autoglo, 'pH_3D' , f3_pH(:,:,:) )
317	CALL iom_get( numrtr, jpdom_autoglo, 'h2CO3_3D', f3_h2co3(:,:,:) )
318	CALL iom_get( numrtr, jpdom_autoglo, 'hCO3_3D' , f3_hco3(:,:,:) )
319	CALL iom_get( numrtr, jpdom_autoglo, 'CO3_3D' , f3_co3(:,:,:) )
320	CALL iom_get( numrtr, jpdom_autoglo, 'omcal_3D', f3_omcal(:,:,:) )
321	CALL iom_get( numrtr, jpdom_autoglo, 'omarg_3D', f3_omarg(:,:,:) )
322	CALL iom_get( numrtr, jpdom_autoglo, 'CCD_CAL' , f2_ccd_cal(:,:) )
323	CALL iom_get( numrtr, jpdom_autoglo, 'CCD_ARG' , f2_ccd_arg(:,:) )
324	!!
325	IF(lwp) WRITE(numout,*) ' MEDUSA averaged Carb-chem stats (min, max, sum) ...'
326	call trc_rst_dia_stat( f3_pH(:,:,1) ,'pH 3D surf')
327	call trc_rst_dia_stat( f3_h2co3(:,:,1),'h2CO3 3D surf')
328	call trc_rst_dia_stat( f3_hco3(:,:,1) ,'hCO3 3D surf' )
329	call trc_rst_dia_stat( f3_co3(:,:,1) ,'CO3 3D surf' )
330	call trc_rst_dia_stat( f3_omcal(:,:,1),'omcal 3D surf')
331	call trc_rst_dia_stat( f3_omarg(:,:,1),'omarg 3D surf')
332	call trc_rst_dia_stat( f2_ccd_cal(:,:),'CCD_CAL')
333	call trc_rst_dia_stat( f2_ccd_arg(:,:),'CCD_ARG')
334
335	ELSE
336	IF(lwp) WRITE(numout,*) 'WARNING : No Carbonate-chem variable in the restart.... '
337	IF(lwp) WRITE(numout,*) 'Is not a problem if start a month, but may be very problematic if not '
338	IF(lwp) WRITE(numout,) 'Check if mod(ktrdt,2592000) == rdt'
339	IF(lwp) WRITE(numout,*) 'Or don t start from uncomplete restart...'
340	ENDIF
341	!
342	IF ( ln_foam_medusa ) THEN
343	!! 2D fields of pCO2 and fCO2 for observation operator on first timestep
344	IF( iom_varid( numrtr, 'PCO2W', ldstop = .FALSE. ) > 0 ) THEN
345	IF(lwp) WRITE(numout,*) ' MEDUSA pCO2 present - reading in ...'
346	CALL iom_get( numrtr, jpdom_autoglo, 'PCO2W', f2_pco2w(:,:) )
347	CALL iom_get( numrtr, jpdom_autoglo, 'FCO2W', f2_fco2w(:,:) )
348	ELSE
349	IF(lwp) WRITE(numout,*) ' MEDUSA pCO2 absent - setting to fill ...'
350	f2_pco2w(:,:) = obfillflt * tmask(:,:,1)
351	f2_fco2w(:,:) = obfillflt * tmask(:,:,1)
352	ENDIF
353	ENDIF
354	# endif
355	IF ( ln_foam_medusa ) THEN
356	!! Fields for ocean colour assimilation on first timestep
357	IF( iom_varid( numrtr, 'pgrow_avg', ldstop = .FALSE. ) > 0 ) THEN
358	IF(lwp) WRITE(numout,*) ' MEDUSA pgrow_avg present - reading in ...'
359	CALL iom_get( numrtr, jpdom_autoglo, 'pgrow_avg', pgrow_avg(:,:) )
360	CALL iom_get( numrtr, jpdom_autoglo, 'ploss_avg', ploss_avg(:,:) )
361	CALL iom_get( numrtr, jpdom_autoglo, 'phyt_avg', phyt_avg(:,:) )
362	CALL iom_get( numrtr, jpdom_autoglo, 'mld_max', mld_max(:,:) )
363	ELSE
364	IF(lwp) WRITE(numout,*) ' MEDUSA pgrow_avg absent - setting to zero ...'
365	pgrow_avg(:,:) = 0.0
366	ploss_avg(:,:) = 0.0
367	phyt_avg(:,:) = 0.0
368	mld_max(:,:) = 0.0
369	ENDIF
370	IF( iom_varid( numrtr, 'pgrow_avg_3d', ldstop = .FALSE. ) > 0 ) THEN
371	IF(lwp) WRITE(numout,*) ' MEDUSA pgrow_avg present - reading in ...'
372	CALL iom_get( numrtr, jpdom_autoglo, 'pgrow_avg_3d', pgrow_avg_3d(:,:,:) )
373	CALL iom_get( numrtr, jpdom_autoglo, 'ploss_avg_3d', ploss_avg_3d(:,:,:) )
374	CALL iom_get( numrtr, jpdom_autoglo, 'phyt_avg_3d', phyt_avg_3d(:,:,:) )
375	ELSE
376	IF(lwp) WRITE(numout,*) ' MEDUSA pgrow_avg_3d absent - setting to zero ...'
377	pgrow_avg_3d(:,:,:) = 0.0
378	ploss_avg_3d(:,:,:) = 0.0
379	phyt_avg_3d(:,:,:) = 0.0
380	ENDIF
381	ENDIF
382
383	#endif
384	!
385	#if defined key_idtra
386	!! JPALM -- 05-01-2016 -- mv idtra and CFC restart reading and
387	!! writting here undre their key.
388	!! problems in CFC restart, maybe because of this...
389	!! and pb in idtra diag or diad-restart writing.
390	!!----------------------------------------------------------------------
391	IF( iom_varid( numrtr, 'qint_IDTRA', ldstop = .FALSE. ) > 0 ) THEN
392	!! YES; in which case read them
393	!!
394	IF(lwp) WRITE(numout,*) ' IDTRA averaged properties present - reading in ...'
395	CALL iom_get( numrtr, jpdom_autoglo, 'qint_IDTRA', qint_idtra(:,:,1) )
396	ELSE
397	!! NO; in which case set them to zero
398	!!
399	IF(lwp) WRITE(numout,*) ' IDTRA averaged properties absent - setting to zero ...'
400	qint_idtra(:,:,1) = 0.0 !! CHN
401	ENDIF
402	!!
403	!! calculate stats on these fields
404	IF(lwp) WRITE(numout,*) ' IDTRA averaged properties stats (min, max, sum) ...'
405	call trc_rst_dia_stat(qint_idtra(:,:,1), 'qint_IDTRA')
406	#endif
407	!
408	#if defined key_cfc
409	DO jl = 1, jp_cfc
410	jn = jp_cfc0 + jl - 1
411	IF( iom_varid( numrtr, 'qint_'//ctrcnm(jn), ldstop = .FALSE. ) > 0 ) THEN
412	!! YES; in which case read them
413	!!
414	IF(lwp) WRITE(numout,*) ' CFC averaged properties present - reading in ...'
415	CALL iom_get( numrtr, jpdom_autoglo, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jl) )
416	ELSE
417	!! NO; in which case set them to zero
418	!!
419	IF(lwp) WRITE(numout,*) ' CFC averaged properties absent - setting to zero ...'
420	qint_cfc(:,:,jl) = 0.0 !! CHN
421	ENDIF
422	!!
423	!! calculate stats on these fields
424	IF(lwp) WRITE(numout,*) ' CFC averaged properties stats (min, max, sum) ...'
425	call trc_rst_dia_stat(qint_cfc(:,:,jl), 'qint_'//ctrcnm(jn))
426	END DO
427	#endif
428	!
429	END SUBROUTINE trc_rst_read
430
431	SUBROUTINE trc_rst_wri( kt )
432	!!----------------------------------------------------------------------
433	!! * trc_rst_wri *
434	!!
435	!! ** purpose : write passive tracer fields in restart files
436	!!----------------------------------------------------------------------
437	INTEGER, INTENT( in ) :: kt ! ocean time-step index
438	!!
439	INTEGER :: jn, jl
440	REAL(wp) :: zarak0
441	!! AXY (05/11/13): temporary variables
442	REAL(wp) :: fq0,fq1,fq2
443	!!----------------------------------------------------------------------
444	!
445	CALL iom_rstput( kt, nitrst, numrtw, 'rdttrc1', rdttrc(1) ) ! surface passive tracer time step
446	! prognostic variables
447	! --------------------
448	DO jn = 1, jptra
449	CALL iom_rstput( kt, nitrst, numrtw, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
450	END DO
451
452	DO jn = 1, jptra
453	CALL iom_rstput( kt, nitrst, numrtw, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
454	END DO
455
456	!! AXY (09/06/14): the ARCHER version of MEDUSA-2 does not include an equivalent
457	!! call to MEDUSA-2 at this point; this suggests that the FCM
458	!! version of NEMO date significantly earlier than the current
459	!! version
460
461	#if defined key_medusa
462	!! AXY (13/01/12): write out "before" and "now" state of seafloor
463	!! sediment pools into restart; this happens
464	!! whether or not the pools are to be used by
465	!! MEDUSA (which is controlled by a switch in the
466	!! namelist_top file)
467	!!
468	IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields - writing out ...'
469	CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_N', zb_sed_n(:,:) )
470	CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_N', zn_sed_n(:,:) )
471	CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_FE', zb_sed_fe(:,:) )
472	CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_FE', zn_sed_fe(:,:) )
473	CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_SI', zb_sed_si(:,:) )
474	CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_SI', zn_sed_si(:,:) )
475	CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_C', zb_sed_c(:,:) )
476	CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_C', zn_sed_c(:,:) )
477	CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_CA', zb_sed_ca(:,:) )
478	CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_CA', zn_sed_ca(:,:) )
479	!!
480	!! calculate stats on these fields
481	IF(lwp) WRITE(numout,*) ' MEDUSA sediment field stats (min, max, sum) ...'
482	call trc_rst_dia_stat(zn_sed_n(:,:), 'Sediment N')
483	call trc_rst_dia_stat(zn_sed_fe(:,:), 'Sediment Fe')
484	call trc_rst_dia_stat(zn_sed_si(:,:), 'Sediment Si')
485	call trc_rst_dia_stat(zn_sed_c(:,:), 'Sediment C')
486	call trc_rst_dia_stat(zn_sed_ca(:,:), 'Sediment Ca')
487	!!
488	!! AXY (07/07/15): write out temporally averaged fields for DMS
489	!! calculations
490	!!
491	IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS - writing out ...'
492	CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_CHN', zb_dms_chn(:,:) )
493	CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_CHN', zn_dms_chn(:,:) )
494	CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_CHD', zb_dms_chd(:,:) )
495	CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_CHD', zn_dms_chd(:,:) )
496	CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_MLD', zb_dms_mld(:,:) )
497	CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_MLD', zn_dms_mld(:,:) )
498	CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_QSR', zb_dms_qsr(:,:) )
499	CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_QSR', zn_dms_qsr(:,:) )
500	CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_DIN', zb_dms_din(:,:) )
501	CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_DIN', zn_dms_din(:,:) )
502	!! JPALM 14-06-2016 -- add CO2 flux and DMS surf through the restart
503	!! -- needed for the coupling with atm
504	CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_srf', zb_dms_srf(:,:) )
505	CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_srf', zn_dms_srf(:,:) )
506	CALL iom_rstput( kt, nitrst, numrtw, 'B_CO2_flx', zb_co2_flx(:,:) )
507	CALL iom_rstput( kt, nitrst, numrtw, 'N_CO2_flx', zn_co2_flx(:,:) )
508	!! JPALM 07-12-2017 -- To make things cleaner, we want to store an
509	!! unscaled Chl field in the restart and only
510	!! scale it when reading it in.
511	CALL iom_rstput( kt, nitrst, numrtw, 'N_CHL_srf', zn_chl_srf(:,:) )
512	!!
513	!! calculate stats on these fields
514	IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS stats (min, max, sum) ...'
515	call trc_rst_dia_stat(zn_dms_chn(:,:), 'DMS, CHN')
516	call trc_rst_dia_stat(zn_dms_chd(:,:), 'DMS, CHD')
517	call trc_rst_dia_stat(zn_dms_mld(:,:), 'DMS, MLD')
518	call trc_rst_dia_stat(zn_dms_qsr(:,:), 'DMS, QSR')
519	call trc_rst_dia_stat(zn_dms_din(:,:), 'DMS, DIN')
520	call trc_rst_dia_stat(zn_dms_srf(:,:), 'DMS surf')
521	call trc_rst_dia_stat(zn_co2_flx(:,:), 'CO2 flux')
522	call trc_rst_dia_stat(zn_chl_srf(:,:), 'unscaled CHL cpl')
523	!!
524	IF(lwp) WRITE(numout,*) ' MEDUSA averaged prop. for dust and iron dep.'
525	call trc_rst_dia_stat(dust(:,:), 'Dust dep')
526	call trc_rst_dia_stat(zirondep(:,:), 'Iron dep')
527	!!
528	!!
529	!! JPALM 14-06-2016 -- add Carbonate chenistry variables through the restart
530	!! -- needed for monthly call of carb-chem routine and better reproducibility
531	# if defined key_roam
532	IF(lwp) WRITE(numout,*) 'Carbonate chem variable - writing out ...'
533	CALL iom_rstput( kt, nitrst, numrtw, 'pH_3D' , f3_pH(:,:,:) )
534	CALL iom_rstput( kt, nitrst, numrtw, 'h2CO3_3D', f3_h2co3(:,:,:) )
535	CALL iom_rstput( kt, nitrst, numrtw, 'hCO3_3D' , f3_hco3(:,:,:) )
536	CALL iom_rstput( kt, nitrst, numrtw, 'CO3_3D' , f3_co3(:,:,:) )
537	CALL iom_rstput( kt, nitrst, numrtw, 'omcal_3D', f3_omcal(:,:,:) )
538	CALL iom_rstput( kt, nitrst, numrtw, 'omarg_3D', f3_omarg(:,:,:) )
539	CALL iom_rstput( kt, nitrst, numrtw, 'CCD_CAL' , f2_ccd_cal(:,:) )
540	CALL iom_rstput( kt, nitrst, numrtw, 'CCD_ARG' , f2_ccd_arg(:,:) )
541	!!
542	IF(lwp) WRITE(numout,*) ' MEDUSA averaged Carb-chem stats (min, max, sum) ...'
543	call trc_rst_dia_stat( f3_pH(:,:,1) ,'pH 3D surf')
544	call trc_rst_dia_stat( f3_h2co3(:,:,1),'h2CO3 3D surf')
545	call trc_rst_dia_stat( f3_hco3(:,:,1) ,'hCO3 3D surf' )
546	call trc_rst_dia_stat( f3_co3(:,:,1) ,'CO3 3D surf' )
547	call trc_rst_dia_stat( f3_omcal(:,:,1),'omcal 3D surf')
548	call trc_rst_dia_stat( f3_omarg(:,:,1),'omarg 3D surf')
549	call trc_rst_dia_stat( f2_ccd_cal(:,:),'CCD_CAL')
550	call trc_rst_dia_stat( f2_ccd_arg(:,:),'CCD_ARG')
551	!!
552	IF ( ln_foam_medusa ) THEN
553	!! Fields for observation operator on first timestep
554	IF(lwp) WRITE(numout,*) ' MEDUSA OBS fields - writing out ...'
555	CALL iom_rstput( kt, nitrst, numrtw, 'PCO2W', f2_pco2w(:,:) )
556	CALL iom_rstput( kt, nitrst, numrtw, 'FCO2W', f2_fco2w(:,:) )
557	ENDIF
558	# endif
559	IF ( ln_foam_medusa ) THEN
560	!! Fields for assimilation on first timestep
561	IF(lwp) WRITE(numout,*) ' MEDUSA ASM fields - writing out ...'
562	CALL iom_rstput( kt, nitrst, numrtw, 'pgrow_avg', pgrow_avg(:,:) )
563	CALL iom_rstput( kt, nitrst, numrtw, 'ploss_avg', ploss_avg(:,:) )
564	CALL iom_rstput( kt, nitrst, numrtw, 'phyt_avg', phyt_avg(:,:) )
565	CALL iom_rstput( kt, nitrst, numrtw, 'mld_max', mld_max(:,:) )
566	CALL iom_rstput( kt, nitrst, numrtw, 'pgrow_avg_3d', pgrow_avg_3d(:,:,:) )
567	CALL iom_rstput( kt, nitrst, numrtw, 'ploss_avg_3d', ploss_avg_3d(:,:,:) )
568	CALL iom_rstput( kt, nitrst, numrtw, 'phyt_avg_3d', phyt_avg_3d(:,:,:) )
569	ENDIF
570	!!
571	#endif
572	!
573	#if defined key_idtra
574	!! JPALM -- 05-01-2016 -- mv idtra and CFC restart reading and
575	!! writting here undre their key.
576	!! problems in CFC restart, maybe because of this...
577	!! and pb in idtra diag or diad-restart writing.
578	!!----------------------------------------------------------------------
579	IF(lwp) WRITE(numout,*) ' IDTRA averaged properties - writing out ...'
580	CALL iom_rstput( kt, nitrst, numrtw, 'qint_IDTRA', qint_idtra(:,:,1) )
581	!!
582	!! calculate stats on these fields
583	IF(lwp) WRITE(numout,*) ' IDTRA averaged properties stats (min, max, sum) ...'
584	call trc_rst_dia_stat(qint_idtra(:,:,1), 'qint_IDTRA')
585	#endif
586	!
587	#if defined key_cfc
588	DO jl = 1, jp_cfc
589	jn = jp_cfc0 + jl - 1
590	IF(lwp) WRITE(numout,*) ' CFC averaged properties - writing out ...'
591	CALL iom_rstput( kt, nitrst, numrtw, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jl) )
592	!!
593	!! calculate stats on these fields
594	IF(lwp) WRITE(numout,*) ' CFC averaged properties stats (min, max, sum) ...'
595	call trc_rst_dia_stat(qint_cfc(:,:,jl), 'qint_'//ctrcnm(jn))
596	END DO
597	#endif
598	!
599
600	IF( kt == nitrst ) THEN
601	CALL trc_rst_stat ! statistics
602	#if defined key_medusa && defined key_roam
603	CALL trc_rst_conserve ! conservation check
604	#endif
605	CALL iom_close( numrtw ) ! close the restart file (only at last time step)
606	#if ! defined key_trdmxl_trc
607	lrst_trc = .FALSE.
608	#endif
609	IF( lk_offline .AND. ln_rst_list ) THEN
610	nrst_lst = nrst_lst + 1
611	nitrst = nstocklist( nrst_lst )
612	ENDIF
613	ENDIF
614	!
615	END SUBROUTINE trc_rst_wri
616
617
618	SUBROUTINE trc_rst_cal( kt, cdrw )
619	!!---------------------------------------------------------------------
620	!! * ROUTINE trc_rst_cal *
621	!!
622	!! ** Purpose : Read or write calendar in restart file:
623	!!
624	!! WRITE(READ) mode:
625	!! kt : number of time step since the begining of the experiment at the
626	!! end of the current(previous) run
627	!! adatrj(0) : number of elapsed days since the begining of the experiment at the
628	!! end of the current(previous) run (REAL -> keep fractions of day)
629	!! ndastp : date at the end of the current(previous) run (coded as yyyymmdd integer)
630	!!
631	!! According to namelist parameter nrstdt,
632	!! nn_rsttr = 0 no control on the date (nittrc000 is arbitrary).
633	!! nn_rsttr = 1 we verify that nittrc000 is equal to the last
634	!! time step of previous run + 1.
635	!! In both those options, the exact duration of the experiment
636	!! since the beginning (cumulated duration of all previous restart runs)
637	!! is not stored in the restart and is assumed to be (nittrc000-1)*rdt.
638	!! This is valid is the time step has remained constant.
639	!!
640	!! nn_rsttr = 2 the duration of the experiment in days (adatrj)
641	!! has been stored in the restart file.
642	!!----------------------------------------------------------------------
643	INTEGER , INTENT(in) :: kt ! ocean time-step
644	CHARACTER(len=*), INTENT(in) :: cdrw ! "READ"/"WRITE" flag
645	!
646	INTEGER :: jlibalt = jprstlib
647	LOGICAL :: llok
648	REAL(wp) :: zkt, zrdttrc1
649	REAL(wp) :: zndastp
650
651	! Time domain : restart
652	! ---------------------
653
654	IF( TRIM(cdrw) == 'READ' ) THEN
655
656	IF(lwp) WRITE(numout,*)
657	IF(lwp) WRITE(numout,*) 'trc_rst_cal : read the TOP restart file for calendar'
658	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
659
660	IF ( jprstlib == jprstdimg ) THEN
661	! eventually read netcdf file (monobloc) for restarting on different number of processors
662	! if {cn_trcrst_in}.nc exists, then set jlibalt to jpnf90
663	INQUIRE( FILE = TRIM(cn_trcrst_indir)//'/'//TRIM(cn_trcrst_in)//'.nc', EXIST = llok )
664	IF ( llok ) THEN ; jlibalt = jpnf90 ; ELSE ; jlibalt = jprstlib ; ENDIF
665	ENDIF
666
667	IF( ln_rsttr ) THEN
668	CALL iom_open( TRIM(cn_trcrst_indir)//'/'//cn_trcrst_in, numrtr, kiolib = jlibalt )
669	CALL iom_get ( numrtr, 'kt', zkt ) ! last time-step of previous run
670
671	IF(lwp) THEN
672	WRITE(numout,) ' ** Info read in restart : '
673	WRITE(numout,*) ' previous time-step : ', NINT( zkt )
674	WRITE(numout,) ' ** restart option'
675	SELECT CASE ( nn_rsttr )
676	CASE ( 0 ) ; WRITE(numout,*) ' nn_rsttr = 0 : no control of nittrc000'
677	CASE ( 1 ) ; WRITE(numout,*) ' nn_rsttr = 1 : no control the date at nittrc000 (use ndate0 read in the namelist)'
678	CASE ( 2 ) ; WRITE(numout,*) ' nn_rsttr = 2 : calendar parameters read in restart'
679	END SELECT
680	WRITE(numout,*)
681	ENDIF
682	! Control of date
683	IF( nittrc000 - NINT( zkt ) /= nn_dttrc .AND. nn_rsttr /= 0 ) &
684	& CALL ctl_stop( ' ===>>>> : problem with nittrc000 for the restart', &
685	& ' verify the restart file or rerun with nn_rsttr = 0 (namelist)' )
686	ENDIF
687	!
688	IF( lk_offline ) THEN
689	! ! set the date in offline mode
690	IF( ln_rsttr .AND. nn_rsttr == 2 ) THEN
691	CALL iom_get( numrtr, 'ndastp', zndastp )
692	ndastp = NINT( zndastp )
693	CALL iom_get( numrtr, 'adatrj', adatrj )
694	ELSE
695	ndastp = ndate0 - 1 ! ndate0 read in the namelist in dom_nam
696	adatrj = ( REAL( nittrc000-1, wp ) * rdttra(1) ) / rday
697	! note this is wrong if time step has changed during run
698	ENDIF
699	!
700	IF(lwp) THEN
701	WRITE(numout,) ' ** Info used values : '
702	WRITE(numout,*) ' date ndastp : ', ndastp
703	WRITE(numout,*) ' number of elapsed days since the begining of run : ', adatrj
704	WRITE(numout,*)
705	ENDIF
706	!
707	IF( ln_rsttr ) THEN ; neuler = 1
708	ELSE ; neuler = 0
709	ENDIF
710	!
711	CALL day_init ! compute calendar
712	!
713	ENDIF
714	!
715	ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
716	!
717	IF( kt == nitrst ) THEN
718	IF(lwp) WRITE(numout,*)
719	IF(lwp) WRITE(numout,*) 'trc_wri : write the TOP restart file (NetCDF) at it= ', kt, ' date= ', ndastp
720	IF(lwp) WRITE(numout,*) '~~~~~~~'
721	ENDIF
722	CALL iom_rstput( kt, nitrst, numrtw, 'kt' , REAL( kt , wp) ) ! time-step
723	CALL iom_rstput( kt, nitrst, numrtw, 'ndastp' , REAL( ndastp, wp) ) ! date
724	CALL iom_rstput( kt, nitrst, numrtw, 'adatrj' , adatrj ) ! number of elapsed days since
725	! ! the begining of the run [s]
726	ENDIF
727
728	END SUBROUTINE trc_rst_cal
729
730
731	SUBROUTINE trc_rst_stat
732	!!----------------------------------------------------------------------
733	!! * trc_rst_stat *
734	!!
735	!! ** purpose : Compute tracers statistics
736	!!----------------------------------------------------------------------
737	INTEGER :: jk, jn
738	REAL(wp) :: ztraf, zmin, zmax, zmean, zdrift
739	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvol
740	!!----------------------------------------------------------------------
741
742	IF( lwp ) THEN
743	WRITE(numout,*)
744	WRITE(numout,*) ' ----TRACER STAT---- '
745	WRITE(numout,*)
746	ENDIF
747	!
748	DO jk = 1, jpk
749	zvol(:,:,jk) = e1e2t(:,:) * fse3t_a(:,:,jk) * tmask(:,:,jk)
750	END DO
751	!
752	DO jn = 1, jptra
753	ztraf = glob_sum( trn(:,:,:,jn) * zvol(:,:,:) )
754	zmin = MINVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
755	zmax = MAXVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
756	IF( lk_mpp ) THEN
757	CALL mpp_min( zmin ) ! min over the global domain
758	CALL mpp_max( zmax ) ! max over the global domain
759	END IF
760	zmean = ztraf / areatot
761	zdrift = ( ( ztraf - trai(jn) ) / ( trai(jn) + 1.e-12 ) ) * 100._wp
762	IF(lwp) WRITE(numout,9000) jn, TRIM( ctrcnm(jn) ), zmean, zmin, zmax, zdrift
763	END DO
764	IF(lwp) WRITE(numout,*)
765	9000 FORMAT(' tracer nb :',i2,' name :',a10,' mean :',e18.10,' min :',e18.10, &
766	& ' max :',e18.10,' drift :',e18.10, ' %')
767	!
768	END SUBROUTINE trc_rst_stat
769
770
771	# if defined key_medusa && defined key_roam
772	SUBROUTINE trc_rst_conserve
773	!!----------------------------------------------------------------------
774	!! * trc_rst_conserve *
775	!!
776	!! ** purpose : Compute tracers conservation statistics
777	!!
778	!! AXY (17/11/2017)
779	!! This routine calculates the "now" inventories of the elemental
780	!! cycles of MEDUSA and compares them to those calculate when the
781	!! model was initialised / restarted; the cycles calculated are:
782	!! nitrogen, silicon, iron, carbon, alkalinity and oxygen
783	!!----------------------------------------------------------------------
784	INTEGER :: ji, jj, jk, jn
785	REAL(wp) :: zsum3d, zsum2d, zinvt, zdelta, zratio
786	REAL(wp), DIMENSION(jpi,jpj,jpk) :: z3d, zvol
787	REAL(wp), DIMENSION(jpi,jpj) :: z2d, zarea
788	REAL(wp), DIMENSION(6) :: loc_cycletot3, loc_cycletot2
789	!!----------------------------------------------------------------------
790	!
791	IF( lwp ) THEN
792	WRITE(numout,*)
793	WRITE(numout,*) ' ----TRACER CONSERVATION---- '
794	WRITE(numout,*)
795	ENDIF
796	!
797	! ocean volume
798	DO jk = 1, jpk
799	zvol(:,:,jk) = e1e2t(:,:) * fse3t_a(:,:,jk) * tmask(:,:,jk)
800	END DO
801	!
802	! ocean area (for sediments)
803	zarea(:,:) = e1e2t(:,:) * tmask(:,:,1)
804	!
805	!----------------------------------------------------------------------
806	! nitrogen
807	z3d(:,:,:) = trn(:,:,:,jpphn) + trn(:,:,:,jpphd) + trn(:,:,:,jpzmi) + &
808	trn(:,:,:,jpzme) + trn(:,:,:,jpdet) + trn(:,:,:,jpdin)
809	z2d(:,:) = zn_sed_n(:,:)
810	zsum3d = glob_sum( z3d(:,:,:) * zvol(:,:,:) )
811	zsum2d = glob_sum( z2d(:,:) * zarea(:,:) )
812	! total tracer, and delta
813	zinvt = zsum3d + zsum2d
814	zdelta = zinvt - cycletot(1)
815	zratio = 1.0e2 * zdelta / cycletot(1)
816	!
817	IF ( lwp ) WRITE(numout,9010) 'nitrogen', zsum3d, zsum2d, zinvt, &
818	cycletot(1), zdelta, zratio
819	IF ( lwp ) WRITE(numout,*)
820	!
821	!----------------------------------------------------------------------
822	! silicon
823	z3d(:,:,:) = trn(:,:,:,jppds) + trn(:,:,:,jpsil)
824	z2d(:,:) = zn_sed_si(:,:)
825	zsum3d = glob_sum( z3d(:,:,:) * zvol(:,:,:) )
826	zsum2d = glob_sum( z2d(:,:) * zarea(:,:) )
827	! total tracer, and delta
828	zinvt = zsum3d + zsum2d
829	zdelta = zinvt - cycletot(2)
830	zratio = 1.0e2 * zdelta / cycletot(2)
831	!
832	IF ( lwp ) WRITE(numout,9010) 'silicon', zsum3d, zsum2d, zinvt, &
833	cycletot(2), zdelta, zratio
834	IF ( lwp ) WRITE(numout,*)
835	!
836	!----------------------------------------------------------------------
837	! iron
838	z3d(:,:,:) = ((trn(:,:,:,jpphn) + trn(:,:,:,jpphd) + trn(:,:,:,jpzmi) + &
839	trn(:,:,:,jpzme) + trn(:,:,:,jpdet)) * xrfn) + trn(:,:,:,jpfer)
840	z2d(:,:) = zn_sed_fe(:,:)
841	zsum3d = glob_sum( z3d(:,:,:) * zvol(:,:,:) )
842	zsum2d = glob_sum( z2d(:,:) * zarea(:,:) )
843	! total tracer, and delta
844	zinvt = zsum3d + zsum2d
845	zdelta = zinvt - cycletot(3)
846	zratio = 1.0e2 * zdelta / cycletot(3)
847	!
848	IF ( lwp ) WRITE(numout,9010) 'iron', zsum3d, zsum2d, zinvt, &
849	cycletot(3), zdelta, zratio
850	IF ( lwp ) WRITE(numout,*)
851	!
852	!----------------------------------------------------------------------
853	! carbon
854	z3d(:,:,:) = (trn(:,:,:,jpphn) * xthetapn) + (trn(:,:,:,jpphd) * xthetapd) + &
855	(trn(:,:,:,jpzmi) * xthetazmi) + (trn(:,:,:,jpzme) * xthetazme) + &
856	trn(:,:,:,jpdtc) + trn(:,:,:,jpdic)
857	z2d(:,:) = zn_sed_c(:,:) + zn_sed_ca(:,:)
858	zsum3d = glob_sum( z3d(:,:,:) * zvol(:,:,:) )
859	zsum2d = glob_sum( z2d(:,:) * zarea(:,:) )
860	! total tracer, and delta
861	zinvt = zsum3d + zsum2d
862	zdelta = zinvt - cycletot(4)
863	zratio = 1.0e2 * zdelta / cycletot(4)
864	!
865	IF ( lwp ) WRITE(numout,9010) 'carbon', zsum3d, zsum2d, zinvt, &
866	cycletot(4), zdelta, zratio
867	IF ( lwp ) WRITE(numout,*)
868	!
869	!----------------------------------------------------------------------
870	! alkalinity
871	z3d(:,:,:) = trn(:,:,:,jpalk)
872	z2d(:,:) = zn_sed_ca(:,:) * 2.0
873	zsum3d = glob_sum( z3d(:,:,:) * zvol(:,:,:) )
874	zsum2d = glob_sum( z2d(:,:) * zarea(:,:) )
875	! total tracer, and delta
876	zinvt = zsum3d + zsum2d
877	zdelta = zinvt - cycletot(5)
878	zratio = 1.0e2 * zdelta / cycletot(5)
879	!
880	IF ( lwp ) WRITE(numout,9010) 'alkalinity', zsum3d, zsum2d, zinvt, &
881	cycletot(5), zdelta, zratio
882	IF ( lwp ) WRITE(numout,*)
883	!
884	!----------------------------------------------------------------------
885	! oxygen
886	z3d(:,:,:) = trn(:,:,:,jpoxy)
887	z2d(:,:) = 0.0
888	zsum3d = glob_sum( z3d(:,:,:) * zvol(:,:,:) )
889	zsum2d = glob_sum( z2d(:,:) * zarea(:,:) )
890	! total tracer, and delta
891	zinvt = zsum3d + zsum2d
892	zdelta = zinvt - cycletot(6)
893	zratio = 1.0e2 * zdelta / cycletot(6)
894	!
895	IF ( lwp ) WRITE(numout,9010) 'oxygen', zsum3d, zsum2d, zinvt, &
896	cycletot(6), zdelta, zratio
897	!
898	!----------------------------------------------------------------------
899	! Check
900	zsum3d = glob_sum( zvol(:,:,:) )
901	zsum2d = glob_sum( zarea(:,:) )
902	IF ( lwp ) THEN
903	WRITE(numout,*)
904	WRITE(numout,*) ' check : cvol : ', zsum3d
905	WRITE(numout,*) ' check : carea : ', zsum2d
906	WRITE(numout,*)
907	ENDIF
908	!
909	9010 FORMAT(' element:',a10, &
910	' 3d sum:',e18.10,' 2d sum:',e18.10, &
911	' total:',e18.10,' initial:',e18.10, &
912	' delta:',e18.10,' %:',e18.10)
913	!
914	END SUBROUTINE trc_rst_conserve
915	# endif
916
917
918	#else
919	!!----------------------------------------------------------------------
920	!! Dummy module : No passive tracer
921	!!----------------------------------------------------------------------
922	CONTAINS
923	SUBROUTINE trc_rst_read ! Empty routines
924	END SUBROUTINE trc_rst_read
925	SUBROUTINE trc_rst_wri( kt )
926	INTEGER, INTENT ( in ) :: kt
927	WRITE(,) 'trc_rst_wri: You should not have seen this print! error?', kt
928	END SUBROUTINE trc_rst_wri
929	#endif
930
931	!!----------------------------------------------------------------------
932	!! NEMO/TOP 3.3 , NEMO Consortium (2010)
933	!! $Id$
934	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
935	!!======================================================================
936	END MODULE trcrst

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

Download in other formats: