/[MITgcm]/MITgcm/pkg/gmredi/gmredi_diagnostics_init.F
ViewVC logotype

Annotation of /MITgcm/pkg/gmredi/gmredi_diagnostics_init.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.11 - (hide annotations) (download)
Fri Jun 21 17:23:31 2013 UTC (10 years, 11 months ago) by m_bates
Branch: MAIN
CVS Tags: checkpoint64j
Changes since 1.10: +276 -1 lines
Added a new eddy diffusivity parameterisation pkg/gmredi.  More detailed description in tag-index.

1 m_bates 1.11 C $Header: /u/gcmpack/MITgcm/pkg/gmredi/gmredi_diagnostics_init.F,v 1.10 2012/03/20 23:36:59 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "GMREDI_OPTIONS.h"
5    
6     CBOP
7     C !ROUTINE: GMREDI_DIAGNOSTICS_INIT
8     C !INTERFACE:
9     SUBROUTINE GMREDI_DIAGNOSTICS_INIT( myThid )
10    
11     C !DESCRIPTION: \bv
12     C *==========================================================*
13     C | SUBROUTINE GMREDI_DIAGNOSTICS_INIT
14     C | o Routine to initialize list of all available diagnostics
15     C | for GM/Redi package
16     C *==========================================================*
17     C \ev
18     C !USES:
19     IMPLICIT NONE
20    
21     C === Global variables ===
22     #include "EEPARAMS.h"
23     c #include "SIZE.h"
24     c #include "PARAMS.h"
25     c #include "GMREDI.h"
26    
27     C !INPUT/OUTPUT PARAMETERS:
28     C === Routine arguments ===
29     C myThid :: my Thread Id number
30     INTEGER myThid
31     CEOP
32    
33     #ifdef ALLOW_DIAGNOSTICS
34     C !LOCAL VARIABLES:
35     C === Local variables ===
36     C diagNum :: diagnostics number in the (long) list of available diag.
37 jmc 1.7 C diagMate :: diag. mate number in the (long) list of available diag.
38 jmc 1.1 C diagName :: local short name (8c) of a diagnostics
39     C diagCode :: local parser field with characteristics of the diagnostics
40     C cf head of S/R DIAGNOSTICS_INIT_EARLY or DIAGNOSTICS_MAIN_INIT
41     C diagUnits :: local string (16c): physical units of a diagnostic field
42     C diagTitle :: local string (80c): description of field in diagnostic
43     INTEGER diagNum
44 jmc 1.7 INTEGER diagMate
45 jmc 1.1 CHARACTER*8 diagName
46     CHARACTER*16 diagCode
47     CHARACTER*16 diagUnits
48     CHARACTER*(80) diagTitle
49    
50     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
51    
52     c IF ( useDiagnotics ) THEN
53    
54     diagName = 'GM_VisbK'
55 jmc 1.7 diagTitle =
56 jmc 1.1 & 'Mixing coefficient from Visbeck etal parameterization'
57     diagUnits = 'm^2/s '
58     diagCode = 'SM P M1 '
59 jmc 1.7 CALL DIAGNOSTICS_ADDTOLIST( diagNum,
60     I diagName, diagCode, diagUnits, diagTitle, 0, myThid )
61 jmc 1.1
62 jmc 1.5 diagName = 'GM_hTrsL'
63     diagTitle = 'Base depth (>0) of the Transition Layer'
64     diagUnits = 'm '
65     diagCode = 'SM P M1 '
66 jmc 1.7 CALL DIAGNOSTICS_ADDTOLIST( diagNum,
67     I diagName, diagCode, diagUnits, diagTitle, 0, myThid )
68 jmc 1.5
69     diagName = 'GM_baseS'
70     diagTitle = 'Slope at the base of the Transition Layer'
71     diagUnits = '1 '
72     diagCode = 'SM P M1 '
73 jmc 1.7 CALL DIAGNOSTICS_ADDTOLIST( diagNum,
74     I diagName, diagCode, diagUnits, diagTitle, 0, myThid )
75 jmc 1.5
76     diagName = 'GM_rLamb'
77     diagTitle =
78     & 'Slope vertical gradient at Trans. Layer Base (=recip.Lambda)'
79     diagUnits = '1/m '
80     diagCode = 'SM P M1 '
81 jmc 1.7 CALL DIAGNOSTICS_ADDTOLIST( diagNum,
82     I diagName, diagCode, diagUnits, diagTitle, 0, myThid )
83 jmc 1.5
84 jmc 1.9 #ifndef GM_EXCLUDE_SUBMESO
85     diagName = 'SubMesLf'
86     diagTitle = 'Sub-Meso horiz. Length Scale (Lf)'
87     diagUnits = 'm '
88     diagCode = 'SM P M1 '
89     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
90     I diagName, diagCode, diagUnits, diagTitle, 0, myThid )
91    
92     diagName = 'SubMpsiX'
93     diagTitle =
94     & 'Sub-Meso transp.stream-funct. magnitude (Psi0): U component'
95     diagUnits = 'm^2/s '
96     diagCode = 'UU M1 '
97     diagMate = diagNum + 2
98     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
99     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
100    
101     diagName = 'SubMpsiY'
102     diagTitle =
103     & 'Sub-Meso transp.stream-funct. magnitude (Psi0): V component'
104     diagUnits = 'm^2/s '
105     diagCode = 'VV M1 '
106     diagMate = diagNum
107     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
108     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
109     #endif
110    
111 jmc 1.1 diagName = 'GM_Kux '
112     diagTitle = 'K_11 element (U.point, X.dir) of GM-Redi tensor'
113     diagUnits = 'm^2/s '
114 jmc 1.7 diagCode = 'UU P MR '
115     diagMate = diagNum + 2
116     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
117     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
118 jmc 1.1
119     diagName = 'GM_Kvy '
120     diagTitle = 'K_22 element (V.point, Y.dir) of GM-Redi tensor'
121     diagUnits = 'm^2/s '
122 jmc 1.7 diagCode = 'VV P MR '
123     diagMate = diagNum
124     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
125     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
126 jmc 1.1
127     diagName = 'GM_Kuz '
128     diagTitle = 'K_13 element (U.point, Z.dir) of GM-Redi tensor'
129     diagUnits = 'm^2/s '
130 jmc 1.7 diagCode = 'UU MR '
131     diagMate = diagNum + 2
132     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
133     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
134 jmc 1.1
135     diagName = 'GM_Kvz '
136     diagTitle = 'K_23 element (V.point, Z.dir) of GM-Redi tensor'
137     diagUnits = 'm^2/s '
138 jmc 1.7 diagCode = 'VV MR '
139     diagMate = diagNum
140     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
141     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
142 jmc 1.1
143     diagName = 'GM_Kwx '
144     diagTitle = 'K_31 element (W.point, X.dir) of GM-Redi tensor'
145     diagUnits = 'm^2/s '
146 jmc 1.7 diagCode = 'UM LR '
147     diagMate = diagNum + 2
148     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
149     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
150 jmc 1.1
151     diagName = 'GM_Kwy '
152     diagTitle = 'K_32 element (W.point, Y.dir) of GM-Redi tensor'
153     diagUnits = 'm^2/s '
154 jmc 1.7 diagCode = 'VM LR '
155     diagMate = diagNum
156     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
157     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
158 jmc 1.1
159     diagName = 'GM_Kwz '
160     diagTitle = 'K_33 element (W.point, Z.dir) of GM-Redi tensor'
161     diagUnits = 'm^2/s '
162     diagCode = 'WM P LR '
163 jmc 1.7 CALL DIAGNOSTICS_ADDTOLIST( diagNum,
164     I diagName, diagCode, diagUnits, diagTitle, 0, myThid )
165 jmc 1.1
166     diagName = 'GM_PsiX '
167 jmc 1.9 diagTitle = 'GM Bolus transport stream-function : U component'
168 jmc 1.1 diagUnits = 'm^2/s '
169 jmc 1.7 diagCode = 'UU LR '
170     diagMate = diagNum + 2
171     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
172     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
173 jmc 1.1
174     diagName = 'GM_PsiY '
175 jmc 1.9 diagTitle = 'GM Bolus transport stream-function : V component'
176 jmc 1.1 diagUnits = 'm^2/s '
177 jmc 1.7 diagCode = 'VV LR '
178     diagMate = diagNum
179     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
180     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
181 jmc 1.1
182 jmc 1.2 diagName = 'GM_KuzTz'
183 dfer 1.4 diagTitle = 'Redi Off-diagonal Temperature flux: X component'
184 jmc 1.2 diagUnits = 'degC.m^3/s '
185 jmc 1.7 diagCode = 'UU MR '
186     diagMate = diagNum + 2
187     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
188     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
189 jmc 1.2
190     diagName = 'GM_KvzTz'
191 dfer 1.4 diagTitle = 'Redi Off-diagonal Temperature flux: Y component'
192 jmc 1.2 diagUnits = 'degC.m^3/s '
193 jmc 1.7 diagCode = 'VV MR '
194     diagMate = diagNum
195     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
196     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
197 jmc 1.2
198 jmc 1.10 diagName = 'GM_KwzTz'
199     diagTitle = 'Redi main-diagonal vertical Temperature flux'
200     diagUnits = 'degC.m^3/s '
201     diagCode = 'WM LR '
202     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
203     I diagName, diagCode, diagUnits, diagTitle, 0, myThid )
204    
205 dfer 1.4 diagName = 'GM_ubT '
206     diagTitle = 'Zonal Mass-Weight Bolus Transp of Pot Temp'
207     diagUnits = 'degC.m^3/s '
208 jmc 1.7 diagCode = 'UUr MR '
209     diagMate = diagNum + 2
210     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
211     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
212 dfer 1.4
213     diagName = 'GM_vbT '
214     diagTitle = 'Meridional Mass-Weight Bolus Transp of Pot Temp'
215     diagUnits = 'degC.m^3/s '
216 jmc 1.7 diagCode = 'VVr MR '
217     diagMate = diagNum
218     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
219     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
220 dfer 1.4
221 gforget 1.6 #ifdef ALLOW_EDDYPSI
222 heimbach 1.3 diagName = 'GMEdTauX'
223     diagTitle = 'eddy-induced stress X-comp. estimated from Kwx'
224     diagUnits = 'N/m^2 '
225 jmc 1.7 diagCode = 'UM LR '
226     diagMate = diagNum + 2
227     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
228     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
229 heimbach 1.3
230     diagName = 'GMEdTauY'
231     diagTitle = 'eddy-induced stress Y-comp. estimated from Kwy'
232     diagUnits = 'N/m^2 '
233 jmc 1.7 diagCode = 'VM LR '
234     diagMate = diagNum
235     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
236     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
237 heimbach 1.3 #endif
238    
239 jmc 1.8 #ifdef GM_BOLUS_BVP
240     diagName = 'GM_BVPcW'
241     diagTitle = 'WKB wave speed (at Western edge location)'
242     diagUnits = 'm/s '
243     diagCode = 'SU P M1 '
244     diagMate = diagNum + 2
245     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
246     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
247    
248     diagName = 'GM_BVPcS'
249     diagTitle = 'WKB wave speed (at Southern edge location)'
250     diagUnits = 'm/s '
251     diagCode = 'SV P M1 '
252     diagMate = diagNum
253     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
254     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
255     #endif
256    
257 m_bates 1.11 #ifdef GM_K3D
258     diagName = 'GM_K3D'
259     diagTitle = '3D diffusivity'
260     diagUnits = 'm**2/s '
261     diagCode = 'SM P MR '
262     diagMate = diagNum
263     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
264     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
265    
266     diagName = 'GM_A3D'
267     diagTitle = '3D lower diagona'
268     diagUnits = '1/m**2 '
269     diagCode = 'SM MR '
270     diagMate = diagNum
271     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
272     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
273    
274     diagName = 'GM_B3D'
275     diagTitle = '3D lower diagona'
276     diagUnits = '1/m**2 '
277     diagCode = 'SM MR '
278     diagMate = diagNum
279     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
280     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
281    
282     diagName = 'GM_C3D'
283     diagTitle = '3D lower diagona'
284     diagUnits = '1/m**2 '
285     diagCode = 'SM MR '
286     diagMate = diagNum
287     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
288     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
289    
290     diagName = 'GM_VAL'
291     diagTitle = 'Eigen value '
292     diagUnits = '1/m**2 '
293     diagCode = 'SM P MR '
294     diagMate = diagNum
295     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
296     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
297    
298     diagName = 'GM_VEC'
299     diagTitle = 'Eigen vector for the first barcolinic mode'
300     diagUnits = 'dimensionless '
301     diagCode = 'SM MR '
302     diagMate = diagNum
303     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
304     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
305    
306     diagName = 'GM_URMS'
307     diagTitle = 'rms Eddy Velocity'
308     diagUnits = 'm/s '
309     diagCode = 'SM P MR '
310     diagMate = diagNum
311     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
312     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
313    
314     diagName = 'GM_UMC '
315     diagTitle = 'ubar-c'
316     diagUnits = 'm/s '
317     diagCode = 'SM MR '
318     diagMate = diagNum
319     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
320     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
321    
322     diagName = 'GM_SFLYR'
323     diagTitle = 'mixed layer depth'
324     diagUnits = 'm '
325     diagCode = 'SM P MR '
326     diagMate = diagNum
327     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
328     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
329    
330     diagName = 'GM_SIGR '
331     diagTitle = 'drho/dr'
332     diagUnits = 'kg/m**4 '
333     diagCode = 'SM MR '
334     diagMate = diagNum
335     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
336     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
337    
338     diagName = 'GM_N2 '
339     diagTitle = 'Square of buoyancy frequency'
340     diagUnits = '1/s**2 '
341     diagCode = 'SM MR '
342     diagMate = diagNum
343     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
344     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
345    
346     diagName = 'GM_SIGX '
347     diagTitle = 'd(rho)/dx'
348     diagUnits = 'kg/m**4 '
349     diagCode = 'UU MR '
350     diagMate = diagNum
351     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
352     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
353    
354     diagName = 'GM_SIGY '
355     diagTitle = 'd(rho)/dy'
356     diagUnits = 'kg/m**4 '
357     diagCode = 'VV MR '
358     diagMate = diagNum
359     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
360     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
361    
362     diagName = 'GM_USTAR'
363     diagTitle = 'u^*'
364     diagUnits = 'm/s '
365     diagCode = 'UU MR '
366     diagMate = diagNum
367     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
368     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
369    
370     diagName = 'GM_VSTAR'
371     diagTitle = 'v^*'
372     diagUnits = 'm/s '
373     diagCode = 'VV MR '
374     diagMate = diagNum
375     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
376     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
377    
378     diagName = 'GM_RDEF'
379     diagTitle = 'Deformation Radius'
380     diagUnits = 'm '
381     diagCode = 'SM P M1 '
382     diagMate = diagNum
383     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
384     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
385    
386     diagName = 'GM_LRHNS'
387     diagTitle = 'Rhines Radius'
388     diagUnits = 'm '
389     diagCode = 'SM P M1 '
390     diagMate = diagNum
391     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
392     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
393    
394     diagName = 'GM_RMIX'
395     diagTitle = 'Unmodulated Mixing Length'
396     diagUnits = 'm '
397     diagCode = 'SM P M1 '
398     diagMate = diagNum
399     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
400     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
401    
402     diagName = 'GM_SUPP'
403     diagTitle = 'Suppression Factor for K3D'
404     diagUnits = 'none '
405     diagCode = 'SM P MR '
406     diagMate = diagNum
407     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
408     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
409    
410     diagName = 'GM_dqdx '
411     diagTitle = 'dq/dx'
412     diagUnits = '1/(m*s) '
413     diagCode = 'UU MR '
414     diagMate = diagNum
415     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
416     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
417    
418     diagName = 'GM_dqdy '
419     diagTitle = 'dq/dy'
420     diagUnits = '1/(m*s) '
421     diagCode = 'VV MR '
422     diagMate = diagNum
423     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
424     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
425    
426     diagName = 'GM_Xix '
427     diagTitle = '-k dq/dx expansion'
428     diagUnits = '1/(m*s) '
429     diagCode = 'UU MR '
430     diagMate = diagNum
431     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
432     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
433    
434     diagName = 'GM_Xiy '
435     diagTitle = '-k dq/dy expansion'
436     diagUnits = '1/(m*s) '
437     diagCode = 'VV MR '
438     diagMate = diagNum
439     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
440     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
441    
442     diagName = 'GM_PSI '
443     diagTitle = 'Eddy induced meridional streamfunction'
444     diagUnits = 'm**3/s '
445     diagCode = 'VV MR '
446     diagMate = diagNum
447     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
448     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
449    
450     diagName = 'GM_Sx '
451     diagTitle = 'Zonal isopycnal slope'
452     diagUnits = 'none '
453     diagCode = 'UU LR '
454     diagMate = diagNum
455     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
456     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
457    
458     diagName = 'GM_Sy '
459     diagTitle = 'Meridional isopycnal slope'
460     diagUnits = 'none '
461     diagCode = 'VV LR '
462     diagMate = diagNum
463     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
464     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
465    
466     diagName = 'GM_TFLXX'
467     diagTitle = 'Zonal thickness flux'
468     diagUnits = '1/(m*s) '
469     diagCode = 'UU MR '
470     diagMate = diagNum
471     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
472     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
473    
474     diagName = 'GM_TFLXY'
475     diagTitle = 'meridional thickness flux'
476     diagUnits = '1/(m*s) '
477     diagCode = 'VV MR '
478     diagMate = diagNum
479     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
480     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
481    
482     diagName = 'GM_C'
483     diagTitle = 'Doppler shifted long Rossby wave speed'
484     diagUnits = 'm/s '
485     diagCode = 'SM M1 '
486     diagMate = diagNum
487     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
488     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
489    
490     diagName = 'GM_UBARO'
491     diagTitle = 'Barotropic velocity'
492     diagUnits = 'm/s '
493     diagCode = 'SM M1 '
494     diagMate = diagNum
495     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
496     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
497    
498     diagName = 'GM_EADY '
499     diagTitle = 'Eady Growth rate'
500     diagUnits = '1/s '
501     diagCode = 'SM M1 '
502     diagMate = diagNum
503     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
504     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
505    
506     diagName = 'GM_UBAR '
507     diagTitle = 'Mean zonal velocity'
508     diagUnits = 'm/s '
509     diagCode = 'SM MR '
510     diagMate = diagNum
511     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
512     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
513    
514     diagName = 'GM_KDEF '
515     diagTitle = 'Rossby wavenumber'
516     diagUnits = '1/m '
517     diagCode = 'SM M1 '
518     diagMate = diagNum
519     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
520     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
521    
522     diagName = 'GM_BVINT'
523     diagTitle = 'vertical integral of the buoyancy frequency'
524     diagUnits = 'm/s '
525     diagCode = 'SM M1 '
526     diagMate = diagNum
527     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
528     I diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
529    
530     #endif
531    
532 jmc 1.1 c ENDIF
533    
534     #endif /* ALLOW_DIAGNOSTICS */
535    
536     RETURN
537     END

  ViewVC Help
Powered by ViewVC 1.1.22