Context Navigation

trdmxl.F90 @ 10774

Last change on this file since 10774 was 10774, checked in by andmirek, 5 years ago
GMED 450 add flush after prints
File size: 44.5 KB

Line
1	MODULE trdmxl
2	!!======================================================================
3	!! * MODULE trdmxl *
4	!! Ocean diagnostics: mixed layer T-S trends
5	!!======================================================================
6	!! History : OPA ! 1995-04 (J. Vialard) Original code
7	!! ! 1997-02 (E. Guilyardi) Adaptation global + base cmo
8	!! ! 1999-09 (E. Guilyardi) Re-writing + netCDF output
9	!! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module
10	!! - ! 2004-08 (C. Talandier) New trends organization
11	!! 2.0 ! 2005-05 (C. Deltel) Diagnose trends of time averaged ML T & S
12	!! 3.5 ! 2012-03 (G. Madec) complete reorganisation + change in the time averaging
13	!!----------------------------------------------------------------------
14
15	!!----------------------------------------------------------------------
16	!! trd_mxl : T and S cumulated trends averaged over the mixed layer
17	!! trd_mxl_zint : T and S trends vertical integration
18	!! trd_mxl_init : initialization step
19	!!----------------------------------------------------------------------
20	USE oce ! ocean dynamics and tracers variables
21	USE dom_oce ! ocean space and time domain variables
22	USE trd_oce ! trends: ocean variables
23	USE trdmxl_oce ! ocean variables trends
24	USE ldftra_oce ! ocean active tracers lateral physics
25	USE zdf_oce ! ocean vertical physics
26	USE in_out_manager ! I/O manager
27	USE phycst ! Define parameters for the routines
28	USE dianam ! build the name of file (routine)
29	USE ldfslp ! iso-neutral slopes
30	USE zdfmxl ! mixed layer depth
31	USE zdfddm ! ocean vertical physics: double diffusion
32	USE ioipsl ! NetCDF library
33	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
34	USE diadimg ! dimg direct access file format output
35	USE trdmxl_rst ! restart for diagnosing the ML trends
36	USE prtctl ! Print control
37	USE restart ! for lrst_oce
38	USE lib_mpp ! MPP library
39	USE wrk_nemo ! Memory allocation
40	USE iom
41
42	IMPLICIT NONE
43	PRIVATE
44
45	PUBLIC trd_mxl ! routine called by step.F90
46	PUBLIC trd_mxl_init ! routine called by opa.F90
47	PUBLIC trd_mxl_zint ! routine called by tracers routines
48
49	INTEGER :: nkstp ! current time step
50
51
52
53	!!gm to be moved from trdmxl_oce
54	! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: hml ! ML depth (sum of e3t over nmln-1 levels) [m]
55	! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: tml , sml ! now ML averaged T & S
56	! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: tmlb_nf, smlb_nf ! not filtered before ML averaged T & S
57	!
58	!
59	! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: hmlb, hmln ! before, now, and after Mixed Layer depths [m]
60	!
61	! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: tb_mlb, tb_mln ! before (not filtered) tracer averaged over before and now ML
62	!
63	! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: tn_mln ! now tracer averaged over now ML
64	!!gm end
65
66	CHARACTER (LEN=40) :: clhstnam ! name of the trends NetCDF file
67	INTEGER :: nh_t, nmoymltrd
68	INTEGER :: nidtrd
69	INTEGER, ALLOCATABLE, SAVE, DIMENSION(:) :: ndextrd1
70	INTEGER :: ndimtrd1
71	INTEGER :: ionce, icount
72
73	!! * Substitutions
74	# include "domzgr_substitute.h90"
75	# include "ldftra_substitute.h90"
76	# include "zdfddm_substitute.h90"
77	!!----------------------------------------------------------------------
78	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
79	!! $Id$
80	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
81	!!----------------------------------------------------------------------
82	CONTAINS
83
84	INTEGER FUNCTION trd_mxl_alloc()
85	!!----------------------------------------------------------------------
86	!! * ROUTINE trd_mxl_alloc *
87	!!----------------------------------------------------------------------
88	ALLOCATE( ndextrd1(jpi*jpj) , STAT=trd_mxl_alloc )
89	!
90	IF( lk_mpp ) CALL mpp_sum ( trd_mxl_alloc )
91	IF( trd_mxl_alloc /= 0 ) CALL ctl_warn('trd_mxl_alloc: failed to allocate array ndextrd1')
92	END FUNCTION trd_mxl_alloc
93
94
95	SUBROUTINE trd_tra_mxl( ptrdx, ptrdy, ktrd, kt, p2dt, kmxln )
96	!!----------------------------------------------------------------------
97	!! * ROUTINE trd_tra_mng *
98	!!
99	!! ** Purpose : Dispatch all trends computation, e.g. 3D output, integral
100	!! constraints, barotropic vorticity, kinetic enrgy,
101	!! potential energy, and/or mixed layer budget.
102	!!----------------------------------------------------------------------
103	REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: ptrdx ! Temperature or U trend
104	REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: ptrdy ! Salinity or V trend
105	INTEGER , INTENT(in ) :: ktrd ! tracer trend index
106	INTEGER , INTENT(in ) :: kt ! time step index
107	REAL(wp) , INTENT(in ) :: p2dt ! time step [s]
108	REAL(wp), DIMENSION(:,:) , INTENT(in ) :: kmxln ! number of t-box for the vertical average
109	!
110	INTEGER :: ji, jj, jk ! dummy loop indices
111	!!----------------------------------------------------------------------
112
113	! !==============================!
114	IF ( kt /= nkstp ) THEN != 1st call at kt time step =!
115	! !==============================!
116	nkstp = kt
117
118
119	! !== reset trend arrays to zero ==!
120	tmltrd(:,:,:) = 0._wp ; smltrd(:,:,:) = 0._wp
121
122
123	!
124	wkx(:,:,:) = 0._wp !== now ML weights for vertical averaging ==!
125	DO jk = 1, jpktrd ! initialize wkx with vertical scale factor in mixed-layer
126	DO jj = 1,jpj
127	DO ji = 1,jpi
128	IF( jk - kmxln(ji,jj) < 0 ) wkx(ji,jj,jk) = fse3t(ji,jj,jk) * tmask(ji,jj,jk)
129	END DO
130	END DO
131	END DO
132	hmxl(:,:) = 0._wp ! NOW mixed-layer depth
133	DO jk = 1, jpktrd
134	hmxl(:,:) = hmxl(:,:) + wkx(:,:,jk)
135	END DO
136	DO jk = 1, jpktrd ! integration weights
137	wkx(:,:,jk) = wkx(:,:,jk) / MAX( 1.e-20_wp, hmxl(:,:) ) * tmask(:,:,1)
138	END DO
139
140
141	!
142	! !== Vertically averaged T and S ==!
143	tml(:,:) = 0._wp ; sml(:,:) = 0._wp
144	DO jk = 1, jpktrd
145	tml(:,:) = tml(:,:) + wkx(:,:,jk) * tsn(:,:,jk,jp_tem)
146	sml(:,:) = sml(:,:) + wkx(:,:,jk) * tsn(:,:,jk,jp_sal)
147	END DO
148	!
149	ENDIF
150
151
152
153	! mean now trends over the now ML
154	tmltrd(:,:,ktrd) = tmltrd(:,:,ktrd) + ptrdx(:,:,jk) * wkx(:,:,jk) ! temperature
155	smltrd(:,:,ktrd) = smltrd(:,:,ktrd) + ptrdy(:,:,jk) * wkx(:,:,jk) ! salinity
156
157
158
159	!!gm to be put juste before the output !
160	! ! Lateral boundary conditions
161	! CALL lbc_lnk( tmltrd(:,:,jl), 'T', 1. )
162	! CALL lbc_lnk( smltrd(:,:,jl), 'T', 1. )
163	!!gm end
164
165
166
167	SELECT CASE( ktrd )
168	CASE( jptra_npc ) ! non-penetrative convection: regrouped with zdf
169	!!gm : to be completed !
170	! IF( ....
171	!!gm end
172	CASE( jptra_zdfp ) ! iso-neutral diffusion: "pure" vertical diffusion
173	! ! regroup iso-neutral diffusion in one term
174	tmltrd(:,:,jpmxl_ldf) = tmltrd(:,:,jpmxl_ldf) + ( tmltrd(:,:,jpmxl_zdf) - tmltrd(:,:,jpmxl_zdfp) )
175	smltrd(:,:,jpmxl_ldf) = smltrd(:,:,jpmxl_ldf) + ( smltrd(:,:,jpmxl_zdf) - smltrd(:,:,jpmxl_zdfp) )
176	! ! put in zdf the dia-neutral diffusion
177	tmltrd(:,:,jpmxl_zdf) = tmltrd(:,:,jpmxl_zdfp)
178	smltrd(:,:,jpmxl_zdf) = smltrd(:,:,jpmxl_zdfp)
179	IF( ln_zdfnpc ) THEN
180	tmltrd(:,:,jpmxl_zdf) = tmltrd(:,:,jpmxl_zdf) + tmltrd(:,:,jpmxl_npc)
181	smltrd(:,:,jpmxl_zdf) = smltrd(:,:,jpmxl_zdf) + smltrd(:,:,jpmxl_npc)
182	ENDIF
183	!
184	CASE( jptra_atf ) ! last trends of the current time step: perform the time averaging & output
185	!
186	! after ML : zhmla NB will be swaped to provide hmln and hmlb
187	!
188	! entrainement ent_1 : tb_mln - tb_mlb ==>> use previous timestep ztn_mla = tb_mln
189	! " " " tn_mln = tb_mlb (unfiltered tb!)
190	! NB: tn_mln itself comes from the 2 time step before (ta_mla)
191	!
192	! atf trend : ztbf_mln - tb_mln ==>> use previous timestep tn_mla = tb_mln
193	! need to compute tbf_mln, using the current tb
194	! which is the before fitered tracer
195	!
196	! entrainement ent_2 : zta_mla - zta_mln ==>> need to compute zta_mla and zta_mln
197	!
198	! time averaging : mean: CALL trd_mean( kt, ptrd, ptrdm )
199	! and out put the starting mean value and the total trends
200	! (i.e. difference between starting and ending values)
201	! hat : CALL trd_hat ( kt, ptrd, ptrdm )
202	! and output the starting hat value and the total hat trends
203	!
204	! swaps : hmlb <== hmln <== zhmla
205	! tb_mlb <== tn_mln <== zta_mla
206	! tb_mln <== ztn_mla ==>> now T over after h, need to be computed here
207	! to be used at next time step (unfiltered before)
208	!
209	END SELECT
210	!
211	END SUBROUTINE trd_tra_mxl
212
213
214	SUBROUTINE trd_mean( kt, ptrd, ptrdm )
215	!!----------------------------------------------------------------------
216	!! * ROUTINE trd_mean *
217	!!
218	!! ** Purpose :
219	!!----------------------------------------------------------------------
220	REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: ptrd ! trend at kt
221	REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: ptrdm ! cumulative trends at kt
222	INTEGER , INTENT(in ) :: kt ! time step index
223	!!----------------------------------------------------------------------
224	!
225	IF ( kt == nn_it000 ) ptrdm(:,:,:) = 0._wp
226	!
227	ptrdm(:,:,:) = ptrdm(:,:,:) + ptrd(:,:,:)
228	!
229	IF ( MOD( kt - nn_it000 + 1, nn_trd ) == 0 ) THEN
230	!
231	! call iom put???? avec en argument le tableau de nom des trends?
232	!
233	ENDIF
234	!
235	END SUBROUTINE trd_mean
236
237
238	SUBROUTINE trd_mxl_zint( pttrdmxl, pstrdmxl, ktrd, ctype )
239	!!----------------------------------------------------------------------
240	!! * ROUTINE trd_mxl_zint *
241	!!
242	!! ** Purpose : Compute the vertical average of the 3D fields given as arguments
243	!! to the subroutine. This vertical average is performed from ocean
244	!! surface down to a chosen control surface.
245	!!
246	!! ** Method/usage :
247	!! The control surface can be either a mixed layer depth (time varying)
248	!! or a fixed surface (jk level or bowl).
249	!! Choose control surface with nn_ctls in namelist NAMTRD :
250	!! nn_ctls = 0 : use mixed layer with density criterion
251	!! nn_ctls = 1 : read index from file 'ctlsurf_idx'
252	!! nn_ctls > 1 : use fixed level surface jk = nn_ctls
253	!! Note: in the remainder of the routine, the volume between the
254	!! surface and the control surface is called "mixed-layer"
255	!!----------------------------------------------------------------------
256	INTEGER , INTENT( in ) :: ktrd ! ocean trend index
257	CHARACTER(len=2) , INTENT( in ) :: ctype ! 2D surface/bottom or 3D interior physics
258	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: pttrdmxl ! temperature trend
259	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: pstrdmxl ! salinity trend
260	!
261	INTEGER :: ji, jj, jk, isum
262	REAL(wp), POINTER, DIMENSION(:,:) :: zvlmsk
263	!!----------------------------------------------------------------------
264
265	CALL wrk_alloc( jpi, jpj, zvlmsk )
266
267	! I. Definition of control surface and associated fields
268	! ------------------------------------------------------
269	! ==> only once per time step <==
270
271	IF( icount == 1 ) THEN
272	!
273
274	!!gm BUG?
275	!!gm CAUTION: double check the definition of nmln it is the nb of w-level, not t-level I guess
276
277
278	! ... Set nmxl(ji,jj) = index of first T point below control surf. or outside mixed-layer
279	IF( nn_ctls == 0 ) THEN ! * control surface = mixed-layer with density criterion
280	nmxl(:,:) = nmln(:,:) ! array nmln computed in zdfmxl.F90
281	ELSEIF( nn_ctls == 1 ) THEN ! * control surface = read index from file
282	nmxl(:,:) = nbol(:,:)
283	ELSEIF( nn_ctls >= 2 ) THEN ! * control surface = model level
284	nn_ctls = MIN( nn_ctls, jpktrd - 1 )
285	nmxl(:,:) = nn_ctls + 1
286	ENDIF
287
288	END IF
289	!
290	CALL wrk_dealloc( jpi, jpj, zvlmsk )
291	!
292	END SUBROUTINE trd_mxl_zint
293
294
295	SUBROUTINE trd_mxl( kt, p2dt )
296	!!----------------------------------------------------------------------
297	!! * ROUTINE trd_mxl *
298	!!
299	!! ** Purpose : Compute and cumulate the mixed layer trends over an analysis
300	!! period, and write NetCDF (or dimg) outputs.
301	!!
302	!! ** Method/usage :
303	!! The stored trends can be chosen twofold (according to the ln_trdmxl_instant
304	!! logical namelist variable) :
305	!! 1) to explain the difference between initial and final
306	!! mixed-layer T & S (where initial and final relate to the
307	!! current analysis window, defined by nn_trd in the namelist)
308	!! 2) to explain the difference between the current and previous
309	!! TIME-AVERAGED mixed-layer T & S (where time-averaging is
310	!! performed over each analysis window).
311	!!
312	!! ** Consistency check :
313	!! If the control surface is fixed ( nn_ctls > 1 ), the residual term (dh/dt
314	!! entrainment) should be zero, at machine accuracy. Note that in the case
315	!! of time-averaged mixed-layer fields, this residual WILL NOT BE ZERO
316	!! over the first two analysis windows (except if restart).
317	!! N.B. For ORCA2_LIM, use e.g. nn_trd=5, rn_ucf=1., nn_ctls=8
318	!! for checking residuals.
319	!! On a NEC-SX5 computer, this typically leads to:
320	!! O(1.e-20) temp. residuals (tml_res) when ln_trdmxl_instant=.false.
321	!! O(1.e-21) temp. residuals (tml_res) when ln_trdmxl_instant=.true.
322	!!
323	!! ** Action :
324	!! At each time step, mixed-layer averaged trends are stored in the
325	!! tmltrd(:,:,jpmxl_xxx) array (see trdmxl_oce.F90 for definitions of jpmxl_xxx).
326	!! This array is known when trd_mxl is called, at the end of the stp subroutine,
327	!! except for the purely vertical K_z diffusion term, which is embedded in the
328	!! lateral diffusion trend.
329	!!
330	!! In I), this K_z term is diagnosed and stored, thus its contribution is removed
331	!! from the lateral diffusion trend.
332	!! In II), the instantaneous mixed-layer T & S are computed, and misc. cumulative
333	!! arrays are updated.
334	!! In III), called only once per analysis window, we compute the total trends,
335	!! along with the residuals and the Asselin correction terms.
336	!! In IV), the appropriate trends are written in the trends NetCDF file.
337	!!
338	!! References : Vialard et al.,2001, JPO.
339	!!----------------------------------------------------------------------
340	INTEGER , INTENT(in ) :: kt ! ocean time-step index
341	REAL(wp), INTENT(in ) :: p2dt ! time step [s]
342	!
343	INTEGER :: ji, jj, jk, jl, ik, it, itmod
344	LOGICAL :: lldebug = .TRUE.
345	REAL(wp) :: zavt, zfn, zfn2
346	! ! z(ts)mltot : dT/dt over the anlysis window (including Asselin)
347	! ! z(ts)mlres : residual = dh/dt entrainment term
348	REAL(wp), POINTER, DIMENSION(:,: ) :: ztmltot , zsmltot , ztmlres , zsmlres , ztmlatf , zsmlatf
349	REAL(wp), POINTER, DIMENSION(:,: ) :: ztmltot2, zsmltot2, ztmlres2, zsmlres2, ztmlatf2, zsmlatf2, ztmltrdm2, zsmltrdm2
350	REAL(wp), POINTER, DIMENSION(:,:,:) :: ztmltrd2, zsmltrd2 ! only needed for mean diagnostics
351	#if defined key_dimgout
352	INTEGER :: iyear,imon,iday
353	CHARACTER(LEN=80) :: cltext, clmode
354	#endif
355	!!----------------------------------------------------------------------
356
357	CALL wrk_alloc( jpi, jpj, ztmltot , zsmltot , ztmlres , zsmlres , ztmlatf , zsmlatf )
358	CALL wrk_alloc( jpi, jpj, ztmltot2, zsmltot2, ztmlres2, zsmlres2, ztmlatf2, zsmlatf2, ztmltrdm2, zsmltrdm2 )
359	CALL wrk_alloc( jpi, jpj, jpltrd, ztmltrd2, zsmltrd2 )
360
361	! ======================================================================
362	! II. Cumulate the trends over the analysis window
363	! ======================================================================
364
365	ztmltrd2(:,:,:) = 0.e0 ; zsmltrd2(:,:,:) = 0.e0 ! <<< reset arrays to zero
366	ztmltot2(:,:) = 0.e0 ; zsmltot2(:,:) = 0.e0
367	ztmlres2(:,:) = 0.e0 ; zsmlres2(:,:) = 0.e0
368	ztmlatf2(:,:) = 0.e0 ; zsmlatf2(:,:) = 0.e0
369
370	! II.1 Set before values of vertically average T and S
371	! ----------------------------------------------------
372	IF( kt > nit000 ) THEN
373	! ... temperature ... ... salinity ...
374	tmlb (:,:) = tml (:,:) ; smlb (:,:) = sml (:,:)
375	tmlatfn(:,:) = tmltrd(:,:,jpmxl_atf) ; smlatfn(:,:) = smltrd(:,:,jpmxl_atf)
376	END IF
377
378
379	! II.3 Initialize mixed-layer "before" arrays for the 1rst analysis window
380	! ------------------------------------------------------------------------
381	IF( kt == 2 ) THEN ! i.e. ( .NOT. ln_rstart ).AND.( kt == nit000 + 1)
382	!
383	! ... temperature ... ... salinity ...
384	tmlbb (:,:) = tmlb (:,:) ; smlbb (:,:) = smlb (:,:)
385	tmlbn (:,:) = tml (:,:) ; smlbn (:,:) = sml (:,:)
386	tmlatfb(:,:) = tmlatfn(:,:) ; smlatfb(:,:) = smlatfn(:,:)
387
388	tmltrd_csum_ub (:,:,:) = 0.e0 ; smltrd_csum_ub (:,:,:) = 0.e0
389	tmltrd_atf_sumb(:,:) = 0.e0 ; smltrd_atf_sumb(:,:) = 0.e0
390
391	hmxlbn(:,:) = hmxl(:,:)
392
393	IF( ln_ctl ) THEN
394	WRITE(numout,*) ' we reach kt == nit000 + 1 = ', nit000+1
395	IF(lflush) CALL flush(numout)
396	CALL prt_ctl(tab2d_1=tmlbb , clinfo1=' tmlbb - : ', mask1=tmask, ovlap=1)
397	CALL prt_ctl(tab2d_1=tmlbn , clinfo1=' tmlbn - : ', mask1=tmask, ovlap=1)
398	CALL prt_ctl(tab2d_1=tmlatfb , clinfo1=' tmlatfb - : ', mask1=tmask, ovlap=1)
399	END IF
400	!
401	END IF
402
403	IF( ( ln_rstart ) .AND. ( kt == nit000 ) .AND. ( ln_ctl ) ) THEN
404	IF( ln_trdmxl_instant ) THEN
405	WRITE(numout,*) ' restart from kt == nit000 = ', nit000
406	IF(lflush) CALL flush(numout)
407	CALL prt_ctl(tab2d_1=tmlbb , clinfo1=' tmlbb - : ', mask1=tmask, ovlap=1)
408	CALL prt_ctl(tab2d_1=tmlbn , clinfo1=' tmlbn - : ', mask1=tmask, ovlap=1)
409	CALL prt_ctl(tab2d_1=tmlatfb , clinfo1=' tmlatfb - : ', mask1=tmask, ovlap=1)
410	ELSE
411	WRITE(numout,*) ' restart from kt == nit000 = ', nit000
412	IF(lflush) CALL flush(numout)
413	CALL prt_ctl(tab2d_1=tmlbn , clinfo1=' tmlbn - : ', mask1=tmask, ovlap=1)
414	CALL prt_ctl(tab2d_1=hmxlbn , clinfo1=' hmxlbn - : ', mask1=tmask, ovlap=1)
415	CALL prt_ctl(tab2d_1=tml_sumb , clinfo1=' tml_sumb - : ', mask1=tmask, ovlap=1)
416	CALL prt_ctl(tab2d_1=tmltrd_atf_sumb, clinfo1=' tmltrd_atf_sumb - : ', mask1=tmask, ovlap=1)
417	CALL prt_ctl(tab3d_1=tmltrd_csum_ub , clinfo1=' tmltrd_csum_ub - : ', mask1=tmask, ovlap=1, kdim=1)
418	END IF
419	END IF
420
421	! II.4 Cumulated trends over the analysis period
422	! ----------------------------------------------
423	!
424	! [ 1rst analysis window ] [ 2nd analysis window ]
425	!
426	! o---[--o-----o-----o-----o--]-[--o-----o-----o-----o-----o--]---o-----o--> time steps
427	! nn_trd 2*nn_trd etc.
428	! 1 2 3 4 =5 e.g. =10
429	!
430	IF( ( kt >= 2 ).OR.( ln_rstart ) ) THEN
431	!
432	nmoymltrd = nmoymltrd + 1
433
434	! ... Cumulate over BOTH physical contributions AND over time steps
435	DO jl = 1, jpltrd
436	tmltrdm(:,:) = tmltrdm(:,:) + tmltrd(:,:,jl)
437	smltrdm(:,:) = smltrdm(:,:) + smltrd(:,:,jl)
438	END DO
439
440	! ... Special handling of the Asselin trend
441	tmlatfm(:,:) = tmlatfm(:,:) + tmlatfn(:,:)
442	smlatfm(:,:) = smlatfm(:,:) + smlatfn(:,:)
443
444	! ... Trends associated with the time mean of the ML T/S
445	tmltrd_sum (:,:,:) = tmltrd_sum (:,:,:) + tmltrd (:,:,:) ! tem
446	tmltrd_csum_ln(:,:,:) = tmltrd_csum_ln(:,:,:) + tmltrd_sum(:,:,:)
447	tml_sum (:,:) = tml_sum (:,:) + tml (:,:)
448	smltrd_sum (:,:,:) = smltrd_sum (:,:,:) + smltrd (:,:,:) ! sal
449	smltrd_csum_ln(:,:,:) = smltrd_csum_ln(:,:,:) + smltrd_sum(:,:,:)
450	sml_sum (:,:) = sml_sum (:,:) + sml (:,:)
451	hmxl_sum (:,:) = hmxl_sum (:,:) + hmxl (:,:) ! rmxl
452	!
453	END IF
454
455	! ======================================================================
456	! III. Prepare fields for output (get here ONCE PER ANALYSIS PERIOD)
457	! ======================================================================
458
459	! Convert to appropriate physical units
460	! N.B. It may be useful to check IOIPSL time averaging with :
461	! tmltrd (:,:,:) = 1. ; smltrd (:,:,:) = 1.
462	tmltrd(:,:,:) = tmltrd(:,:,:) * rn_ucf ! (actually needed for 1:jpltrd-1, but trdmxl(:,:,jpltrd)
463	smltrd(:,:,:) = smltrd(:,:,:) * rn_ucf ! is no longer used, and is reset to 0. at next time step)
464
465	! define time axis
466	it = kt
467	itmod = kt - nit000 + 1
468
469	MODULO_NTRD : IF( MOD( itmod, nn_trd ) == 0 ) THEN ! nitend MUST be multiple of nn_trd
470	!
471	ztmltot (:,:) = 0.e0 ; zsmltot (:,:) = 0.e0 ! reset arrays to zero
472	ztmlres (:,:) = 0.e0 ; zsmlres (:,:) = 0.e0
473	ztmltot2(:,:) = 0.e0 ; zsmltot2(:,:) = 0.e0
474	ztmlres2(:,:) = 0.e0 ; zsmlres2(:,:) = 0.e0
475
476	zfn = REAL( nmoymltrd, wp ) ; zfn2 = zfn * zfn
477
478	! III.1 Prepare fields for output ("instantaneous" diagnostics)
479	! -------------------------------------------------------------
480
481	!-- Compute total trends
482	ztmltot(:,:) = ( tml(:,:) - tmlbn(:,:) + tmlb(:,:) - tmlbb(:,:) ) / p2dt
483	zsmltot(:,:) = ( sml(:,:) - smlbn(:,:) + smlb(:,:) - smlbb(:,:) ) / p2dt
484
485	!-- Compute residuals
486	ztmlres(:,:) = ztmltot(:,:) - ( tmltrdm(:,:) - tmlatfn(:,:) + tmlatfb(:,:) )
487	zsmlres(:,:) = zsmltot(:,:) - ( smltrdm(:,:) - smlatfn(:,:) + smlatfb(:,:) )
488
489	!-- Diagnose Asselin trend over the analysis window
490	ztmlatf(:,:) = tmlatfm(:,:) - tmlatfn(:,:) + tmlatfb(:,:)
491	zsmlatf(:,:) = smlatfm(:,:) - smlatfn(:,:) + smlatfb(:,:)
492
493	!-- Lateral boundary conditions
494	! ... temperature ... ... salinity ...
495	CALL lbc_lnk( ztmltot , 'T', 1. ) ; CALL lbc_lnk( zsmltot , 'T', 1. )
496	CALL lbc_lnk( ztmlres , 'T', 1. ) ; CALL lbc_lnk( zsmlres , 'T', 1. )
497	CALL lbc_lnk( ztmlatf , 'T', 1. ) ; CALL lbc_lnk( zsmlatf , 'T', 1. )
498
499
500	! III.2 Prepare fields for output ("mean" diagnostics)
501	! ----------------------------------------------------
502
503	!-- Update the ML depth time sum (to build the Leap-Frog time mean)
504	hmxl_sum(:,:) = hmxlbn(:,:) + 2 * ( hmxl_sum(:,:) - hmxl(:,:) ) + hmxl(:,:)
505
506	!-- Compute temperature total trends
507	tml_sum (:,:) = tmlbn(:,:) + 2 * ( tml_sum(:,:) - tml(:,:) ) + tml(:,:)
508	ztmltot2(:,:) = ( tml_sum(:,:) - tml_sumb(:,:) ) / p2dt ! now in degC/s
509
510	!-- Compute salinity total trends
511	sml_sum (:,:) = smlbn(:,:) + 2 * ( sml_sum(:,:) - sml(:,:) ) + sml(:,:)
512	zsmltot2(:,:) = ( sml_sum(:,:) - sml_sumb(:,:) ) / p2dt ! now in psu/s
513
514	!-- Compute temperature residuals
515	DO jl = 1, jpltrd
516	ztmltrd2(:,:,jl) = tmltrd_csum_ub(:,:,jl) + tmltrd_csum_ln(:,:,jl)
517	END DO
518
519	ztmltrdm2(:,:) = 0.e0
520	DO jl = 1, jpltrd
521	ztmltrdm2(:,:) = ztmltrdm2(:,:) + ztmltrd2(:,:,jl)
522	END DO
523
524	ztmlres2(:,:) = ztmltot2(:,:) - &
525	( ztmltrdm2(:,:) - tmltrd_sum(:,:,jpmxl_atf) + tmltrd_atf_sumb(:,:) )
526
527	!-- Compute salinity residuals
528	DO jl = 1, jpltrd
529	zsmltrd2(:,:,jl) = smltrd_csum_ub(:,:,jl) + smltrd_csum_ln(:,:,jl)
530	END DO
531
532	zsmltrdm2(:,:) = 0.
533	DO jl = 1, jpltrd
534	zsmltrdm2(:,:) = zsmltrdm2(:,:) + zsmltrd2(:,:,jl)
535	END DO
536
537	zsmlres2(:,:) = zsmltot2(:,:) - &
538	( zsmltrdm2(:,:) - smltrd_sum(:,:,jpmxl_atf) + smltrd_atf_sumb(:,:) )
539
540	!-- Diagnose Asselin trend over the analysis window
541	ztmlatf2(:,:) = ztmltrd2(:,:,jpmxl_atf) - tmltrd_sum(:,:,jpmxl_atf) + tmltrd_atf_sumb(:,:)
542	zsmlatf2(:,:) = zsmltrd2(:,:,jpmxl_atf) - smltrd_sum(:,:,jpmxl_atf) + smltrd_atf_sumb(:,:)
543
544	!-- Lateral boundary conditions
545	! ... temperature ... ... salinity ...
546	CALL lbc_lnk( ztmltot2, 'T', 1. ) ; CALL lbc_lnk( zsmltot2, 'T', 1. )
547	CALL lbc_lnk( ztmlres2, 'T', 1. ) ; CALL lbc_lnk( zsmlres2, 'T', 1. )
548	DO jl = 1, jpltrd
549	CALL lbc_lnk( ztmltrd2(:,:,jl), 'T', 1. ) ! \ these will be output
550	CALL lbc_lnk( zsmltrd2(:,:,jl), 'T', 1. ) ! / in the NetCDF trends file
551	END DO
552
553	! III.3 Time evolution array swap
554	! -------------------------------
555
556	! For T/S instantaneous diagnostics
557	! ... temperature ... ... salinity ...
558	tmlbb (:,:) = tmlb (:,:) ; smlbb (:,:) = smlb (:,:)
559	tmlbn (:,:) = tml (:,:) ; smlbn (:,:) = sml (:,:)
560	tmlatfb(:,:) = tmlatfn(:,:) ; smlatfb(:,:) = smlatfn(:,:)
561
562	! For T mean diagnostics
563	tmltrd_csum_ub (:,:,:) = zfn * tmltrd_sum(:,:,:) - tmltrd_csum_ln(:,:,:)
564	tml_sumb (:,:) = tml_sum(:,:)
565	tmltrd_atf_sumb(:,:) = tmltrd_sum(:,:,jpmxl_atf)
566
567	! For S mean diagnostics
568	smltrd_csum_ub (:,:,:) = zfn * smltrd_sum(:,:,:) - smltrd_csum_ln(:,:,:)
569	sml_sumb (:,:) = sml_sum(:,:)
570	smltrd_atf_sumb(:,:) = smltrd_sum(:,:,jpmxl_atf)
571
572	! ML depth
573	hmxlbn (:,:) = hmxl (:,:)
574
575	IF( ln_ctl ) THEN
576	IF( ln_trdmxl_instant ) THEN
577	CALL prt_ctl(tab2d_1=tmlbb , clinfo1=' tmlbb - : ', mask1=tmask, ovlap=1)
578	CALL prt_ctl(tab2d_1=tmlbn , clinfo1=' tmlbn - : ', mask1=tmask, ovlap=1)
579	CALL prt_ctl(tab2d_1=tmlatfb , clinfo1=' tmlatfb - : ', mask1=tmask, ovlap=1)
580	ELSE
581	CALL prt_ctl(tab2d_1=tmlbn , clinfo1=' tmlbn - : ', mask1=tmask, ovlap=1)
582	CALL prt_ctl(tab2d_1=hmxlbn , clinfo1=' hmxlbn - : ', mask1=tmask, ovlap=1)
583	CALL prt_ctl(tab2d_1=tml_sumb , clinfo1=' tml_sumb - : ', mask1=tmask, ovlap=1)
584	CALL prt_ctl(tab2d_1=tmltrd_atf_sumb, clinfo1=' tmltrd_atf_sumb - : ', mask1=tmask, ovlap=1)
585	CALL prt_ctl(tab3d_1=tmltrd_csum_ub , clinfo1=' tmltrd_csum_ub - : ', mask1=tmask, ovlap=1, kdim=1)
586	END IF
587	END IF
588
589	! III.4 Convert to appropriate physical units
590	! -------------------------------------------
591
592	! ... temperature ... ... salinity ...
593	ztmltot (:,:) = ztmltot(:,:) * rn_ucf/zfn ; zsmltot (:,:) = zsmltot(:,:) * rn_ucf/zfn
594	ztmlres (:,:) = ztmlres(:,:) * rn_ucf/zfn ; zsmlres (:,:) = zsmlres(:,:) * rn_ucf/zfn
595	ztmlatf (:,:) = ztmlatf(:,:) * rn_ucf/zfn ; zsmlatf (:,:) = zsmlatf(:,:) * rn_ucf/zfn
596
597	tml_sum (:,:) = tml_sum (:,:) / (2zfn) ; sml_sum (:,:) = sml_sum (:,:) / (2zfn)
598	ztmltot2(:,:) = ztmltot2(:,:) * rn_ucf/zfn2 ; zsmltot2(:,:) = zsmltot2(:,:) * rn_ucf/zfn2
599	ztmltrd2(:,:,:) = ztmltrd2(:,:,:)* rn_ucf/zfn2 ; zsmltrd2(:,:,:) = zsmltrd2(:,:,:)* rn_ucf/zfn2
600	ztmlatf2(:,:) = ztmlatf2(:,:) * rn_ucf/zfn2 ; zsmlatf2(:,:) = zsmlatf2(:,:) * rn_ucf/zfn2
601	ztmlres2(:,:) = ztmlres2(:,:) * rn_ucf/zfn2 ; zsmlres2(:,:) = zsmlres2(:,:) * rn_ucf/zfn2
602
603	hmxl_sum(:,:) = hmxl_sum(:,:) / (2*zfn) ! similar to tml_sum and sml_sum
604
605	! * Debugging information *
606	IF( lldebug ) THEN
607	!
608	WRITE(numout,*)
609	WRITE(numout,*) 'trd_mxl : write trends in the Mixed Layer for debugging process:'
610	WRITE(numout,*) '~~~~~~~ '
611	WRITE(numout,*) ' TRA kt = ', kt, 'nmoymltrd = ', nmoymltrd
612	WRITE(numout,*)
613	WRITE(numout,*) ' >>>>>>>>>>>>>>>>>> TRA TEMPERATURE <<<<<<<<<<<<<<<<<<'
614	WRITE(numout,*) ' TRA ztmlres : ', SUM(ztmlres(:,:))
615	WRITE(numout,*) ' TRA ztmltot : ', SUM(ztmltot(:,:))
616	WRITE(numout,*) ' TRA tmltrdm : ', SUM(tmltrdm(:,:))
617	WRITE(numout,*) ' TRA tmlatfb : ', SUM(tmlatfb(:,:))
618	WRITE(numout,*) ' TRA tmlatfn : ', SUM(tmlatfn(:,:))
619	DO jl = 1, jpltrd
620	WRITE(numout,) ' TRA TREND INDEX jpmxl_xxx = jl = ', jl, &
621	& ' tmltrd : ', SUM(tmltrd(:,:,jl))
622	END DO
623	WRITE(numout,*) ' TRA ztmlres (jpi/2,jpj/2) : ', ztmlres (jpi/2,jpj/2)
624	WRITE(numout,*) ' TRA ztmlres2(jpi/2,jpj/2) : ', ztmlres2(jpi/2,jpj/2)
625	WRITE(numout,*)
626	WRITE(numout,*) ' >>>>>>>>>>>>>>>>>> TRA SALINITY <<<<<<<<<<<<<<<<<<'
627	WRITE(numout,*) ' TRA zsmlres : ', SUM(zsmlres(:,:))
628	WRITE(numout,*) ' TRA zsmltot : ', SUM(zsmltot(:,:))
629	WRITE(numout,*) ' TRA smltrdm : ', SUM(smltrdm(:,:))
630	WRITE(numout,*) ' TRA smlatfb : ', SUM(smlatfb(:,:))
631	WRITE(numout,*) ' TRA smlatfn : ', SUM(smlatfn(:,:))
632	DO jl = 1, jpltrd
633	WRITE(numout,) ' TRA TREND INDEX jpmxl_xxx = jl = ', jl, &
634	& ' smltrd : ', SUM(smltrd(:,:,jl))
635	END DO
636	WRITE(numout,*) ' TRA zsmlres (jpi/2,jpj/2) : ', zsmlres (jpi/2,jpj/2)
637	WRITE(numout,*) ' TRA zsmlres2(jpi/2,jpj/2) : ', zsmlres2(jpi/2,jpj/2)
638	IF(lflush) CALL flush(numout)
639	!
640	END IF
641	!
642	END IF MODULO_NTRD
643
644	! ======================================================================
645	! IV. Write trends in the NetCDF file
646	! ======================================================================
647
648	!-- Write the trends for T/S instantaneous diagnostics
649
650	IF( ln_trdmxl_instant ) THEN
651
652	CALL iom_put( "mxl_depth", hmxl(:,:) )
653
654	!................................. ( ML temperature ) ...................................
655
656	!-- Output the fields
657	CALL iom_put( "tml" , tml (:,:) )
658	CALL iom_put( "tml_tot" , ztmltot(:,:) )
659	CALL iom_put( "tml_res" , ztmlres(:,:) )
660
661	DO jl = 1, jpltrd - 1
662	CALL iom_put( trim("tml"//ctrd(jl,2)), tmltrd (:,:,jl) )
663	END DO
664
665	CALL iom_put( trim("tml"//ctrd(jpmxl_atf,2)), ztmlatf(:,:) )
666
667	!.................................. ( ML salinity ) .....................................
668
669	!-- Output the fields
670	CALL iom_put( "sml" , sml (:,:) )
671	CALL iom_put( "sml_tot", zsmltot(:,:) )
672	CALL iom_put( "sml_res", zsmlres(:,:) )
673
674	DO jl = 1, jpltrd - 1
675	CALL iom_put( trim("sml"//ctrd(jl,2)), smltrd(:,:,jl) )
676	END DO
677
678	CALL iom_put( trim("sml"//ctrd(jpmxl_atf,2)), zsmlatf(:,:) )
679
680
681
682	ELSE !-- Write the trends for T/S mean diagnostics
683
684	CALL iom_put( "mxl_depth", hmxl_sum(:,:) )
685
686	!................................. ( ML temperature ) ...................................
687
688	!-- Output the fields
689	CALL iom_put( "tml" , tml_sum (:,:) )
690	CALL iom_put( "tml_tot" , ztmltot2(:,:) )
691	CALL iom_put( "tml_res" , ztmlres2(:,:) )
692
693	DO jl = 1, jpltrd - 1
694	CALL iom_put( trim("tml"//ctrd(jl,2)), ztmltrd2(:,:,jl) )
695	END DO
696
697	CALL iom_put( trim("tml"//ctrd(jpmxl_atf,2)), ztmlatf2(:,:) )
698
699	!.................................. ( ML salinity ) .....................................
700
701	!-- Output the fields
702	CALL iom_put( "sml" , sml_sum (:,:) )
703	CALL iom_put( "sml_tot", zsmltot2(:,:) )
704	CALL iom_put( "sml_res", zsmlres2(:,:) )
705
706	DO jl = 1, jpltrd - 1
707	CALL iom_put( trim("sml"//ctrd(jl,2)), zsmltrd2(:,:,jl) )
708	END DO
709
710	CALL iom_put( trim("sml"//ctrd(jpmxl_atf,2)), zsmlatf2(:,:) )
711	!
712	END IF
713	!
714
715	IF( MOD( itmod, nn_trd ) == 0 ) THEN
716	!
717	! III.5 Reset cumulative arrays to zero
718	! -------------------------------------
719	nmoymltrd = 0
720
721	! ... temperature ... ... salinity ...
722	tmltrdm (:,:) = 0.e0 ; smltrdm (:,:) = 0.e0
723	tmlatfm (:,:) = 0.e0 ; smlatfm (:,:) = 0.e0
724	tml_sum (:,:) = 0.e0 ; sml_sum (:,:) = 0.e0
725	tmltrd_csum_ln (:,:,:) = 0.e0 ; smltrd_csum_ln (:,:,:) = 0.e0
726	tmltrd_sum (:,:,:) = 0.e0 ; smltrd_sum (:,:,:) = 0.e0
727
728	hmxl_sum (:,:) = 0.e0
729	!
730	END IF
731
732	! ======================================================================
733	! V. Write restart file
734	! ======================================================================
735
736	IF( lrst_oce ) CALL trd_mxl_rst_write( kt )
737
738	CALL wrk_dealloc( jpi, jpj, ztmltot , zsmltot , ztmlres , zsmlres , ztmlatf , zsmlatf )
739	CALL wrk_dealloc( jpi, jpj, ztmltot2, zsmltot2, ztmlres2, zsmlres2, ztmlatf2, zsmlatf2, ztmltrdm2, zsmltrdm2 )
740	CALL wrk_dealloc( jpi, jpj, jpltrd, ztmltrd2, zsmltrd2 )
741	!
742	END SUBROUTINE trd_mxl
743
744
745	SUBROUTINE trd_mxl_init
746	!!----------------------------------------------------------------------
747	!! * ROUTINE trd_mxl_init *
748	!!
749	!! ** Purpose : computation of vertically integrated T and S budgets
750	!! from ocean surface down to control surface (NetCDF output)
751	!!----------------------------------------------------------------------
752	INTEGER :: jl ! dummy loop indices
753	INTEGER :: inum ! logical unit
754	INTEGER :: ios ! local integer
755	REAL(wp) :: zjulian, zsto, zout
756	CHARACTER (LEN=40) :: clop
757	CHARACTER (LEN=12) :: clmxl, cltu, clsu
758	!!
759	NAMELIST/namtrd_mxl/ nn_trd , cn_trdrst_in , ln_trdmxl_restart, &
760	& nn_ctls, cn_trdrst_out, ln_trdmxl_instant, rn_ucf, rn_rho_c
761	!!----------------------------------------------------------------------
762	!
763	REWIND( numnam_ref ) ! Namelist namtrd_mxl in reference namelist : mixed layer trends diagnostic
764	READ ( numnam_ref, namtrd_mxl, IOSTAT = ios, ERR = 901 )
765	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrd_mxl in reference namelist', lwp )
766
767	REWIND( numnam_cfg ) ! Namelist namtrd_mxl in configuration namelist : mixed layer trends diagnostic
768	READ ( numnam_cfg, namtrd_mxl, IOSTAT = ios, ERR = 902 )
769	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrd_mxl in configuration namelist', lwp )
770	IF(lwm .AND. nprint > 2) WRITE( numond, namtrd_mxl )
771	!
772	IF(lwp) THEN ! control print
773	WRITE(numout,*)
774	WRITE(numout,*) ' trd_mxl_init : Mixed-layer trends'
775	WRITE(numout,*) ' ~~~~~~~~~~'
776	WRITE(numout,*) ' Namelist namtrd : set trends parameters'
777	WRITE(numout,*) ' frequency of trends diagnostics (glo) nn_trd = ', nn_trd
778	WRITE(numout,*) ' density criteria used to defined the MLD rn_rho_c = ', rn_rho_c
779	WRITE(numout,*) ' control surface type (mld) nn_ctls = ', nn_ctls
780	WRITE(numout,*) ' restart for ML diagnostics ln_trdmxl_restart = ', ln_trdmxl_restart
781	WRITE(numout,*) ' instantaneous or mean ML T/S ln_trdmxl_instant = ', ln_trdmxl_instant
782	WRITE(numout,*) ' unit conversion factor rn_ucf = ', rn_ucf
783	WRITE(numout,*) ' criteria to compute the MLD rn_rho_c = ', rn_rho_c
784	IF(lflush) CALL flush(numout)
785	ENDIF
786
787
788
789	! I.1 Check consistency of user defined preferences
790	! -------------------------------------------------
791
792	IF ( rn_rho_c /= rho_c ) CALL ctl_warn( 'Unless you have good reason to do so, you should use the value ', &
793	& 'defined in zdfmxl.F90 module to calculate the mixed layer depth' )
794
795	IF( MOD( nitend, nn_trd ) /= 0 ) THEN
796	WRITE(numout,cform_err)
797	WRITE(numout,*) ' Your nitend parameter, nitend = ', nitend
798	WRITE(numout,*) ' is no multiple of the trends diagnostics frequency '
799	WRITE(numout,*) ' you defined, nn_trd = ', nn_trd
800	WRITE(numout,*) ' This will not allow you to restart from this simulation. '
801	WRITE(numout,*) ' You should reconsider this choice. '
802	WRITE(numout,*)
803	WRITE(numout,*) ' N.B. the nitend parameter is also constrained to be a '
804	WRITE(numout,*) ' multiple of the nn_fsbc parameter '
805	nstop = nstop + 1
806	END IF
807
808	IF( nn_cla == 1 ) CALL ctl_warn( ' You set n_cla = 1. Note that the Mixed-Layer diagnostics ', &
809	& ' are not exact along the corresponding straits. ')
810
811	! ! allocate trdmxl arrays
812	IF( trd_mxl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'trd_mxl_init : unable to allocate trdmxl arrays' )
813	IF( trdmxl_oce_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'trd_mxl_init : unable to allocate trdmxl_oce arrays' )
814
815
816
817
818	nkstp = nit000 - 1 ! current time step indicator initialization
819
820
821
822
823	! I.2 Initialize arrays to zero or read a restart file
824	! ----------------------------------------------------
825
826	nmoymltrd = 0
827
828	! ... temperature ... ... salinity ...
829	tml (:,:) = 0.e0 ; sml (:,:) = 0.e0 ! inst.
830	tmltrdm (:,:) = 0.e0 ; smltrdm (:,:) = 0.e0
831	tmlatfm (:,:) = 0.e0 ; smlatfm (:,:) = 0.e0
832	tml_sum (:,:) = 0.e0 ; sml_sum (:,:) = 0.e0 ! mean
833	tmltrd_sum (:,:,:) = 0.e0 ; smltrd_sum (:,:,:) = 0.e0
834	tmltrd_csum_ln (:,:,:) = 0.e0 ; smltrd_csum_ln (:,:,:) = 0.e0
835
836	hmxl (:,:) = 0.e0
837	hmxl_sum (:,:) = 0.e0
838
839	IF( ln_rstart .AND. ln_trdmxl_restart ) THEN
840	CALL trd_mxl_rst_read
841	ELSE
842	! ... temperature ... ... salinity ...
843	tmlb (:,:) = 0.e0 ; smlb (:,:) = 0.e0 ! inst.
844	tmlbb (:,:) = 0.e0 ; smlbb (:,:) = 0.e0
845	tmlbn (:,:) = 0.e0 ; smlbn (:,:) = 0.e0
846	tml_sumb (:,:) = 0.e0 ; sml_sumb (:,:) = 0.e0 ! mean
847	tmltrd_csum_ub (:,:,:) = 0.e0 ; smltrd_csum_ub (:,:,:) = 0.e0
848	tmltrd_atf_sumb(:,:) = 0.e0 ; smltrd_atf_sumb(:,:) = 0.e0
849	END IF
850
851	icount = 1 ; ionce = 1 ! open specifier
852
853	! I.3 Read control surface from file ctlsurf_idx
854	! ----------------------------------------------
855
856	IF( nn_ctls == 1 ) THEN
857	CALL ctl_opn( inum, 'ctlsurf_idx', 'OLD', 'UNFORMATTED', 'SEQUENTIAL', -1, numout, lwp )
858	READ ( inum, * ) nbol
859	CLOSE( inum )
860	END IF
861
862	! ======================================================================
863	! II. netCDF output initialization
864	! ======================================================================
865
866	! clmxl = legend root for netCDF output
867	IF( nn_ctls == 0 ) THEN ! control surface = mixed-layer with density criterion
868	clmxl = 'Mixed Layer ' ! (array nmln computed in zdfmxl.F90)
869	ELSE IF( nn_ctls == 1 ) THEN ! control surface = read index from file
870	clmxl = ' Bowl '
871	ELSE IF( nn_ctls >= 2 ) THEN ! control surface = model level
872	WRITE(clmxl,'(A10,I2,1X)') 'Levels 1 -', nn_ctls
873	END IF
874
875
876
877	! II.3 Define the T grid trend file (nidtrd)
878	! ------------------------------------------
879	!-- Define long and short names for the NetCDF output variables
880	! ==> choose them according to trdmxl_oce.F90 <==
881
882	ctrd(jpmxl_xad,1) = " Zonal advection" ; ctrd(jpmxl_xad,2) = "_xad"
883	ctrd(jpmxl_yad,1) = " Meridional advection" ; ctrd(jpmxl_yad,2) = "_yad"
884	ctrd(jpmxl_zad,1) = " Vertical advection" ; ctrd(jpmxl_zad,2) = "_zad"
885	ctrd(jpmxl_ldf,1) = " Lateral diffusion" ; ctrd(jpmxl_ldf,2) = "_ldf"
886	ctrd(jpmxl_for,1) = " Forcing" ; ctrd(jpmxl_for,2) = "_for"
887	ctrd(jpmxl_zdf,1) = " Vertical diff. (Kz)" ; ctrd(jpmxl_zdf,2) = "_zdf"
888	ctrd(jpmxl_bbc,1) = " Geothermal flux" ; ctrd(jpmxl_bbc,2) = "_bbc"
889	ctrd(jpmxl_bbl,1) = " Adv/diff. Bottom boundary layer" ; ctrd(jpmxl_bbl,2) = "_bbl"
890	ctrd(jpmxl_dmp,1) = " Tracer damping" ; ctrd(jpmxl_dmp,2) = "_dmp"
891	ctrd(jpmxl_npc,1) = " Non penetrative convec. adjust." ; ctrd(jpmxl_npc,2) = "_npc"
892	ctrd(jpmxl_atf,1) = " Asselin time filter" ; ctrd(jpmxl_atf,2) = "_atf"
893
894
895	!-- Define physical units
896	IF ( rn_ucf == 1. ) THEN ; cltu = "degC/s" ; clsu = "p.s.u./s"
897	ELSEIF ( rn_ucf == 3600.*24.) THEN ; cltu = "degC/day" ; clsu = "p.s.u./day"
898	ELSE ; cltu = "unknown?" ; clsu = "unknown?"
899	END IF
900	!
901	END SUBROUTINE trd_mxl_init
902
903	!!======================================================================
904	END MODULE trdmxl

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

New URL for NEMO forge! http://forge.nemo-ocean.eu

Context Navigation

source: branches/UKMO/dev_r5518_GO6_package_text_diagnostics/NEMOGCM/NEMO/OPA_SRC/TRD/trdmxl.F90 @ 10774

Download in other formats: