/[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.46 - (show annotations) (download)
Fri Feb 2 21:04:47 2001 UTC (23 years, 3 months ago) by adcroft
Branch: MAIN
Changes since 1.45: +23 -36 lines
File MIME type: text/plain
Merged changes from branch "branch-atmos-merge" into MAIN (checkpoint34)
 - substantial modifications to algorithm sequence (dynamics.F)
 - packaged OBCS, Shapiro filter, Zonal filter, Atmospheric Physics

1 C $Header: /u/gcmpack/models/MITgcmUV/model/inc/PARAMS.h,v 1.45.2.6 2001/01/30 21:02:58 adcroft Exp $
2 C
3 C /==========================================================\
4 C | PARAMS.h |
5 C | o Header file defining model "parameters". |
6 C |==========================================================|
7 C | The values from the model standard input file are |
8 C | stored into the variables held here. Notes describing |
9 C | the parameters can also be found here. |
10 C \==========================================================/
11
12 C Macros for special grid options
13 #include "PARAMS_MACROS.h"
14
15 C-- Contants
16 C Useful physical values
17 Real*8 PI
18 PARAMETER ( PI = 3.14159265358979323844D0 )
19 Real*8 deg2rad
20 PARAMETER ( deg2rad = 2.D0*PI/360.D0 )
21
22 C Symbolic values
23 C precXXXX - Used to indicate what precision to use for
24 C dumping model state.
25 INTEGER precFloat32
26 PARAMETER ( precFloat32 = 32 )
27 INTEGER precFloat64
28 PARAMETER ( precFloat64 = 64 )
29 C UNSET_xxx - Used to indicate variables that have not been given a value
30 Real*8 UNSET_FLOAT8
31 PARAMETER ( UNSET_FLOAT8 = 1.234567D5 )
32 Real*4 UNSET_FLOAT4
33 PARAMETER ( UNSET_FLOAT4 = 1.234567E5 )
34 _RL UNSET_RL
35 PARAMETER ( UNSET_RL = 1.234567D5 )
36 _RS UNSET_RS
37 PARAMETER ( UNSET_RS = 1.234567E5 )
38 INTEGER UNSET_I
39 PARAMETER ( UNSET_I = 123456789 )
40
41 C Checkpoint data
42 INTEGER maxNoChkptLev
43 PARAMETER ( maxNoChkptLev = 2 )
44
45 C-- COMMON /PARM_C/ Character valued parameters used by the model.
46 C checkPtSuff - List of checkpoint file suffices
47 C bathyFile - File containing bathymetry. If not defined bathymetry
48 C is taken from inline function.
49 C hydrogThetaFile - File containing initial hydrographic data for potential
50 C temperature.
51 C hydrogSaltFile - File containing initial hydrographic data for salinity.
52 C zonalWindFile - File containing zonal wind data
53 C meridWindFile - File containing meridional wind data
54 C thetaClimFile - File containing theta climataology used
55 C in relaxation term -lambda(theta-theta*)
56 C saltClimFile - File containing salt climataology used
57 C in relaxation term -lambda(salt-salt*)
58 C surfQfile - File containing surface heat flux
59 C surfQswfile - File containing surface shortwave radiation
60 C EmPmRfile - File containing surface fresh water flux
61 C buoyancyRelation - Flag used to indicate which relation to use to
62 C get buoyancy.
63 COMMON /PARM_C/ checkPtSuff,
64 & bathyFile, hydrogThetaFile, hydrogSaltFile,
65 & zonalWindFile, meridWindFile, thetaClimFile,
66 & saltClimFile, buoyancyRelation,
67 & EmPmRfile, surfQfile, surfQswfile,
68 & uVelInitFile, vVelInitFile, pSurfInitFile
69 CHARACTER*(5) checkPtSuff(maxNoChkptLev)
70 CHARACTER*(MAX_LEN_FNAM) bathyFile
71 CHARACTER*(MAX_LEN_FNAM) hydrogThetaFile
72 CHARACTER*(MAX_LEN_FNAM) hydrogSaltFile
73 CHARACTER*(MAX_LEN_FNAM) zonalWindFile
74 CHARACTER*(MAX_LEN_FNAM) meridWindFile
75 CHARACTER*(MAX_LEN_FNAM) thetaClimFile
76 CHARACTER*(MAX_LEN_FNAM) saltClimFile
77 CHARACTER*(MAX_LEN_FNAM) surfQfile
78 CHARACTER*(MAX_LEN_FNAM) surfQswfile
79 CHARACTER*(MAX_LEN_FNAM) EmPmRfile
80 CHARACTER*(MAX_LEN_FNAM) buoyancyRelation
81 CHARACTER*(MAX_LEN_FNAM) uVelInitFile
82 CHARACTER*(MAX_LEN_FNAM) vVelInitFile
83 CHARACTER*(MAX_LEN_FNAM) pSurfInitFile
84
85 C-- COMMON /PARM_I/ Integer valued parameters used by the model.
86 C cg2dMaxIters - Maximum number of iterations in the
87 C two-dimensional con. grad solver.
88 C cg2dChkResFreq - Frequency with which to check residual
89 C in con. grad solver.
90 C cg3dMaxIters - Maximum number of iterations in the
91 C three-dimensional con. grad solver.
92 C cg3dChkResFreq - Frequency with which to check residual
93 C in con. grad solver.
94 C nIter0 - Start time-step number of for this run
95 C nTimeSteps - Number of timesteps to execute
96 C numStepsPerPickup - For offline setup. Frequency of pickup
97 C of flow fields.
98 C writeStatePrec - Precision used for writing model state.
99 C writeBinaryPrec - Precision used for writing binary files
100 C readBinaryPrec - Precision used for reading binary files
101 C nCheckLev - Holds current checkpoint level
102 C nShap - "Order" of shapiro filter to apply to
103 C model prognositic fields.
104 C nShap == 1 => del2
105 C nShap == 2 => del4
106 C nShap == 3 => del6
107 C etc...
108
109 COMMON /PARM_I/
110 & cg2dMaxIters,
111 & cg2dChkResFreq,
112 & cg3dMaxIters,
113 & cg3dChkResFreq,
114 & nIter0, nTimeSteps, nEndIter,
115 & numStepsPerPickup,
116 & writeStatePrec, nCheckLev,
117 & writeBinaryPrec, readBinaryPrec,
118 & nShap, zonal_filt_sinpow, zonal_filt_cospow
119 INTEGER cg2dMaxIters
120 INTEGER cg2dChkResFreq
121 INTEGER cg3dMaxIters
122 INTEGER cg3dChkResFreq
123 INTEGER nIter0
124 INTEGER nTimeSteps
125 INTEGER nEndIter
126 INTEGER numStepsPerPickup
127 INTEGER writeStatePrec
128 INTEGER writeBinaryPrec
129 INTEGER readBinaryPrec
130 INTEGER nCheckLev
131 INTEGER nShap
132 INTEGER zonal_filt_sinpow
133 INTEGER zonal_filt_cospow
134
135 C-- COMMON /PARM_L/ Logical valued parameters used by the model.
136 C usingCartesianGrid - If TRUE grid generation will be in a cartesian
137 C coordinate frame.
138 C usingSphericalPolarGrid - If TRUE grid generation will be in a
139 C spherical polar frame.
140 C no_slip_sides - Impose "no-slip" at lateral boundaries.
141 C no_slip_bottom- Impose "no-slip" at bottom boundary.
142 C staggerTimeStep - enable a Stagger time stepping T,S Rho then U,V
143 C momViscosity - Flag which turns momentum friction terms on and off.
144 C momAdvection - Flag which turns advection of momentum on and off.
145 C momForcing - Flag which turns external forcing of momentum on
146 C and off.
147 C momPressureForcing - Flag which turns pressure term in momentum equation
148 C on and off.
149 C metricTerms - Flag which turns metric terms on or off.
150 C usingSphericalPolarMTerms - If TRUE use spherical polar metric terms.
151 C useCoriolis - Flag which turns the coriolis terms on and off.
152 C tempDiffusion - Flag which turns diffusion of temperature on
153 C and off.
154 C tempAdvection - Flag which turns advection of temperature on
155 C and off.
156 C tempForcing - Flag which turns external forcing of temperature on
157 C and off.
158 C saltDiffusion - Flag which turns diffusion of salinit on
159 C and off.
160 C saltAdvection - Flag which turns advection of salinit on
161 C and off.
162 C saltForcing - Flag which turns external forcing of salinit on
163 C and off.
164 C implicitFreeSurface - Set to true to use implcit free surface
165 C rigidLid - Set to true to use rigid lid
166 C momStepping - Turns momentum equation time-stepping off
167 C tempStepping - Turns temperature equation time-stepping off
168 C saltStepping - Turns salinity equation time-stepping off
169 C useConstantF - Coriolis parameter set to f0
170 C useBetaPlaneF - Coriolis parameter set to f0 + beta.y
171 C useSphereF - Coriolis parameter set to 2.omega.sin(phi)
172 C implicitDiffusion - Turns implicit vertical diffusion on
173 C implicitViscosity - Turns implicit vertical viscosity on
174 C doThetaClimRelax - Set true if relaxation to temperature
175 C climatology is required.
176 C doSaltClimRelax - Set true if relaxation to salinity
177 C climatology is required.
178 C periodicExternalForcing - Set true if forcing is time-dependant
179 C usingPCoords - Set to indicate that we are working in pressure
180 C coords.
181 C usingZCoords - Set to indicate that we are working in height
182 C coords.
183 C nonHydrostatic - Using non-hydrostatic terms
184 C globalFiles - Selects between "global" and "tiled" files
185 C allowFreezing - Allows water to freeze and form ice
186 C groundAtK1 - put the surface(k=1) at the Lower Boundary (=ground)
187 COMMON /PARM_L/ usingCartesianGrid, usingSphericalPolarGrid,
188 & no_slip_sides,no_slip_bottom,
189 & staggerTimeStep,
190 & momViscosity, momAdvection, momForcing, useCoriolis,
191 & momPressureForcing,tempDiffusion, tempAdvection, tempForcing,
192 & saltDiffusion, saltAdvection, saltForcing,
193 & implicitFreeSurface, rigidLid,
194 & momStepping, tempStepping, saltStepping,
195 & metricTerms, usingSphericalPolarMTerms,
196 & useConstantF, useBetaPlaneF, useSphereF,
197 & implicitDiffusion, implicitViscosity,
198 & doThetaClimRelax, doSaltClimRelax,
199 & periodicExternalForcing, usingPCoords, usingZCoords,
200 & nonHydrostatic, globalFiles,
201 & allowFreezing, groundAtK1
202 LOGICAL usingCartesianGrid
203 LOGICAL usingSphericalPolarGrid
204 LOGICAL usingSphericalPolarMTerms
205 LOGICAL no_slip_sides
206 LOGICAL no_slip_bottom
207 LOGICAL staggerTimeStep
208 LOGICAL momViscosity
209 LOGICAL momAdvection
210 LOGICAL momForcing
211 LOGICAL momPressureForcing
212 LOGICAL useCoriolis
213 LOGICAL tempDiffusion
214 LOGICAL tempAdvection
215 LOGICAL tempForcing
216 LOGICAL saltDiffusion
217 LOGICAL saltAdvection
218 LOGICAL saltForcing
219 LOGICAL implicitFreeSurface
220 LOGICAL rigidLid
221 LOGICAL momStepping
222 LOGICAL tempStepping
223 LOGICAL saltStepping
224 LOGICAL metricTerms
225 LOGICAL useConstantF
226 LOGICAL useBetaPlaneF
227 LOGICAL useSphereF
228 LOGICAL implicitDiffusion
229 LOGICAL implicitViscosity
230 LOGICAL doThetaClimRelax
231 LOGICAL doSaltClimRelax
232 LOGICAL periodicExternalForcing
233 LOGICAL usingPCoords
234 LOGICAL usingZCoords
235 LOGICAL nonHydrostatic
236 LOGICAL globalFiles
237 LOGICAL allowFreezing
238 LOGICAL groundAtK1
239
240 C-- COMMON /PARM_R/ "Real" valued parameters used by the model.
241 C cg2dTargetResidual
242 C - Target residual for cg2d solver.
243 C cg3dTargetResidual
244 C - Target residual for cg3d solver.
245 C cg2dpcOffDFac - Averaging weight for preconditioner off-diagonal.
246 C Note. 20th May 1998
247 C I made a weird discovery! In the model paper we argue
248 C for the form of the preconditioner used here ( see
249 C A Finite-volume, Incompressible Navier-Stokes Model
250 C ...., Marshall et. al ). The algebra gives a simple
251 C 0.5 factor for the averaging of ac and aCw to get a
252 C symmettric pre-conditioner. By using a factor of 0.51
253 C i.e. scaling the off-diagonal terms in the
254 C preconditioner down slightly I managed to get the
255 C number of iterations for convergence in a test case to
256 C drop form 192 -> 134! Need to investigate this further!
257 C For now I have introduced a parameter cg2dpcOffDFac which
258 C defaults to 0.51 but can be set at runtime.
259 C delP - Vertical grid spacing ( Pa ).
260 C delZ - Vertical grid spacing ( m ).
261 C delR - Vertical grid spacing ( units of r ).
262 C delX - Separation between cell faces (m) or (deg), depending
263 C delY on input flags.
264 C gravity - Accel. due to gravity ( m/s^2 )
265 C recip_gravity and its inverse
266 C gBaro - Accel. due to gravity used in barotropic equation ( m/s^2 )
267 C ronil - Reference density
268 C rhoConst - Vertically constant reference density
269 C startTime - Start time for model ( s )
270 C phiMin - Latitude of southern most cell face.
271 C thetaMin - Longitude of western most cell face (this
272 C is an "inert" parameter but it is included
273 C to make geographical references simple.)
274 C rSphere - Radius of sphere for a spherical polar grid ( m ).
275 C recip_RSphere - Reciprocal radius of sphere ( m ).
276 C f0 - Reference coriolis parameter ( 1/s )
277 C ( Southern edge f for beta plane )
278 C beta - df/dy ( s^-1.m^-1 )
279 C omega - Angular velocity ( rad/s )
280 C viscAh - Eddy viscosity coeff. for mixing of
281 C momentum laterally ( m^2/s )
282 C viscAz - Eddy viscosity coeff. for mixing of
283 C momentum vertically ( m^2/s )
284 C viscAp - Eddy viscosity coeff. for mixing of
285 C momentum vertically ( Pa^2/s )
286 C viscAr - Eddy viscosity coeff. for mixing of
287 C momentum vertically ( units of r^2/s )
288 C viscA4 - Biharmonic viscosity coeff. for mixing of
289 C momentum laterally ( m^4/s )
290 C diffKhT - Laplacian diffusion coeff. for mixing of
291 C heat laterally ( m^2/s )
292 C diffKzT - Laplacian diffusion coeff. for mixing of
293 C heat vertically ( m^2/s )
294 C diffKpT - Laplacian diffusion coeff. for mixing of
295 C heat vertically ( Pa^2/s )
296 C diffKrT - Laplacian diffusion coeff. for mixing of
297 C heat vertically ( units of r^2/s )
298 C diffK4T - Biharmonic diffusion coeff. for mixing of
299 C heat laterally ( m^4/s )
300 C diffKhS - Laplacian diffusion coeff. for mixing of
301 C salt laterally ( m^2/s )
302 C diffKzS - Laplacian diffusion coeff. for mixing of
303 C salt vertically ( m^2/s )
304 C diffKpS - Laplacian diffusion coeff. for mixing of
305 C salt vertically ( Pa^2/s )
306 C diffKrS - Laplacian diffusion coeff. for mixing of
307 C salt vertically ( units of r^2/s )
308 C diffK4S - Biharmonic diffusion coeff. for mixing of
309 C salt laterally ( m^4/s )
310 C deltaT - Default timestep ( s )
311 C deltaTClock - Timestep used as model "clock". This determines the
312 C IO frequencies and is used in tagging output. It can
313 C be totally different to the dynamical time. Typically
314 C it will be the deep-water timestep for accelerated runs.
315 C Frequency of checkpointing and dumping of the model state
316 C are referenced to this clock. ( s )
317 C deltaTMom - Timestep for momemtum equations ( s )
318 C deltaTtracer - Timestep for tracer equations ( s )
319 C freesurfFac - Parameter to turn implicit free surface term on or off
320 C freesurfac = 1. uses implicit free surface
321 C freesurfac = 0. uses rigid lid
322 C hFacMin - Minimum fraction size of a cell (affects hFacC etc...)
323 C hFacMinDz - Minimum dimesional size of a cell (affects hFacC etc..., m)
324 C hFacMinDp - Minimum dimesional size of a cell (affects hFacC etc..., Pa)
325 C hFacMinDr - Minimum dimesional size of a cell (affects hFacC etc..., units of r)
326 C tauCD - CD scheme coupling timescale ( 1/s )
327 C rCD - CD scheme normalised coupling parameter ( 0-1 )
328 C startTime - Starting time for this integration ( s ).
329 C endTime - Ending time for this integration ( s ).
330 C chkPtFreq - Frequency of rolling check pointing ( s ).
331 C pChkPtFreq - Frequency of permanent check pointing ( s ).
332 C dumpFreq - Frequency with which model state is written to
333 C post-processing files ( s ).
334 C afFacMom - Advection of momentum term tracer parameter
335 C vfFacMom - Momentum viscosity tracer parameter
336 C pfFacMom - Momentum pressure forcing tracer parameter
337 C cfFacMom - Coriolis term tracer parameter
338 C foFacMom - Momentum forcing tracer parameter
339 C mtFacMom - Metric terms tracer parameter
340 C cosPower - Power of cosine of latitude to multiply viscosity
341 C cAdjFreq - Frequency of convective adjustment
342 C
343 C taveFreq - Frequency with which time-averaged model state is written to
344 C post-processing files ( s ).
345 C tauThetaClimRelax - Relaxation to climatology time scale ( s ).
346 C lambdaThetaClimRelax - Inverse time scale for relaxation ( 1/s ).
347 C tauSaltClimRelax - Relaxation to climatology time scale ( s ).
348 C lambdaSaltClimRelax - Inverse time scale for relaxation ( 1/s ).
349 C externForcingPeriod - Is the period of which forcing varies (eg. 1 month)
350 C externForcingCycle - Is the repeat time of the forcing (eg. 1 year)
351 C (note: externForcingCycle must be an integer
352 C number times externForcingPeriod)
353 C horiVertRatio - Ratio on units in vertical to units in horizontal.
354 C recip_horiVertRatio ( 1 if horiz in m and vertical in m ).
355 C ( g*rho if horiz in m and vertical in Pa ).
356 C latFFTFiltLo - Low latitude for FFT filtering of latitude
357 C circles ( see filter*.F )
358 C Ro_SeaLevel - standard position of Sea-Level in "R" coordinate, used as
359 C starting value (k=1) for vertical coordinate (rf(1)=Ro_SeaLevel)
360 C bottomDragLinear - Drag coefficient built in to core dynamics
361 C " Quadratic ( linear: 1/s, quadratic: 1/m )
362 COMMON /PARM_R/ cg2dTargetResidual, cg2dpcOffDFac,
363 & cg3dTargetResidual,
364 & delP, delZ, delR, delX, delY,
365 & deltaT,deltaTmom, deltaTtracer, deltaTClock,abeps, startTime,
366 & phiMin, thetaMin, rSphere, recip_RSphere, f0, fCori, beta,
367 & viscAh, viscAz, viscA4, viscAr,
368 & diffKhT, diffKzT, diffK4T, diffKrT,
369 & diffKhS, diffKzS, diffK4S, diffKrS,
370 & delT, tauCD, rCD, freeSurfFac, hFacMin, hFacMinDz,
371 & gravity, recip_Gravity, gBaro, rhonil, recip_rhonil,
372 & recip_rhoConst, rhoConst, tRef, sRef,
373 & endTime, chkPtFreq, pchkPtFreq, dumpFreq, taveFreq,
374 & afFacMom, vfFacMom, pfFacMom, cfFacMom, foFacMom, mtFacMom,
375 & cosPower,
376 & cAdjFreq, omega, tauThetaClimRelax, lambdaThetaClimRelax,
377 & tauSaltClimRelax, lambdaSaltClimRelax,
378 & externForcingCycle, externForcingPeriod,
379 & viscAp, diffKpT, diffKpS, hFacMinDr, hFacMinDp,
380 & theta_S, specVol_S, horiVertRatio, recip_horiVertRatio,
381 & latFFTFiltLo, ivdc_kappa, Ro_SeaLevel, zonal_filt_lat,
382 & bottomDragLinear,bottomDragQuadratic
383
384 _RL cg2dTargetResidual
385 _RL cg3dTargetResidual
386 _RL cg2dpcOffDFac
387 _RL delZ(Nr)
388 _RL delP(Nr)
389 _RL delR(Nr)
390 _RL delX(Nx)
391 _RL delY(Ny)
392 _RL deltaT
393 _RL deltaTClock
394 _RL deltaTmom
395 _RL deltaTtracer
396 _RL abeps
397 _RL phiMin
398 _RL thetaMin
399 _RL rSphere
400 _RL recip_RSphere
401 _RL f0
402 _RL freeSurfFac
403 _RL hFacMin
404 _RL hFacMinDz
405 _RL hFacMinDp
406 _RL hFacMinDr
407 _RL beta
408 _RL viscAh
409 _RL viscAz
410 _RL viscAp
411 _RL viscAr
412 _RL viscA4
413 _RL diffKhT
414 _RL diffKrT
415 _RL diffKzT
416 _RL diffKpT
417 _RL diffK4T
418 _RL diffKhS
419 _RL diffKrS
420 _RL diffKzS
421 _RL diffKpS
422 _RL diffK4S
423 _RL delt
424 _RL tauCD
425 _RL rCD
426 _RL gravity
427 _RL recip_gravity
428 _RL gBaro
429 _RL rhonil
430 _RL recip_rhonil
431 _RL rhoConst
432 _RL recip_rhoConst
433 _RL specVol_S(Nr)
434 _RL tRef(Nr)
435 _RL theta_S(Nr)
436 _RL sRef(Nr)
437 _RS Fcori(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
438 _RL startTime
439 _RL endTime
440 _RL chkPtFreq
441 _RL pChkPtFreq
442 _RL dumpFreq
443 _RL taveFreq
444 _RL afFacMom
445 _RL vfFacMom
446 _RL pfFacMom
447 _RL cfFacMom
448 _RL foFacMom
449 _RL mTFacMom
450 _RL cosPower
451 _RL cAdjFreq
452 _RL omega
453 _RL tauThetaClimRelax
454 _RL lambdaThetaClimRelax
455 _RL tauSaltClimRelax
456 _RL lambdaSaltClimRelax
457 _RL externForcingCycle
458 _RL externForcingPeriod
459 _RL horiVertRatio
460 _RL recip_horiVertRatio
461 _RL latFFTFiltLo
462 _RL ivdc_kappa
463 _RL Ro_SeaLevel
464 _RL zonal_filt_lat
465 _RL bottomDragLinear
466 _RL bottomDragQuadratic
467
468 COMMON /PARM_A/ HeatCapacity_Cp,recip_Cp,
469 & Lamba_theta
470 _RL HeatCapacity_Cp
471 _RL Lamba_theta
472 _RL recip_Cp
473
474 C Equation of State (polynomial coeffients)
475 COMMON /PARM_EOS_NL/ eosC,eosSig0,eosRefT,eosRefS
476 _RL eosC(9,Nr+1),eosSig0(Nr+1),eosRefT(Nr+1),eosRefS(Nr+1)
477 C Linear equation of state
478 C tAlpha - Linear EOS thermal expansion coefficient ( 1/degree ).
479 C sBeta - Linear EOS haline contraction coefficient.
480 COMMON /PARM_EOS_LIN/ tAlpha,sBeta,eosType
481 _RL tAlpha
482 _RL sBeta
483 character*(6) eosType
484
485 C Logical flags for selecting packages
486 LOGICAL useKPP
487 LOGICAL useGMRedi
488 LOGICAL useOBCS
489 LOGICAL useECCO
490 COMMON /PARM_PACKAGES/
491 & useKPP, useGMRedi, useOBCS, useECCO

  ViewVC Help
Powered by ViewVC 1.1.22