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