/[MITgcm]/MITgcm/model/inc/PARAMS.h
ViewVC logotype

Contents of /MITgcm/model/inc/PARAMS.h

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


Revision 1.289 - (show annotations) (download)
Thu Nov 2 17:43:26 2017 UTC (6 years, 4 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint66o, checkpoint66n, checkpoint66m, HEAD
Changes since 1.288: +24 -3 lines
File MIME type: text/plain
- new options 1) to account for true vertical distance (including hFac)
     in vertical viscous flux and diffusive flux ;
  2) to increase vertical mixing near surface and/or bottom where partial
     cell is too thin ;
- for now, both additions above are within: #ifndef EXCLUDE_PCELL_MIX_CODE ;
- also re-order options in CPP_OPTIONS.h (little bit more logical)

1 C $Header: /u/gcmpack/MITgcm/model/inc/PARAMS.h,v 1.288 2017/10/04 20:30:04 jmc Exp $
2 C $Name: $
3 C
4
5 CBOP
6 C !ROUTINE: PARAMS.h
7 C !INTERFACE:
8 C #include PARAMS.h
9
10 C !DESCRIPTION:
11 C Header file defining model "parameters". The values from the
12 C model standard input file are stored into the variables held
13 C here. Notes describing the parameters can also be found here.
14
15 CEOP
16
17 C-- Contants
18 C Useful physical values
19 Real*8 PI
20 PARAMETER ( PI = 3.14159265358979323844D0 )
21 Real*8 deg2rad
22 PARAMETER ( deg2rad = 2.D0*PI/360.D0 )
23
24 C-- COMMON /PARM_C/ Character valued parameters used by the model.
25 C buoyancyRelation :: Flag used to indicate which relation to use to
26 C get buoyancy.
27 C eosType :: choose the equation of state:
28 C LINEAR, POLY3, UNESCO, JMD95Z, JMD95P, MDJWF, IDEALGAS
29 C pickupSuff :: force to start from pickup files (even if nIter0=0)
30 C and read pickup files with this suffix (max 10 Char.)
31 C mdsioLocalDir :: read-write tiled file from/to this directory name
32 C (+ 4 digits Processor-Rank) instead of current dir.
33 C adTapeDir :: read-write checkpointing tape files from/to this
34 C directory name instead of current dir. Conflicts
35 C mdsioLocalDir, so only one of the two can be set.
36 C tRefFile :: File containing reference Potential Temperat. tRef (1.D)
37 C sRefFile :: File containing reference salinity/spec.humid. sRef (1.D)
38 C rhoRefFile :: File containing reference density profile rhoRef (1.D)
39 C gravityFile :: File containing gravity vertical profile (1.D)
40 C delRFile :: File containing vertical grid spacing delR (1.D array)
41 C delRcFile :: File containing vertical grid spacing delRc (1.D array)
42 C hybSigmFile :: File containing hybrid-sigma vertical coord. coeff. (2x 1.D)
43 C delXFile :: File containing X-spacing grid definition (1.D array)
44 C delYFile :: File containing Y-spacing grid definition (1.D array)
45 C horizGridFile :: File containing horizontal-grid definition
46 C (only when using curvilinear_grid)
47 C bathyFile :: File containing bathymetry. If not defined bathymetry
48 C is taken from inline function.
49 C topoFile :: File containing the topography of the surface (unit=m)
50 C (mainly used for the atmosphere = ground height).
51 C addWwallFile :: File containing 2-D additional Western cell-edge wall
52 C addSwallFile :: File containing 2-D additional Southern cell-edge wall
53 C (e.g., to add "thin-wall" where it is =1)
54 C hydrogThetaFile :: File containing initial hydrographic data (3-D)
55 C for potential temperature.
56 C hydrogSaltFile :: File containing initial hydrographic data (3-D)
57 C for salinity.
58 C diffKrFile :: File containing 3D specification of vertical diffusivity
59 C viscAhDfile :: File containing 3D specification of horizontal viscosity
60 C viscAhZfile :: File containing 3D specification of horizontal viscosity
61 C viscA4Dfile :: File containing 3D specification of horizontal viscosity
62 C viscA4Zfile :: File containing 3D specification of horizontal viscosity
63 C zonalWindFile :: File containing zonal wind data
64 C meridWindFile :: File containing meridional wind data
65 C thetaClimFile :: File containing surface theta climataology used
66 C in relaxation term -lambda(theta-theta*)
67 C saltClimFile :: File containing surface salt climataology used
68 C in relaxation term -lambda(salt-salt*)
69 C surfQfile :: File containing surface heat flux, excluding SW
70 C (old version, kept for backward compatibility)
71 C surfQnetFile :: File containing surface net heat flux
72 C surfQswFile :: File containing surface shortwave radiation
73 C EmPmRfile :: File containing surface fresh water flux
74 C NOTE: for backward compatibility EmPmRfile is specified in
75 C m/s when using external_fields_load.F. It is converted
76 C to kg/m2/s by multiplying by rhoConstFresh.
77 C saltFluxFile :: File containing surface salt flux
78 C pLoadFile :: File containing pressure loading
79 C addMassFile :: File containing source/sink of fluid in the interior
80 C eddyPsiXFile :: File containing zonal Eddy streamfunction data
81 C eddyPsiYFile :: File containing meridional Eddy streamfunction data
82 C the_run_name :: string identifying the name of the model "run"
83 COMMON /PARM_C/
84 & buoyancyRelation, eosType,
85 & pickupSuff, mdsioLocalDir, adTapeDir,
86 & tRefFile, sRefFile, rhoRefFile, gravityFile,
87 & delRFile, delRcFile, hybSigmFile,
88 & delXFile, delYFile, horizGridFile,
89 & bathyFile, topoFile, addWwallFile, addSwallFile,
90 & viscAhDfile, viscAhZfile,
91 & viscA4Dfile, viscA4Zfile,
92 & hydrogThetaFile, hydrogSaltFile, diffKrFile,
93 & zonalWindFile, meridWindFile, thetaClimFile,
94 & saltClimFile,
95 & EmPmRfile, saltFluxFile,
96 & surfQfile, surfQnetFile, surfQswFile,
97 & lambdaThetaFile, lambdaSaltFile,
98 & uVelInitFile, vVelInitFile, pSurfInitFile,
99 & pLoadFile, addMassFile,
100 & eddyPsiXFile, eddyPsiYFile, geothermalFile,
101 & the_run_name
102 CHARACTER*(MAX_LEN_FNAM) buoyancyRelation
103 CHARACTER*(6) eosType
104 CHARACTER*(10) pickupSuff
105 CHARACTER*(MAX_LEN_FNAM) mdsioLocalDir
106 CHARACTER*(MAX_LEN_FNAM) adTapeDir
107 CHARACTER*(MAX_LEN_FNAM) tRefFile
108 CHARACTER*(MAX_LEN_FNAM) sRefFile
109 CHARACTER*(MAX_LEN_FNAM) rhoRefFile
110 CHARACTER*(MAX_LEN_FNAM) gravityFile
111 CHARACTER*(MAX_LEN_FNAM) delRFile
112 CHARACTER*(MAX_LEN_FNAM) delRcFile
113 CHARACTER*(MAX_LEN_FNAM) hybSigmFile
114 CHARACTER*(MAX_LEN_FNAM) delXFile
115 CHARACTER*(MAX_LEN_FNAM) delYFile
116 CHARACTER*(MAX_LEN_FNAM) horizGridFile
117 CHARACTER*(MAX_LEN_FNAM) bathyFile, topoFile
118 CHARACTER*(MAX_LEN_FNAM) addWwallFile, addSwallFile
119 CHARACTER*(MAX_LEN_FNAM) hydrogThetaFile, hydrogSaltFile
120 CHARACTER*(MAX_LEN_FNAM) diffKrFile
121 CHARACTER*(MAX_LEN_FNAM) viscAhDfile
122 CHARACTER*(MAX_LEN_FNAM) viscAhZfile
123 CHARACTER*(MAX_LEN_FNAM) viscA4Dfile
124 CHARACTER*(MAX_LEN_FNAM) viscA4Zfile
125 CHARACTER*(MAX_LEN_FNAM) zonalWindFile
126 CHARACTER*(MAX_LEN_FNAM) meridWindFile
127 CHARACTER*(MAX_LEN_FNAM) thetaClimFile
128 CHARACTER*(MAX_LEN_FNAM) saltClimFile
129 CHARACTER*(MAX_LEN_FNAM) surfQfile
130 CHARACTER*(MAX_LEN_FNAM) surfQnetFile
131 CHARACTER*(MAX_LEN_FNAM) surfQswFile
132 CHARACTER*(MAX_LEN_FNAM) EmPmRfile
133 CHARACTER*(MAX_LEN_FNAM) saltFluxFile
134 CHARACTER*(MAX_LEN_FNAM) uVelInitFile
135 CHARACTER*(MAX_LEN_FNAM) vVelInitFile
136 CHARACTER*(MAX_LEN_FNAM) pSurfInitFile
137 CHARACTER*(MAX_LEN_FNAM) pLoadFile
138 CHARACTER*(MAX_LEN_FNAM) addMassFile
139 CHARACTER*(MAX_LEN_FNAM) eddyPsiXFile
140 CHARACTER*(MAX_LEN_FNAM) eddyPsiYFile
141 CHARACTER*(MAX_LEN_FNAM) geothermalFile
142 CHARACTER*(MAX_LEN_FNAM) lambdaThetaFile
143 CHARACTER*(MAX_LEN_FNAM) lambdaSaltFile
144 CHARACTER*(MAX_LEN_PREC/2) the_run_name
145
146 C-- COMMON /PARM_I/ Integer valued parameters used by the model.
147 C cg2dMaxIters :: Maximum number of iterations in the
148 C two-dimensional con. grad solver.
149 C cg2dChkResFreq :: Frequency with which to check residual
150 C in con. grad solver.
151 C cg2dPreCondFreq :: Frequency for updating cg2d preconditioner
152 C (non-linear free-surf.)
153 C cg2dUseMinResSol :: =0 : use last-iteration/converged solution
154 C =1 : use solver minimum-residual solution
155 C cg3dMaxIters :: Maximum number of iterations in the
156 C three-dimensional con. grad solver.
157 C cg3dChkResFreq :: Frequency with which to check residual
158 C in con. grad solver.
159 C printResidualFreq :: Frequency for printing residual in CG iterations
160 C nIter0 :: Start time-step number of for this run
161 C nTimeSteps :: Number of timesteps to execute
162 C nTimeSteps_l2 :: Number of inner timesteps to execute per timestep
163 C selectCoriMap :: select setting of Coriolis parameter map:
164 C =0 f-Plane (Constant Coriolis, = f0)
165 C =1 Beta-Plane Coriolis (= f0 + beta.y)
166 C =2 Spherical Coriolis (= 2.omega.sin(phi))
167 C =3 Read Coriolis 2-d fields from files.
168 C selectSigmaCoord :: option related to sigma vertical coordinate
169 C nonlinFreeSurf :: option related to non-linear free surface
170 C =0 Linear free surface ; >0 Non-linear
171 C select_rStar :: option related to r* vertical coordinate
172 C =0 (default) use r coord. ; > 0 use r*
173 C selectNHfreeSurf :: option for Non-Hydrostatic (free-)Surface formulation:
174 C =0 (default) hydrostatic surf. ; > 0 add NH effects.
175 C selectP_inEOS_Zc :: select which pressure to use in EOS (for z-coords)
176 C =0: simply: -g*rhoConst*z
177 C =1: use pRef = integral{-g*rho(Tref,Sref,pRef)*dz}
178 C =2: use hydrostatic dynamical pressure
179 C =3: use full (Hyd+NH) dynamical pressure
180 C selectAddFluid :: option to add mass source/sink of fluid in the interior
181 C (3-D generalisation of oceanic real-fresh water flux)
182 C =0 off ; =1 add fluid ; =-1 virtual flux (no mass added)
183 C selectImplicitDrag :: select Implicit treatment of bottom/top drag
184 C = 0: fully explicit
185 C = 1: implicit on provisional velocity
186 C (i.e., before grad.Eta increment)
187 C = 2: fully implicit (combined with Impl Surf.Press)
188 C momForcingOutAB :: =1: take momentum forcing contribution
189 C out of (=0: in) Adams-Bashforth time stepping.
190 C tracForcingOutAB :: =1: take tracer (Temp,Salt,pTracers) forcing contribution
191 C out of (=0: in) Adams-Bashforth time stepping.
192 C tempAdvScheme :: Temp. Horiz.Advection scheme selector
193 C tempVertAdvScheme :: Temp. Vert. Advection scheme selector
194 C saltAdvScheme :: Salt. Horiz.advection scheme selector
195 C saltVertAdvScheme :: Salt. Vert. Advection scheme selector
196 C selectKEscheme :: Kinetic Energy scheme selector (Vector Inv.)
197 C selectVortScheme :: Scheme selector for Vorticity term (Vector Inv.)
198 C selectBotDragQuadr :: quadratic bottom drag discretisation option:
199 C =0: average KE from grid center to U & V location
200 C =1: use local velocity norm @ U & V location
201 C =2: same with wet-point averaging of other component
202 C pCellMix_select :: select option to enhance mixing near surface & bottom
203 C unit digit: near bottom ; tens digit: near surface
204 C with digit =0 : disable ;
205 C = 1 : increases mixing linearly with recip_hFac
206 C = 2,3,4 : increases mixing by recip_hFac^(2,3,4)
207 C readBinaryPrec :: Precision used for reading binary files
208 C writeStatePrec :: Precision used for writing model state.
209 C writeBinaryPrec :: Precision used for writing binary files
210 C rwSuffixType :: controls the format of the mds file suffix.
211 C =0 (default): use iteration number (myIter, I10.10);
212 C =1: 100*myTime (100th sec); =2: myTime (seconds);
213 C =3: myTime/360 (10th of hr); =4: myTime/3600 (hours).
214 C monitorSelect :: select group of variables to monitor
215 C =1 : dynvars ; =2 : + vort ; =3 : + surface
216 C- debugLevel :: controls printing of algorithm intermediate results
217 C and statistics ; higher -> more writing
218 C- plotLevel :: controls printing of field maps ; higher -> more flds
219
220 COMMON /PARM_I/
221 & cg2dMaxIters, cg2dChkResFreq,
222 & cg2dPreCondFreq, cg2dUseMinResSol,
223 & cg3dMaxIters, cg3dChkResFreq,
224 & printResidualFreq,
225 & nIter0, nTimeSteps, nTimeSteps_l2, nEndIter,
226 & selectCoriMap,
227 & selectSigmaCoord,
228 & nonlinFreeSurf, select_rStar,
229 & selectNHfreeSurf, selectP_inEOS_Zc,
230 & selectAddFluid, selectImplicitDrag,
231 & momForcingOutAB, tracForcingOutAB,
232 & tempAdvScheme, tempVertAdvScheme,
233 & saltAdvScheme, saltVertAdvScheme,
234 & selectKEscheme, selectVortScheme,
235 & selectBotDragQuadr, pCellMix_select,
236 & readBinaryPrec, writeBinaryPrec, writeStatePrec,
237 & rwSuffixType, monitorSelect, debugLevel, plotLevel
238 INTEGER cg2dMaxIters
239 INTEGER cg2dChkResFreq
240 INTEGER cg2dPreCondFreq
241 INTEGER cg2dUseMinResSol
242 INTEGER cg3dMaxIters
243 INTEGER cg3dChkResFreq
244 INTEGER printResidualFreq
245 INTEGER nIter0
246 INTEGER nTimeSteps
247 INTEGER nTimeSteps_l2
248 INTEGER nEndIter
249 INTEGER selectCoriMap
250 INTEGER selectSigmaCoord
251 INTEGER nonlinFreeSurf
252 INTEGER select_rStar
253 INTEGER selectNHfreeSurf
254 INTEGER selectP_inEOS_Zc
255 INTEGER selectAddFluid
256 INTEGER selectImplicitDrag
257 INTEGER momForcingOutAB, tracForcingOutAB
258 INTEGER tempAdvScheme, tempVertAdvScheme
259 INTEGER saltAdvScheme, saltVertAdvScheme
260 INTEGER selectKEscheme
261 INTEGER selectVortScheme
262 INTEGER selectBotDragQuadr
263 INTEGER pCellMix_select
264 INTEGER readBinaryPrec
265 INTEGER writeStatePrec
266 INTEGER writeBinaryPrec
267 INTEGER rwSuffixType
268 INTEGER monitorSelect
269 INTEGER debugLevel
270 INTEGER plotLevel
271
272 C-- COMMON /PARM_L/ Logical valued parameters used by the model.
273 C- Coordinate + Grid params:
274 C fluidIsAir :: Set to indicate that the fluid major constituent
275 C is air
276 C fluidIsWater :: Set to indicate that the fluid major constituent
277 C is water
278 C usingPCoords :: Set to indicate that we are working in a pressure
279 C type coordinate (p or p*).
280 C usingZCoords :: Set to indicate that we are working in a height
281 C type coordinate (z or z*)
282 C usingCartesianGrid :: If TRUE grid generation will be in a cartesian
283 C coordinate frame.
284 C usingSphericalPolarGrid :: If TRUE grid generation will be in a
285 C spherical polar frame.
286 C rotateGrid :: rotate grid coordinates to geographical coordinates
287 C according to Euler angles phiEuler, thetaEuler, psiEuler
288 C usingCylindricalGrid :: If TRUE grid generation will be Cylindrical
289 C usingCurvilinearGrid :: If TRUE, use a curvilinear grid (to be provided)
290 C hasWetCSCorners :: domain contains CS-type corners where dynamics is solved
291 C deepAtmosphere :: deep model (drop the shallow-atmosphere approximation)
292 C setInterFDr :: set Interface depth (put cell-Center at the middle)
293 C setCenterDr :: set cell-Center depth (put Interface at the middle)
294 C useMin4hFacEdges :: set hFacW,hFacS as minimum of adjacent hFacC factor
295 C interViscAr_pCell :: account for partial-cell in interior vert. viscosity
296 C interDiffKr_pCell :: account for partial-cell in interior vert. diffusion
297 C- Momentum params:
298 C no_slip_sides :: Impose "no-slip" at lateral boundaries.
299 C no_slip_bottom :: Impose "no-slip" at bottom boundary.
300 C bottomVisc_pCell :: account for partial-cell in bottom visc. (no-slip BC)
301 C useSmag3D :: Use isotropic 3-D Smagorinsky
302 C useFullLeith :: Set to true to use full Leith viscosity(may be unstable
303 C on irregular grids)
304 C useStrainTensionVisc:: Set to true to use Strain-Tension viscous terms
305 C useAreaViscLength :: Set to true to use old scaling for viscous lengths,
306 C e.g., L2=Raz. May be preferable for cube sphere.
307 C momViscosity :: Flag which turns momentum friction terms on and off.
308 C momAdvection :: Flag which turns advection of momentum on and off.
309 C momForcing :: Flag which turns external forcing of momentum on and off.
310 C momTidalForcing :: Flag which turns tidal forcing on and off.
311 C momPressureForcing :: Flag which turns pressure term in momentum equation
312 C on and off.
313 C metricTerms :: Flag which turns metric terms on or off.
314 C useNHMTerms :: If TRUE use non-hydrostatic metric terms.
315 C useCoriolis :: Flag which turns the coriolis terms on and off.
316 C use3dCoriolis :: Turns the 3-D coriolis terms (in Omega.cos Phi) on - off
317 C useCDscheme :: use CD-scheme to calculate Coriolis terms.
318 C vectorInvariantMomentum :: use Vector-Invariant form (mom_vecinv package)
319 C (default = F = use mom_fluxform package)
320 C useJamartWetPoints :: Use wet-point method for Coriolis (Jamart & Ozer 1986)
321 C useJamartMomAdv :: Use wet-point method for V.I. non-linear term
322 C upwindVorticity :: bias interpolation of vorticity in the Coriolis term
323 C highOrderVorticity :: use 3rd/4th order interp. of vorticity (V.I., advection)
324 C useAbsVorticity :: work with f+zeta in Coriolis terms
325 C upwindShear :: use 1rst order upwind interp. (V.I., vertical advection)
326 C momStepping :: Turns momentum equation time-stepping off
327 C calc_wVelocity :: Turns vertical velocity calculation off
328 C- Temp. & Salt params:
329 C tempStepping :: Turns temperature equation time-stepping on/off
330 C saltStepping :: Turns salinity equation time-stepping on/off
331 C addFrictionHeating :: account for frictional heating
332 C tempAdvection :: Flag which turns advection of temperature on and off.
333 C tempVertDiff4 :: use vertical bi-harmonic diffusion for temperature
334 C tempIsActiveTr :: Pot.Temp. is a dynamically active tracer
335 C tempForcing :: Flag which turns external forcing of temperature on/off
336 C saltAdvection :: Flag which turns advection of salinity on and off.
337 C saltVertDiff4 :: use vertical bi-harmonic diffusion for salinity
338 C saltIsActiveTr :: Salinity is a dynamically active tracer
339 C saltForcing :: Flag which turns external forcing of salinity on/off
340 C maskIniTemp :: apply mask to initial Pot.Temp.
341 C maskIniSalt :: apply mask to initial salinity
342 C checkIniTemp :: check for points with identically zero initial Pot.Temp.
343 C checkIniSalt :: check for points with identically zero initial salinity
344 C- Pressure solver related parameters (PARM02)
345 C useSRCGSolver :: Set to true to use conjugate gradient
346 C solver with single reduction (only one call of
347 C s/r mpi_allreduce), default is false
348 C- Time-stepping & free-surface params:
349 C rigidLid :: Set to true to use rigid lid
350 C implicitFreeSurface :: Set to true to use implicit free surface
351 C uniformLin_PhiSurf :: Set to true to use a uniform Bo_surf in the
352 C linear relation Phi_surf = Bo_surf*eta
353 C uniformFreeSurfLev :: TRUE if free-surface level-index is uniform (=1)
354 C exactConserv :: Set to true to conserve exactly the total Volume
355 C linFSConserveTr :: Set to true to correct source/sink of tracer
356 C at the surface due to Linear Free Surface
357 C useRealFreshWaterFlux :: if True (=Natural BCS), treats P+R-E flux
358 C as a real Fresh Water (=> changes the Sea Level)
359 C if F, converts P+R-E to salt flux (no SL effect)
360 C storePhiHyd4Phys :: store hydrostatic potential for use in Physics/EOS
361 C this requires specific code for restart & exchange
362 C quasiHydrostatic :: Using non-hydrostatic terms in hydrostatic algorithm
363 C nonHydrostatic :: Using non-hydrostatic algorithm
364 C use3Dsolver :: set to true to use 3-D pressure solver
365 C implicitIntGravWave :: treat Internal Gravity Wave implicitly
366 C staggerTimeStep :: enable a Stagger time stepping U,V (& W) then T,S
367 C applyExchUV_early :: Apply EXCH to U,V earlier, just before integr_continuity
368 C doResetHFactors :: Do reset thickness factors @ beginning of each time-step
369 C implicitDiffusion :: Turns implicit vertical diffusion on
370 C implicitViscosity :: Turns implicit vertical viscosity on
371 C tempImplVertAdv :: Turns on implicit vertical advection for Temperature
372 C saltImplVertAdv :: Turns on implicit vertical advection for Salinity
373 C momImplVertAdv :: Turns on implicit vertical advection for Momentum
374 C multiDimAdvection :: Flag that enable multi-dimension advection
375 C useMultiDimAdvec :: True if multi-dim advection is used at least once
376 C momDissip_In_AB :: if False, put Dissipation tendency contribution
377 C out off Adams-Bashforth time stepping.
378 C doAB_onGtGs :: if the Adams-Bashforth time stepping is used, always
379 C apply AB on tracer tendencies (rather than on Tracer)
380 C- Other forcing params -
381 C balanceEmPmR :: substract global mean of EmPmR at every time step
382 C balanceQnet :: substract global mean of Qnet at every time step
383 C balancePrintMean:: print substracted global means to STDOUT
384 C doThetaClimRelax :: Set true if relaxation to temperature
385 C climatology is required.
386 C doSaltClimRelax :: Set true if relaxation to salinity
387 C climatology is required.
388 C balanceThetaClimRelax :: substract global mean effect at every time step
389 C balanceSaltClimRelax :: substract global mean effect at every time step
390 C allowFreezing :: Allows surface water to freeze and form ice
391 C periodicExternalForcing :: Set true if forcing is time-dependant
392 C- I/O parameters -
393 C globalFiles :: Selects between "global" and "tiled" files.
394 C On some platforms with MPI, option globalFiles is either
395 C slow or does not work. Use useSingleCpuIO instead.
396 C useSingleCpuIO :: moved to EEPARAMS.h
397 C pickupStrictlyMatch :: check and stop if pickup-file do not stricly match
398 C startFromPickupAB2 :: with AB-3 code, start from an AB-2 pickup
399 C usePickupBeforeC54 :: start from old-pickup files, generated with code from
400 C before checkpoint-54a, Jul 06, 2004.
401 C pickup_write_mdsio :: use mdsio to write pickups
402 C pickup_read_mdsio :: use mdsio to read pickups
403 C pickup_write_immed :: echo the pickup immediately (for conversion)
404 C writePickupAtEnd :: write pickup at the last timestep
405 C timeave_mdsio :: use mdsio for timeave output
406 C snapshot_mdsio :: use mdsio for "snapshot" (dumpfreq/diagfreq) output
407 C monitor_stdio :: use stdio for monitor output
408 C dumpInitAndLast :: dumps model state to files at Initial (nIter0)
409 C & Last iteration, in addition multiple of dumpFreq iter.
410
411 COMMON /PARM_L/
412 & fluidIsAir, fluidIsWater,
413 & usingPCoords, usingZCoords,
414 & usingCartesianGrid, usingSphericalPolarGrid, rotateGrid,
415 & usingCylindricalGrid, usingCurvilinearGrid, hasWetCSCorners,
416 & deepAtmosphere, setInterFDr, setCenterDr, useMin4hFacEdges,
417 & interViscAr_pCell, interDiffKr_pCell,
418 & no_slip_sides, no_slip_bottom, bottomVisc_pCell, useSmag3D,
419 & useFullLeith, useStrainTensionVisc, useAreaViscLength,
420 & momViscosity, momAdvection, momForcing, momTidalForcing,
421 & momPressureForcing, metricTerms, useNHMTerms,
422 & useCoriolis, use3dCoriolis,
423 & useCDscheme, vectorInvariantMomentum,
424 & useEnergyConservingCoriolis, useJamartWetPoints, useJamartMomAdv,
425 & upwindVorticity, highOrderVorticity,
426 & useAbsVorticity, upwindShear,
427 & momStepping, calc_wVelocity, tempStepping, saltStepping,
428 & addFrictionHeating,
429 & tempAdvection, tempVertDiff4, tempIsActiveTr, tempForcing,
430 & saltAdvection, saltVertDiff4, saltIsActiveTr, saltForcing,
431 & maskIniTemp, maskIniSalt, checkIniTemp, checkIniSalt,
432 & useSRCGSolver,
433 & rigidLid, implicitFreeSurface,
434 & uniformLin_PhiSurf, uniformFreeSurfLev,
435 & exactConserv, linFSConserveTr, useRealFreshWaterFlux,
436 & storePhiHyd4Phys, quasiHydrostatic, nonHydrostatic,
437 & use3Dsolver, implicitIntGravWave, staggerTimeStep,
438 & applyExchUV_early, doResetHFactors,
439 & implicitDiffusion, implicitViscosity,
440 & tempImplVertAdv, saltImplVertAdv, momImplVertAdv,
441 & multiDimAdvection, useMultiDimAdvec,
442 & momDissip_In_AB, doAB_onGtGs,
443 & balanceEmPmR, balanceQnet, balancePrintMean,
444 & balanceThetaClimRelax, balanceSaltClimRelax,
445 & doThetaClimRelax, doSaltClimRelax,
446 & allowFreezing,
447 & periodicExternalForcing,
448 & globalFiles,
449 & pickupStrictlyMatch, usePickupBeforeC54, startFromPickupAB2,
450 & pickup_read_mdsio, pickup_write_mdsio, pickup_write_immed,
451 & writePickupAtEnd,
452 & timeave_mdsio, snapshot_mdsio, monitor_stdio,
453 & outputTypesInclusive, dumpInitAndLast
454
455 LOGICAL fluidIsAir
456 LOGICAL fluidIsWater
457 LOGICAL usingPCoords
458 LOGICAL usingZCoords
459 LOGICAL usingCartesianGrid
460 LOGICAL usingSphericalPolarGrid, rotateGrid
461 LOGICAL usingCylindricalGrid
462 LOGICAL usingCurvilinearGrid, hasWetCSCorners
463 LOGICAL deepAtmosphere
464 LOGICAL setInterFDr
465 LOGICAL setCenterDr
466 LOGICAL useMin4hFacEdges
467 LOGICAL interViscAr_pCell
468 LOGICAL interDiffKr_pCell
469
470 LOGICAL no_slip_sides
471 LOGICAL no_slip_bottom
472 LOGICAL bottomVisc_pCell
473 LOGICAL useSmag3D
474 LOGICAL useFullLeith
475 LOGICAL useStrainTensionVisc
476 LOGICAL useAreaViscLength
477 LOGICAL momViscosity
478 LOGICAL momAdvection
479 LOGICAL momForcing
480 LOGICAL momTidalForcing
481 LOGICAL momPressureForcing
482 LOGICAL metricTerms
483 LOGICAL useNHMTerms
484
485 LOGICAL useCoriolis
486 LOGICAL use3dCoriolis
487 LOGICAL useCDscheme
488 LOGICAL vectorInvariantMomentum
489 LOGICAL useEnergyConservingCoriolis
490 LOGICAL useJamartWetPoints
491 LOGICAL useJamartMomAdv
492 LOGICAL upwindVorticity
493 LOGICAL highOrderVorticity
494 LOGICAL useAbsVorticity
495 LOGICAL upwindShear
496 LOGICAL momStepping
497 LOGICAL calc_wVelocity
498 LOGICAL tempStepping
499 LOGICAL saltStepping
500 LOGICAL addFrictionHeating
501 LOGICAL tempAdvection
502 LOGICAL tempVertDiff4
503 LOGICAL tempIsActiveTr
504 LOGICAL tempForcing
505 LOGICAL saltAdvection
506 LOGICAL saltVertDiff4
507 LOGICAL saltIsActiveTr
508 LOGICAL saltForcing
509 LOGICAL maskIniTemp
510 LOGICAL maskIniSalt
511 LOGICAL checkIniTemp
512 LOGICAL checkIniSalt
513 LOGICAL useSRCGSolver
514 LOGICAL rigidLid
515 LOGICAL implicitFreeSurface
516 LOGICAL uniformLin_PhiSurf
517 LOGICAL uniformFreeSurfLev
518 LOGICAL exactConserv
519 LOGICAL linFSConserveTr
520 LOGICAL useRealFreshWaterFlux
521 LOGICAL storePhiHyd4Phys
522 LOGICAL quasiHydrostatic
523 LOGICAL nonHydrostatic
524 LOGICAL use3Dsolver
525 LOGICAL implicitIntGravWave
526 LOGICAL staggerTimeStep
527 LOGICAL applyExchUV_early
528 LOGICAL doResetHFactors
529 LOGICAL implicitDiffusion
530 LOGICAL implicitViscosity
531 LOGICAL tempImplVertAdv
532 LOGICAL saltImplVertAdv
533 LOGICAL momImplVertAdv
534 LOGICAL multiDimAdvection
535 LOGICAL useMultiDimAdvec
536 LOGICAL momDissip_In_AB
537 LOGICAL doAB_onGtGs
538 LOGICAL balanceEmPmR
539 LOGICAL balanceQnet
540 LOGICAL balancePrintMean
541 LOGICAL doThetaClimRelax
542 LOGICAL doSaltClimRelax
543 LOGICAL balanceThetaClimRelax
544 LOGICAL balanceSaltClimRelax
545 LOGICAL allowFreezing
546 LOGICAL periodicExternalForcing
547 LOGICAL globalFiles
548 LOGICAL pickupStrictlyMatch
549 LOGICAL usePickupBeforeC54
550 LOGICAL startFromPickupAB2
551 LOGICAL pickup_read_mdsio, pickup_write_mdsio
552 LOGICAL pickup_write_immed, writePickupAtEnd
553 LOGICAL timeave_mdsio, snapshot_mdsio, monitor_stdio
554 LOGICAL outputTypesInclusive
555 LOGICAL dumpInitAndLast
556
557 C-- COMMON /PARM_R/ "Real" valued parameters used by the model.
558 C cg2dTargetResidual
559 C :: Target residual for cg2d solver; no unit (RHS normalisation)
560 C cg2dTargetResWunit
561 C :: Target residual for cg2d solver; W unit (No RHS normalisation)
562 C cg3dTargetResidual
563 C :: Target residual for cg3d solver.
564 C cg2dpcOffDFac :: Averaging weight for preconditioner off-diagonal.
565 C Note. 20th May 1998
566 C I made a weird discovery! In the model paper we argue
567 C for the form of the preconditioner used here ( see
568 C A Finite-volume, Incompressible Navier-Stokes Model
569 C ...., Marshall et. al ). The algebra gives a simple
570 C 0.5 factor for the averaging of ac and aCw to get a
571 C symmettric pre-conditioner. By using a factor of 0.51
572 C i.e. scaling the off-diagonal terms in the
573 C preconditioner down slightly I managed to get the
574 C number of iterations for convergence in a test case to
575 C drop form 192 -> 134! Need to investigate this further!
576 C For now I have introduced a parameter cg2dpcOffDFac which
577 C defaults to 0.51 but can be set at runtime.
578 C delR :: Vertical grid spacing ( units of r ).
579 C delRc :: Vertical grid spacing between cell centers (r unit).
580 C delX :: Separation between cell faces (m) or (deg), depending
581 C delY on input flags. Note: moved to header file SET_GRID.h
582 C xgOrigin :: Origin of the X-axis (Cartesian Grid) / Longitude of Western
583 C :: most cell face (Lat-Lon grid) (Note: this is an "inert"
584 C :: parameter but it makes geographical references simple.)
585 C ygOrigin :: Origin of the Y-axis (Cartesian Grid) / Latitude of Southern
586 C :: most face (Lat-Lon grid).
587 C rSphere :: Radius of sphere for a spherical polar grid ( m ).
588 C recip_rSphere :: Reciprocal radius of sphere ( m^-1 ).
589 C radius_fromHorizGrid :: sphere Radius of input horiz. grid (Curvilinear Grid)
590 C seaLev_Z :: the reference height of sea-level (usually zero)
591 C top_Pres :: pressure (P-Coords) or reference pressure (Z-Coords) at the top
592 C rSigmaBnd :: vertical position (in r-unit) of r/sigma transition (Hybrid-Sigma)
593 C gravity :: Acceleration due to constant gravity ( m/s^2 )
594 C recip_gravity :: Reciprocal gravity acceleration ( s^2/m )
595 C gBaro :: Accel. due to gravity used in barotropic equation ( m/s^2 )
596 C gravFacC :: gravity factor (vs surf. gravity) vert. profile at cell-Center
597 C gravFacF :: gravity factor (vs surf. gravity) vert. profile at cell-interF
598 C rhoNil :: Reference density for the linear equation of state
599 C rhoConst :: Vertically constant reference density (Boussinesq)
600 C rho1Ref :: reference vertical profile for density (anelastic)
601 C rhoFacC :: normalized (by rhoConst) reference density at cell-Center
602 C rhoFacF :: normalized (by rhoConst) reference density at cell-interFace
603 C rhoConstFresh :: Constant reference density for fresh water (rain)
604 C thetaConst :: Constant reference for potential temperature
605 C tRef :: reference vertical profile for potential temperature
606 C sRef :: reference vertical profile for salinity/specific humidity
607 C pRef4EOS :: reference pressure used in EOS (case selectP_inEOS_Zc=1)
608 C phiRef :: reference potential (press/rho, geopot) profile (m^2/s^2)
609 C dBdrRef :: vertical gradient of reference buoyancy [(m/s/r)^2]:
610 C :: z-coord: = N^2_ref = Brunt-Vaissala frequency [s^-2]
611 C :: p-coord: = -(d.alpha/dp)_ref [(m^2.s/kg)^2]
612 C rVel2wUnit :: units conversion factor (Non-Hydrostatic code),
613 C :: from r-coordinate vertical velocity to vertical velocity [m/s].
614 C :: z-coord: = 1 ; p-coord: wSpeed [m/s] = rVel [Pa/s] * rVel2wUnit
615 C wUnit2rVel :: units conversion factor (Non-Hydrostatic code),
616 C :: from vertical velocity [m/s] to r-coordinate vertical velocity.
617 C :: z-coord: = 1 ; p-coord: rVel [Pa/s] = wSpeed [m/s] * wUnit2rVel
618 C mass2rUnit :: units conversion factor (surface forcing),
619 C :: from mass per unit area [kg/m2] to vertical r-coordinate unit.
620 C :: z-coord: = 1/rhoConst ( [kg/m2] / rho = [m] ) ;
621 C :: p-coord: = gravity ( [kg/m2] * g = [Pa] ) ;
622 C rUnit2mass :: units conversion factor (surface forcing),
623 C :: from vertical r-coordinate unit to mass per unit area [kg/m2].
624 C :: z-coord: = rhoConst ( [m] * rho = [kg/m2] ) ;
625 C :: p-coord: = 1/gravity ( [Pa] / g = [kg/m2] ) ;
626 C f0 :: Reference coriolis parameter ( 1/s )
627 C ( Southern edge f for beta plane )
628 C beta :: df/dy ( s^-1.m^-1 )
629 C fPrime :: Second Coriolis parameter ( 1/s ), related to Y-component
630 C of rotation (reference value = 2.Omega.Cos(Phi))
631 C omega :: Angular velocity ( rad/s )
632 C rotationPeriod :: Rotation period (s) (= 2.pi/omega)
633 C viscArNr :: vertical profile of Eddy viscosity coeff.
634 C for vertical mixing of momentum ( units of r^2/s )
635 C viscAh :: Eddy viscosity coeff. for mixing of
636 C momentum laterally ( m^2/s )
637 C viscAhW :: Eddy viscosity coeff. for mixing of vertical
638 C momentum laterally, no effect for hydrostatic
639 C model, defaults to viscAhD if unset ( m^2/s )
640 C Not used if variable horiz. viscosity is used.
641 C viscA4 :: Biharmonic viscosity coeff. for mixing of
642 C momentum laterally ( m^4/s )
643 C viscA4W :: Biharmonic viscosity coeff. for mixing of vertical
644 C momentum laterally, no effect for hydrostatic
645 C model, defaults to viscA4D if unset ( m^2/s )
646 C Not used if variable horiz. viscosity is used.
647 C viscAhD :: Eddy viscosity coeff. for mixing of momentum laterally
648 C (act on Divergence part) ( m^2/s )
649 C viscAhZ :: Eddy viscosity coeff. for mixing of momentum laterally
650 C (act on Vorticity part) ( m^2/s )
651 C viscA4D :: Biharmonic viscosity coeff. for mixing of momentum laterally
652 C (act on Divergence part) ( m^4/s )
653 C viscA4Z :: Biharmonic viscosity coeff. for mixing of momentum laterally
654 C (act on Vorticity part) ( m^4/s )
655 C smag3D_coeff :: Isotropic 3-D Smagorinsky coefficient (-)
656 C viscC2leith :: Leith non-dimensional viscosity factor (grad(vort))
657 C viscC2leithD :: Modified Leith non-dimensional visc. factor (grad(div))
658 C viscC4leith :: Leith non-dimensional viscosity factor (grad(vort))
659 C viscC4leithD :: Modified Leith non-dimensional viscosity factor (grad(div))
660 C viscC2smag :: Smagorinsky non-dimensional viscosity factor (harmonic)
661 C viscC4smag :: Smagorinsky non-dimensional viscosity factor (biharmonic)
662 C viscAhMax :: Maximum eddy viscosity coeff. for mixing of
663 C momentum laterally ( m^2/s )
664 C viscAhReMax :: Maximum gridscale Reynolds number for eddy viscosity
665 C coeff. for mixing of momentum laterally (non-dim)
666 C viscAhGrid :: non-dimensional grid-size dependent viscosity
667 C viscAhGridMax:: maximum and minimum harmonic viscosity coefficients ...
668 C viscAhGridMin:: in terms of non-dimensional grid-size dependent visc.
669 C viscA4Max :: Maximum biharmonic viscosity coeff. for mixing of
670 C momentum laterally ( m^4/s )
671 C viscA4ReMax :: Maximum Gridscale Reynolds number for
672 C biharmonic viscosity coeff. momentum laterally (non-dim)
673 C viscA4Grid :: non-dimensional grid-size dependent bi-harmonic viscosity
674 C viscA4GridMax:: maximum and minimum biharmonic viscosity coefficients ...
675 C viscA4GridMin:: in terms of non-dimensional grid-size dependent viscosity
676 C diffKhT :: Laplacian diffusion coeff. for mixing of
677 C heat laterally ( m^2/s )
678 C diffK4T :: Biharmonic diffusion coeff. for mixing of
679 C heat laterally ( m^4/s )
680 C diffKrNrT :: vertical profile of Laplacian diffusion coeff.
681 C for mixing of heat vertically ( units of r^2/s )
682 C diffKr4T :: vertical profile of Biharmonic diffusion coeff.
683 C for mixing of heat vertically ( units of r^4/s )
684 C diffKhS :: Laplacian diffusion coeff. for mixing of
685 C salt laterally ( m^2/s )
686 C diffK4S :: Biharmonic diffusion coeff. for mixing of
687 C salt laterally ( m^4/s )
688 C diffKrNrS :: vertical profile of Laplacian diffusion coeff.
689 C for mixing of salt vertically ( units of r^2/s ),
690 C diffKr4S :: vertical profile of Biharmonic diffusion coeff.
691 C for mixing of salt vertically ( units of r^4/s )
692 C diffKrBL79surf :: T/S surface diffusivity (m^2/s) Bryan and Lewis, 1979
693 C diffKrBL79deep :: T/S deep diffusivity (m^2/s) Bryan and Lewis, 1979
694 C diffKrBL79scl :: depth scale for arctan fn (m) Bryan and Lewis, 1979
695 C diffKrBL79Ho :: depth offset for arctan fn (m) Bryan and Lewis, 1979
696 C BL79LatVary :: polarwise of this latitude diffKrBL79 is applied with
697 C gradual transition to diffKrBLEQ towards Equator
698 C diffKrBLEQsurf :: same as diffKrBL79surf but at Equator
699 C diffKrBLEQdeep :: same as diffKrBL79deep but at Equator
700 C diffKrBLEQscl :: same as diffKrBL79scl but at Equator
701 C diffKrBLEQHo :: same as diffKrBL79Ho but at Equator
702 C pCellMix_maxFac :: maximum enhanced mixing factor for thin partial-cell
703 C pCellMix_delR :: thickness criteria for too thin partial-cell
704 C pCellMix_viscAr :: vertical viscosity for too thin partial-cell
705 C pCellMix_diffKr :: vertical diffusivity for too thin partial-cell
706 C deltaT :: Default timestep ( s )
707 C deltaTClock :: Timestep used as model "clock". This determines the
708 C IO frequencies and is used in tagging output. It can
709 C be totally different to the dynamical time. Typically
710 C it will be the deep-water timestep for accelerated runs.
711 C Frequency of checkpointing and dumping of the model state
712 C are referenced to this clock. ( s )
713 C deltaTMom :: Timestep for momemtum equations ( s )
714 C dTtracerLev :: Timestep for tracer equations ( s ), function of level k
715 C deltaTFreeSurf :: Timestep for free-surface equation ( s )
716 C freeSurfFac :: Parameter to turn implicit free surface term on or off
717 C freeSurFac = 1. uses implicit free surface
718 C freeSurFac = 0. uses rigid lid
719 C abEps :: Adams-Bashforth-2 stabilizing weight
720 C alph_AB :: Adams-Bashforth-3 primary factor
721 C beta_AB :: Adams-Bashforth-3 secondary factor
722 C implicSurfPress :: parameter of the Crank-Nickelson time stepping :
723 C Implicit part of Surface Pressure Gradient ( 0-1 )
724 C implicDiv2DFlow :: parameter of the Crank-Nickelson time stepping :
725 C Implicit part of barotropic flow Divergence ( 0-1 )
726 C implicitNHPress :: parameter of the Crank-Nickelson time stepping :
727 C Implicit part of Non-Hydrostatic Pressure Gradient ( 0-1 )
728 C hFacMin :: Minimum fraction size of a cell (affects hFacC etc...)
729 C hFacMinDz :: Minimum dimensional size of a cell (affects hFacC etc..., m)
730 C hFacMinDp :: Minimum dimensional size of a cell (affects hFacC etc..., Pa)
731 C hFacMinDr :: Minimum dimensional size of a cell (-> hFacC etc..., r units)
732 C hFacInf :: Threshold (inf and sup) for fraction size of surface cell
733 C hFacSup that control vanishing and creating levels
734 C tauCD :: CD scheme coupling timescale ( s )
735 C rCD :: CD scheme normalised coupling parameter (= 1 - deltaT/tauCD)
736 C epsAB_CD :: Adams-Bashforth-2 stabilizing weight used in CD scheme
737 C baseTime :: model base time (time origin) = time @ iteration zero
738 C startTime :: Starting time for this integration ( s ).
739 C endTime :: Ending time for this integration ( s ).
740 C chkPtFreq :: Frequency of rolling check pointing ( s ).
741 C pChkPtFreq :: Frequency of permanent check pointing ( s ).
742 C dumpFreq :: Frequency with which model state is written to
743 C post-processing files ( s ).
744 C diagFreq :: Frequency with which model writes diagnostic output
745 C of intermediate quantities.
746 C afFacMom :: Advection of momentum term tracer parameter
747 C vfFacMom :: Momentum viscosity tracer parameter
748 C pfFacMom :: Momentum pressure forcing tracer parameter
749 C cfFacMom :: Coriolis term tracer parameter
750 C foFacMom :: Momentum forcing tracer parameter
751 C mtFacMom :: Metric terms tracer parameter
752 C cosPower :: Power of cosine of latitude to multiply viscosity
753 C cAdjFreq :: Frequency of convective adjustment
754 C
755 C taveFreq :: Frequency with which time-averaged model state
756 C is written to post-processing files ( s ).
757 C tave_lastIter :: (for state variable only) fraction of the last time
758 C step (of each taveFreq period) put in the time average.
759 C (fraction for 1rst iter = 1 - tave_lastIter)
760 C tauThetaClimRelax :: Relaxation to climatology time scale ( s ).
761 C tauSaltClimRelax :: Relaxation to climatology time scale ( s ).
762 C latBandClimRelax :: latitude band where Relaxation to Clim. is applied,
763 C i.e. where |yC| <= latBandClimRelax
764 C externForcingPeriod :: Is the period of which forcing varies (eg. 1 month)
765 C externForcingCycle :: Is the repeat time of the forcing (eg. 1 year)
766 C (note: externForcingCycle must be an integer
767 C number times externForcingPeriod)
768 C convertFW2Salt :: salinity, used to convert Fresh-Water Flux to Salt Flux
769 C (use model surface (local) value if set to -1)
770 C temp_EvPrRn :: temperature of Rain & Evap.
771 C salt_EvPrRn :: salinity of Rain & Evap.
772 C temp_addMass :: temperature of addMass array
773 C salt_addMass :: salinity of addMass array
774 C (notes: a) tracer content of Rain/Evap only used if both
775 C NonLin_FrSurf & useRealFreshWater are set.
776 C b) use model surface (local) value if set to UNSET_RL)
777 C hMixCriteria:: criteria for mixed-layer diagnostic
778 C dRhoSmall :: parameter for mixed-layer diagnostic
779 C hMixSmooth :: Smoothing parameter for mixed-layer diag (default=0=no smoothing)
780 C ivdc_kappa :: implicit vertical diffusivity for convection [m^2/s]
781 C sideDragFactor :: side-drag scaling factor (used only if no_slip_sides)
782 C (default=2: full drag ; =1: gives half-slip BC)
783 C bottomDragLinear :: Linear bottom-drag coefficient (units of [r]/s)
784 C bottomDragQuadratic :: Quadratic bottom-drag coefficient (units of [r]/m)
785 C (if using zcoordinate, units becomes linear: m/s, quadratic: [-])
786 C smoothAbsFuncRange :: 1/2 of interval around zero, for which FORTRAN ABS
787 C is to be replace by a smoother function
788 C (affects myabs, mymin, mymax)
789 C nh_Am2 :: scales the non-hydrostatic terms and changes internal scales
790 C (i.e. allows convection at different Rayleigh numbers)
791 C tCylIn :: Temperature of the cylinder inner boundary
792 C tCylOut :: Temperature of the cylinder outer boundary
793 C phiEuler :: Euler angle, rotation about original z-axis
794 C thetaEuler :: Euler angle, rotation about new x-axis
795 C psiEuler :: Euler angle, rotation about new z-axis
796 COMMON /PARM_R/ cg2dTargetResidual, cg2dTargetResWunit,
797 & cg2dpcOffDFac, cg3dTargetResidual,
798 & delR, delRc, xgOrigin, ygOrigin, rSphere, recip_rSphere,
799 & radius_fromHorizGrid, seaLev_Z, top_Pres, rSigmaBnd,
800 & deltaT, deltaTMom, dTtracerLev, deltaTFreeSurf, deltaTClock,
801 & abEps, alph_AB, beta_AB,
802 & f0, beta, fPrime, omega, rotationPeriod,
803 & viscFacAdj, viscAh, viscAhW, smag3D_coeff,
804 & viscAhMax, viscAhGrid, viscAhGridMax, viscAhGridMin,
805 & viscC2leith, viscC2leithD,
806 & viscC2smag, viscC4smag,
807 & viscAhD, viscAhZ, viscA4D, viscA4Z,
808 & viscA4, viscA4W, viscA4Max,
809 & viscA4Grid, viscA4GridMax, viscA4GridMin,
810 & viscAhReMax, viscA4ReMax,
811 & viscC4leith, viscC4leithD, viscArNr,
812 & diffKhT, diffK4T, diffKrNrT, diffKr4T,
813 & diffKhS, diffK4S, diffKrNrS, diffKr4S,
814 & diffKrBL79surf, diffKrBL79deep, diffKrBL79scl, diffKrBL79Ho,
815 & BL79LatVary,
816 & diffKrBLEQsurf, diffKrBLEQdeep, diffKrBLEQscl, diffKrBLEQHo,
817 & pCellMix_maxFac, pCellMix_delR, pCellMix_viscAr, pCellMix_diffKr,
818 & tauCD, rCD, epsAB_CD,
819 & freeSurfFac, implicSurfPress, implicDiv2DFlow, implicitNHPress,
820 & hFacMin, hFacMinDz, hFacInf, hFacSup,
821 & gravity, recip_gravity, gBaro,
822 & gravFacC, recip_gravFacC, gravFacF, recip_gravFacF,
823 & rhoNil, rhoConst, recip_rhoConst, rho1Ref,
824 & rhoFacC, recip_rhoFacC, rhoFacF, recip_rhoFacF, rhoConstFresh,
825 & thetaConst, tRef, sRef, pRef4EOS, phiRef, dBdrRef,
826 & rVel2wUnit, wUnit2rVel, mass2rUnit, rUnit2mass,
827 & baseTime, startTime, endTime,
828 & chkPtFreq, pChkPtFreq, dumpFreq, adjDumpFreq,
829 & diagFreq, taveFreq, tave_lastIter, monitorFreq, adjMonitorFreq,
830 & afFacMom, vfFacMom, pfFacMom, cfFacMom, foFacMom, mtFacMom,
831 & cosPower, cAdjFreq,
832 & tauThetaClimRelax, tauSaltClimRelax, latBandClimRelax,
833 & externForcingCycle, externForcingPeriod,
834 & convertFW2Salt, temp_EvPrRn, salt_EvPrRn,
835 & temp_addMass, salt_addMass, hFacMinDr, hFacMinDp,
836 & ivdc_kappa, hMixCriteria, dRhoSmall, hMixSmooth,
837 & sideDragFactor, bottomDragLinear, bottomDragQuadratic, nh_Am2,
838 & smoothAbsFuncRange,
839 & tCylIn, tCylOut,
840 & phiEuler, thetaEuler, psiEuler
841
842 _RL cg2dTargetResidual
843 _RL cg2dTargetResWunit
844 _RL cg3dTargetResidual
845 _RL cg2dpcOffDFac
846 _RL delR(Nr)
847 _RL delRc(Nr+1)
848 _RL xgOrigin
849 _RL ygOrigin
850 _RL rSphere
851 _RL recip_rSphere
852 _RL radius_fromHorizGrid
853 _RL seaLev_Z
854 _RL top_Pres
855 _RL rSigmaBnd
856 _RL deltaT
857 _RL deltaTClock
858 _RL deltaTMom
859 _RL dTtracerLev(Nr)
860 _RL deltaTFreeSurf
861 _RL abEps, alph_AB, beta_AB
862 _RL f0
863 _RL beta
864 _RL fPrime
865 _RL omega
866 _RL rotationPeriod
867 _RL freeSurfFac
868 _RL implicSurfPress
869 _RL implicDiv2DFlow
870 _RL implicitNHPress
871 _RL hFacMin
872 _RL hFacMinDz
873 _RL hFacMinDp
874 _RL hFacMinDr
875 _RL hFacInf
876 _RL hFacSup
877 _RL viscArNr(Nr)
878 _RL viscFacAdj
879 _RL viscAh
880 _RL viscAhW
881 _RL viscAhD
882 _RL viscAhZ
883 _RL smag3D_coeff
884 _RL viscAhMax
885 _RL viscAhReMax
886 _RL viscAhGrid, viscAhGridMax, viscAhGridMin
887 _RL viscC2leith
888 _RL viscC2leithD
889 _RL viscC2smag
890 _RL viscA4
891 _RL viscA4W
892 _RL viscA4D
893 _RL viscA4Z
894 _RL viscA4Max
895 _RL viscA4ReMax
896 _RL viscA4Grid, viscA4GridMax, viscA4GridMin
897 _RL viscC4leith
898 _RL viscC4leithD
899 _RL viscC4smag
900 _RL diffKhT
901 _RL diffK4T
902 _RL diffKrNrT(Nr)
903 _RL diffKr4T(Nr)
904 _RL diffKhS
905 _RL diffK4S
906 _RL diffKrNrS(Nr)
907 _RL diffKr4S(Nr)
908 _RL diffKrBL79surf
909 _RL diffKrBL79deep
910 _RL diffKrBL79scl
911 _RL diffKrBL79Ho
912 _RL BL79LatVary
913 _RL diffKrBLEQsurf
914 _RL diffKrBLEQdeep
915 _RL diffKrBLEQscl
916 _RL diffKrBLEQHo
917 _RL pCellMix_maxFac
918 _RL pCellMix_delR
919 _RL pCellMix_viscAr(Nr)
920 _RL pCellMix_diffKr(Nr)
921 _RL tauCD, rCD, epsAB_CD
922 _RL gravity, recip_gravity
923 _RL gBaro
924 _RL gravFacC(Nr), recip_gravFacC(Nr)
925 _RL gravFacF(Nr+1), recip_gravFacF(Nr+1)
926 _RL rhoNil
927 _RL rhoConst, recip_rhoConst
928 _RL rho1Ref(Nr)
929 _RL rhoFacC(Nr), recip_rhoFacC(Nr)
930 _RL rhoFacF(Nr+1), recip_rhoFacF(Nr+1)
931 _RL rhoConstFresh
932 _RL thetaConst
933 _RL tRef(Nr)
934 _RL sRef(Nr)
935 _RL pRef4EOS(Nr)
936 _RL phiRef(2*Nr+1)
937 _RL dBdrRef(Nr)
938 _RL rVel2wUnit(Nr+1), wUnit2rVel(Nr+1)
939 _RL mass2rUnit, rUnit2mass
940 _RL baseTime
941 _RL startTime
942 _RL endTime
943 _RL chkPtFreq
944 _RL pChkPtFreq
945 _RL dumpFreq
946 _RL adjDumpFreq
947 _RL diagFreq
948 _RL taveFreq
949 _RL tave_lastIter
950 _RL monitorFreq
951 _RL adjMonitorFreq
952 _RL afFacMom
953 _RL vfFacMom
954 _RL pfFacMom
955 _RL cfFacMom
956 _RL foFacMom
957 _RL mtFacMom
958 _RL cosPower
959 _RL cAdjFreq
960 _RL tauThetaClimRelax
961 _RL tauSaltClimRelax
962 _RL latBandClimRelax
963 _RL externForcingCycle
964 _RL externForcingPeriod
965 _RL convertFW2Salt
966 _RL temp_EvPrRn
967 _RL salt_EvPrRn
968 _RL temp_addMass
969 _RL salt_addMass
970 _RL ivdc_kappa
971 _RL hMixCriteria
972 _RL dRhoSmall
973 _RL hMixSmooth
974 _RL sideDragFactor
975 _RL bottomDragLinear
976 _RL bottomDragQuadratic
977 _RL smoothAbsFuncRange
978 _RL nh_Am2
979 _RL tCylIn, tCylOut
980 _RL phiEuler, thetaEuler, psiEuler
981
982 C-- COMMON /PARM_A/ Thermodynamics constants ?
983 COMMON /PARM_A/ HeatCapacity_Cp
984 _RL HeatCapacity_Cp
985
986 C-- COMMON /PARM_ATM/ Atmospheric physical parameters (Ideal Gas EOS, ...)
987 C celsius2K :: convert centigrade (Celsius) degree to Kelvin
988 C atm_Po :: standard reference pressure
989 C atm_Cp :: specific heat (Cp) of the (dry) air at constant pressure
990 C atm_Rd :: gas constant for dry air
991 C atm_kappa :: kappa = R/Cp (R: constant of Ideal Gas EOS)
992 C atm_Rq :: water vapour specific volume anomaly relative to dry air
993 C (e.g. typical value = (29/18 -1) 10^-3 with q [g/kg])
994 C integr_GeoPot :: option to select the way we integrate the geopotential
995 C (still a subject of discussions ...)
996 C selectFindRoSurf :: select the way surf. ref. pressure (=Ro_surf) is
997 C derived from the orography. Implemented: 0,1 (see INI_P_GROUND)
998 COMMON /PARM_ATM/
999 & celsius2K,
1000 & atm_Cp, atm_Rd, atm_kappa, atm_Rq, atm_Po,
1001 & integr_GeoPot, selectFindRoSurf
1002 _RL celsius2K
1003 _RL atm_Po, atm_Cp, atm_Rd, atm_kappa, atm_Rq
1004 INTEGER integr_GeoPot, selectFindRoSurf
1005
1006 C----------------------------------------
1007 C-- Logical flags for selecting packages
1008 LOGICAL useGAD
1009 LOGICAL useOBCS
1010 LOGICAL useSHAP_FILT
1011 LOGICAL useZONAL_FILT
1012 LOGICAL useOPPS
1013 LOGICAL usePP81
1014 LOGICAL useKL10
1015 LOGICAL useMY82
1016 LOGICAL useGGL90
1017 LOGICAL useKPP
1018 LOGICAL useGMRedi
1019 LOGICAL useDOWN_SLOPE
1020 LOGICAL useBBL
1021 LOGICAL useCAL
1022 LOGICAL useEXF
1023 LOGICAL useBulkForce
1024 LOGICAL useEBM
1025 LOGICAL useCheapAML
1026 LOGICAL useAUTODIFF
1027 LOGICAL useGrdchk
1028 LOGICAL useSMOOTH
1029 LOGICAL usePROFILES
1030 LOGICAL useECCO
1031 LOGICAL useCTRL
1032 LOGICAL useSBO
1033 LOGICAL useFLT
1034 LOGICAL usePTRACERS
1035 LOGICAL useGCHEM
1036 LOGICAL useRBCS
1037 LOGICAL useOffLine
1038 LOGICAL useMATRIX
1039 LOGICAL useFRAZIL
1040 LOGICAL useSEAICE
1041 LOGICAL useSALT_PLUME
1042 LOGICAL useShelfIce
1043 LOGICAL useStreamIce
1044 LOGICAL useICEFRONT
1045 LOGICAL useThSIce
1046 LOGICAL useLand
1047 LOGICAL useATM2d
1048 LOGICAL useAIM
1049 LOGICAL useAtm_Phys
1050 LOGICAL useFizhi
1051 LOGICAL useGridAlt
1052 LOGICAL useDiagnostics
1053 LOGICAL useREGRID
1054 LOGICAL useLayers
1055 LOGICAL useMNC
1056 LOGICAL useRunClock
1057 LOGICAL useEMBED_FILES
1058 LOGICAL useMYPACKAGE
1059 COMMON /PARM_PACKAGES/
1060 & useGAD, useOBCS, useSHAP_FILT, useZONAL_FILT,
1061 & useOPPS, usePP81, useKL10, useMY82, useGGL90, useKPP,
1062 & useGMRedi, useBBL, useDOWN_SLOPE,
1063 & useCAL, useEXF, useBulkForce, useEBM, useCheapAML,
1064 & useGrdchk, useSMOOTH, usePROFILES, useECCO, useCTRL,
1065 & useSBO, useFLT, useAUTODIFF,
1066 & usePTRACERS, useGCHEM, useRBCS, useOffLine, useMATRIX,
1067 & useFRAZIL, useSEAICE, useSALT_PLUME, useShelfIce,
1068 & useStreamIce, useICEFRONT, useThSIce, useLand,
1069 & useATM2D, useAIM, useAtm_Phys, useFizhi, useGridAlt,
1070 & useDiagnostics, useREGRID, useLayers, useMNC,
1071 & useRunClock, useEMBED_FILES,
1072 & useMYPACKAGE
1073
1074 CEH3 ;;; Local Variables: ***
1075 CEH3 ;;; mode:fortran ***
1076 CEH3 ;;; End: ***

  ViewVC Help
Powered by ViewVC 1.1.22