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

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

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


Revision 1.33 - (hide annotations) (download)
Wed Oct 28 03:11:36 1998 UTC (25 years, 6 months ago) by cnh
Branch: MAIN
CVS Tags: checkpoint17, checkpoint16
Changes since 1.32: +18 -6 lines
File MIME type: text/plain
Changes to support
 - g77 compilation under Linux
 - LR(1) form of 64-bit is D or E for constants
 - Modified adjoint of exch with adjoint variables
   acuumulated.

1 cnh 1.33 C $Header: /u/gcmpack/models/MITgcmUV/model/inc/PARAMS.h,v 1.32 1998/09/09 15:19:07 cnh Exp $
2 cnh 1.1 C
3     C /==========================================================\
4     C | PARAMS.h |
5     C | o Header file defining model "parameters". |
6     C |==========================================================|
7 cnh 1.33 C | The values from the model standard input file are |
8 cnh 1.1 C | stored into the variables held here. Notes describing |
9     C | the parameters can also be found here. |
10     C \==========================================================/
11 cnh 1.13
12     C Macros for special grid options
13     #include "PARAMS_MACROS.h"
14 cnh 1.1
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 cnh 1.7 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 = 0 )
27     INTEGER precFloat64
28     PARAMETER ( precFloat64 = 1 )
29 cnh 1.28 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 cnh 1.7
41     C Checkpoint data
42     INTEGER maxNoChkptLev
43     PARAMETER ( maxNoChkptLev = 2 )
44    
45 cnh 1.1 C-- COMMON /PARM_C/ Character valued parameters used by the model.
46 cnh 1.7 C checkPtSuff - List of checkpoint file suffices
47 cnh 1.16 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 cnh 1.18 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 cnh 1.26 C buoyancyRelation - Flag used to indicate which relation to use to
59     C get buoyancy.
60 cnh 1.16 COMMON /PARM_C/ checkPtSuff,
61     & bathyFile, hydrogThetaFile, hydrogSaltFile,
62 cnh 1.18 & zonalWindFile, meridWindFile, thetaClimFile,
63 cnh 1.26 & saltClimFile, buoyancyRelation
64 cnh 1.7 CHARACTER*(5) checkPtSuff(maxNoChkptLev)
65 cnh 1.16 CHARACTER*(MAX_LEN_FNAM) bathyFile
66     CHARACTER*(MAX_LEN_FNAM) hydrogThetaFile
67     CHARACTER*(MAX_LEN_FNAM) hydrogSaltFile
68     CHARACTER*(MAX_LEN_FNAM) zonalWindFile
69     CHARACTER*(MAX_LEN_FNAM) meridWindFile
70 cnh 1.18 CHARACTER*(MAX_LEN_FNAM) thetaClimFile
71     CHARACTER*(MAX_LEN_FNAM) saltClimFile
72 cnh 1.26 CHARACTER*(MAX_LEN_FNAM) buoyancyRelation
73 cnh 1.7
74 cnh 1.1 C-- COMMON /PARM_I/ Integer valued parameters used by the model.
75     C cg2dMaxIters - Maximum number of iterations in the
76     C two-dimensional con. grad solver.
77     C cg2dChkResFreq - Frequency with which to check residual
78     C in con. grad solver.
79 cnh 1.33 C cg3dMaxIters - Maximum number of iterations in the
80     C three-dimensional con. grad solver.
81     C cg3dChkResFreq - Frequency with which to check residual
82     C in con. grad solver.
83 cnh 1.1 C nIter0 - Start time-step number of for this run
84     C nTimeSteps - Number of timesteps to execute
85     C numStepsPerPickup - For offline setup. Frequency of pickup
86     C of flow fields.
87 cnh 1.7 C writeStatePrec - Precision used for writing model state.
88     C writeBinaryPrec - Precision used for writing binary files
89     C readBinaryPrec - Precision used for reading binary files
90     C nCheckLev - Holds current checkpoint level
91 cnh 1.1 COMMON /PARM_I/
92     & cg2dMaxIters,
93     & cg2dChkResFreq,
94 cnh 1.33 & cg3dMaxIters,
95     & cg3dChkResFreq,
96 cnh 1.1 & nIter0, nTimeSteps,
97 cnh 1.7 & numStepsPerPickup,
98     & writeStatePrec, nCheckLev,
99     & writeBinaryPrec, readBinaryPrec
100 cnh 1.1 INTEGER cg2dMaxIters
101     INTEGER cg2dChkResFreq
102 cnh 1.33 INTEGER cg3dMaxIters
103     INTEGER cg3dChkResFreq
104 cnh 1.1 INTEGER nIter0
105     INTEGER nTimeSteps
106     INTEGER numStepsPerPickup
107 cnh 1.7 INTEGER writeStatePrec
108     INTEGER writeBinaryPrec
109     INTEGER readBinaryPrec
110     INTEGER nCheckLev
111 cnh 1.1
112     C-- COMMON /PARM_L/ Logical valued parameters used by the model.
113     C usingCartesianGrid - If TRUE grid generation will be in a cartesian
114     C coordinate frame.
115     C usingSphericalPolarGrid - If TRUE grid generation will be in a
116     C spherical polar frame.
117     C momViscosity - Flag which turns momentum friction terms on and off.
118     C momAdvection - Flag which turns advection of momentum on and off.
119     C momForcing - Flag which turns external forcing of momentum on
120     C and off.
121 cnh 1.9 C momPressureForcing - Flag which turns pressure term in momentum equation
122     C on and off.
123 cnh 1.15 C metricTerms - Flag which turns metric terms on or off.
124     C usingSphericalPolarMTerms - If TRUE use spherical polar metric terms.
125 cnh 1.1 C useCoriolis - Flag which turns the coriolis terms on and off.
126     C tempDiffusion - Flag which turns diffusion of temperature on
127     C and off.
128     C tempAdvection - Flag which turns advection of temperature on
129     C and off.
130     C tempForcing - Flag which turns external forcing of temperature on
131     C and off.
132     C saltDiffusion - Flag which turns diffusion of salinit on
133     C and off.
134     C saltAdvection - Flag which turns advection of salinit on
135     C and off.
136     C saltForcing - Flag which turns external forcing of salinit on
137     C and off.
138 cnh 1.8 C implicitFreeSurface - Set to true to use implcit free surface
139     C rigidLid - Set to true to use rigid lid
140 cnh 1.10 C momStepping - Turns momentum equation time-stepping off
141 cnh 1.15 C tempStepping - Turns temperature equation time-stepping off
142 adcroft 1.17 C saltStepping - Turns salinity equation time-stepping off
143 cnh 1.15 C useConstantF - Coriolis parameter set to f0
144     C useBetaPlaneF - Coriolis parameter set to f0 + beta.y
145     C useSphereF - Coriolis parameter set to 2.omega.sin(phi)
146 adcroft 1.14 C implicitDiffusion - Turns implicit vertical diffusion on
147 cnh 1.18 C doThetaClimRelax - Set true if relaxation to temperature
148     C climatology is required.
149     C doSaltClimRelax - Set true if relaxation to salinity
150     C climatology is required.
151 adcroft 1.19 C periodicExternalForcing - Set true if forcing is time-dependant
152 cnh 1.29 C usingPCoords - Set to indicate that we are working in pressure
153     C coords.
154     C usingZCoords - Set to indicate that we are working in height
155     C coords.
156 cnh 1.1 COMMON /PARM_L/ usingCartesianGrid, usingSphericalPolarGrid,
157 cnh 1.33 & momViscosity, momAdvection, momForcing, useCoriolis,
158     & momPressureForcing,tempDiffusion, tempAdvection, tempForcing,
159 cnh 1.8 & saltDiffusion, saltAdvection, saltForcing,
160 adcroft 1.17 & implicitFreeSurface, rigidLid,
161     & momStepping, tempStepping, saltStepping,
162 cnh 1.15 & metricTerms, usingSphericalPolarMTerms,
163     & useConstantF, useBetaPlaneF, useSphereF,
164 adcroft 1.19 & implicitDiffusion, doThetaClimRelax, doSaltClimRelax,
165 cnh 1.29 & periodicExternalForcing, usingPCoords, usingZCoords
166 cnh 1.1 LOGICAL usingCartesianGrid
167     LOGICAL usingSphericalPolarGrid
168 cnh 1.15 LOGICAL usingSphericalPolarMTerms
169 cnh 1.1 LOGICAL momViscosity
170     LOGICAL momAdvection
171     LOGICAL momForcing
172 cnh 1.9 LOGICAL momPressureForcing
173 cnh 1.1 LOGICAL useCoriolis
174     LOGICAL tempDiffusion
175     LOGICAL tempAdvection
176     LOGICAL tempForcing
177     LOGICAL saltDiffusion
178     LOGICAL saltAdvection
179     LOGICAL saltForcing
180 cnh 1.8 LOGICAL implicitFreeSurface
181     LOGICAL rigidLid
182 cnh 1.10 LOGICAL momStepping
183     LOGICAL tempStepping
184 adcroft 1.17 LOGICAL saltStepping
185 cnh 1.15 LOGICAL metricTerms
186     LOGICAL useConstantF
187     LOGICAL useBetaPlaneF
188     LOGICAL useSphereF
189 adcroft 1.14 LOGICAL implicitDiffusion
190 cnh 1.18 LOGICAL doThetaClimRelax
191     LOGICAL doSaltClimRelax
192 adcroft 1.19 LOGICAL periodicExternalForcing
193 cnh 1.29 LOGICAL usingPCoords
194     LOGICAL usingZCoords
195 cnh 1.1
196     C-- COMMON /PARM_R/ "Real" valued parameters used by the model.
197     C cg2dTargetResidual
198     C - Target residual for cg2d solver.
199 cnh 1.33 C cg3dTargetResidual
200     C - Target residual for cg3d solver.
201 cnh 1.7 C cg2dpcOffDFac - Averaging weight for preconditioner off-diagonal.
202     C Note. 20th May 1998
203     C I made a weird discovery! In the model paper we argue
204     C for the form of the preconditioner used here ( see
205     C A Finite-volume, Incompressible Navier-Stokes Model
206     C ...., Marshall et. al ). The algebra gives a simple
207     C 0.5 factor for the averaging of ac and aCw to get a
208     C symmettric pre-conditioner. By using a factor of 0.51
209     C i.e. scaling the off-diagonal terms in the
210     C preconditioner down slightly I managed to get the
211     C number of iterations for convergence in a test case to
212     C drop form 192 -> 134! Need to investigate this further!
213     C For now I have introduced a parameter cg2dpcOffDFac which
214     C defaults to 0.51 but can be set at runtime.
215 cnh 1.24 C delP - Vertical grid spacing ( Pa ).
216     C delZ - Vertical grid spacing ( m ).
217 cnh 1.23 C delR - Vertical grid spacing ( units of r ).
218 cnh 1.1 C delX - Separation between cell faces (m) or (deg), depending
219     C delY on input flags.
220     C gravity - Accel. due to gravity ( m/s^2 )
221 cnh 1.32 C recip_gravity and its inverse
222 cnh 1.8 C gBaro - Accel. due to gravity used in barotropic equation ( m/s^2 )
223 cnh 1.1 C ronil - Reference density
224 cnh 1.24 C rhoConst - Vertically constant reference density
225 cnh 1.3 C startTime - Start time for model ( s )
226 cnh 1.1 C phiMin - Latitude of southern most cell face.
227     C thetaMin - Longitude of western most cell face (this
228     C is an "inert" parameter but it is included
229     C to make geographical references simple.)
230 cnh 1.3 C rSphere - Radius of sphere for a spherical polar grid ( m ).
231 cnh 1.25 C recip_RSphere - Reciprocal radius of sphere ( m ).
232 cnh 1.1 C f0 - Reference coriolis parameter ( 1/s )
233     C ( Southern edge f for beta plane )
234     C beta - df/dy ( s^-1.m^-1 )
235 cnh 1.15 C omega - Angular velocity ( rad/s )
236 cnh 1.1 C viscAh - Eddy viscosity coeff. for mixing of
237     C momentum laterally ( m^2/s )
238     C viscAz - Eddy viscosity coeff. for mixing of
239     C momentum vertically ( m^2/s )
240 cnh 1.27 C viscAp - Eddy viscosity coeff. for mixing of
241     C momentum vertically ( Pa^2/s )
242 cnh 1.24 C viscAr - Eddy viscosity coeff. for mixing of
243     C momentum vertically ( units of r^2/s )
244 cnh 1.1 C viscA4 - Biharmonic viscosity coeff. for mixing of
245     C momentum laterally ( m^4/s )
246     C diffKhT - Laplacian diffusion coeff. for mixing of
247     C heat laterally ( m^2/s )
248     C diffKzT - Laplacian diffusion coeff. for mixing of
249     C heat vertically ( m^2/s )
250 cnh 1.27 C diffKpT - Laplacian diffusion coeff. for mixing of
251     C heat vertically ( Pa^2/s )
252 cnh 1.24 C diffKrT - Laplacian diffusion coeff. for mixing of
253     C heat vertically ( units of r^2/s )
254 cnh 1.1 C diffK4T - Biharmonic diffusion coeff. for mixing of
255     C heat laterally ( m^4/s )
256     C diffKhS - Laplacian diffusion coeff. for mixing of
257     C salt laterally ( m^2/s )
258     C diffKzS - Laplacian diffusion coeff. for mixing of
259     C salt vertically ( m^2/s )
260 cnh 1.27 C diffKpS - Laplacian diffusion coeff. for mixing of
261     C salt vertically ( Pa^2/s )
262 cnh 1.24 C diffKrS - Laplacian diffusion coeff. for mixing of
263     C salt vertically ( units of r^2/s )
264 cnh 1.1 C diffK4S - Biharmonic diffusion coeff. for mixing of
265     C salt laterally ( m^4/s )
266 cnh 1.3 C deltaT - Default timestep ( s )
267 cnh 1.7 C deltaTClock - Timestep used as model "clock". This determines the
268     C IO frequencies and is used in tagging output. It can
269     C be totally different to the dynamical time. Typically
270     C it will be the deep-water timestep for accelerated runs.
271     C Frequency of checkpointing and dumping of the model state
272     C are referenced to this clock. ( s )
273 cnh 1.3 C deltaTMom - Timestep for momemtum equations ( s )
274     C deltaTtracer - Timestep for tracer equations ( s )
275 cnh 1.8 C freesurfFac - Parameter to turn implicit free surface term on or off
276     C freesurfac = 1. uses implicit free surface
277     C freesurfac = 0. uses rigid lid
278 cnh 1.24 C hFacMin - Minimum fraction size of a cell (affects hFacC etc...)
279     C hFacMinDz - Minimum dimesional size of a cell (affects hFacC etc..., m)
280 cnh 1.28 C hFacMinDp - Minimum dimesional size of a cell (affects hFacC etc..., Pa)
281 cnh 1.24 C hFacMinDr - Minimum dimesional size of a cell (affects hFacC etc..., units of r)
282     C tauCD - CD scheme coupling timescale ( 1/s )
283     C rCD - CD scheme normalised coupling parameter ( 0-1 )
284 cnh 1.27 C GMmaxslope - max. slope allowed in GM/Redi tensor
285     C GMlength - Length to use in Visbeck et al. formula for K (m)
286     C GMalpha - alpha to use in Visbeck et al. formula for K
287     C GMdepth - Depth over which to integrate Richardson # (Visbeck et al.)
288 cnh 1.24 C GMkbackground - background value of GM/Redi diffusion coefficient ( m^2/s )
289     C GMmaxval - max. value of KapGM allowed in GM/Redi scheme ( m^2/s )
290     C startTime - Starting time for this integration ( s ).
291     C endTime - Ending time for this integration ( s ).
292     C chkPtFreq - Frequency of rolling check pointing ( s ).
293     C pChkPtFreq - Frequency of permanent check pointing ( s ).
294     C dumpFreq - Frequency with which model state is written to
295     C post-processing files ( s ).
296     C afFacMom - Advection of momentum term tracer parameter
297     C vfFacMom - Momentum viscosity tracer parameter
298     C pfFacMom - Momentum pressure forcing tracer parameter
299     C cfFacMom - Coriolis term tracer parameter
300     C foFacMom - Momentum forcing tracer parameter
301     C mtFacMom - Metric terms tracer parameter
302     C cAdjFreq - Frequency of convective adjustment
303     C
304 cnh 1.27 C taveFreq - Frequency with which time-averaged model state is written to
305     C post-processing files ( s ).
306 cnh 1.18 C tauThetaClimRelax - Relaxation to climatology time scale ( s ).
307     C lambdaThetaClimRelax - Inverse time scale for relaxation ( 1/s ).
308     C tauSaltClimRelax - Relaxation to climatology time scale ( s ).
309     C lambdaSaltClimRelax - Inverse time scale for relaxation ( 1/s ).
310 adcroft 1.19 C externForcingPeriod - Is the period of which forcing varies (eg. 1 month)
311     C externForcingCycle - Is the repeat time of the forcing (eg. 1 year)
312     C (note: externForcingCycle must be an integer
313     C number times externForcingPeriod)
314 cnh 1.30 C horiVertRatio - Ratio on units in vertical to units in horizontal.
315 cnh 1.31 C recip_horiVertRatio ( 1 if horiz in m and vertical in m ).
316 cnh 1.30 C ( g*rho if horiz in m and vertical in Pa ).
317 cnh 1.24 COMMON /PARM_R/ cg2dTargetResidual, cg2dpcOffDFac,
318 cnh 1.33 & cg3dTargetResidual,
319 cnh 1.24 & delP, delZ, delR, delX, delY,
320 cnh 1.33 & deltaT,deltaTmom, deltaTtracer, deltaTClock,abeps, startTime,
321     & phiMin, thetaMin, rSphere, recip_RSphere, f0, fCori, beta,
322 cnh 1.24 & viscAh, viscAz, viscA4, viscAr,
323     & diffKhT, diffKzT, diffK4T, diffKrT,
324 cnh 1.25 & diffKhS, diffKzS, diffK4S, diffKrS,
325 cnh 1.27 & delT, tauCD, rCD, freeSurfFac, hFacMin, hFacMinDz,
326 adcroft 1.21 & GMmaxslope,GMlength,GMalpha,GMdepth,GMkbackground,GMmaxval,
327 cnh 1.32 & gravity, recip_Gravity, gBaro, rhonil, recip_rhonil,
328 cnh 1.27 & recip_rhoConst, rhoConst, tRef, sRef,
329 adcroft 1.20 & endTime, chkPtFreq, pchkPtFreq, dumpFreq, taveFreq,
330 cnh 1.15 & afFacMom, vfFacMom, pfFacMom, cfFacMom, foFacMom, mtFacMom,
331 cnh 1.18 & cAdjFreq, omega, tauThetaClimRelax, lambdaThetaClimRelax,
332 adcroft 1.19 & tauSaltClimRelax, lambdaSaltClimRelax,
333 cnh 1.27 & externForcingCycle, externForcingPeriod,
334 cnh 1.28 & viscAp, diffKpT, diffKpS, hFacMinDr, hFacMinDp,
335 cnh 1.31 & theta_S, specVol_S, horiVertRatio, recip_horiVertRatio
336 cnh 1.27
337 adcroft 1.6 _RL cg2dTargetResidual
338 cnh 1.33 _RL cg3dTargetResidual
339 cnh 1.7 _RL cg2dpcOffDFac
340 cnh 1.25 _RL delZ(Nr)
341     _RL delP(Nr)
342     _RL delR(Nr)
343 adcroft 1.6 _RL delX(Nx)
344     _RL delY(Ny)
345     _RL deltaT
346 cnh 1.7 _RL deltaTClock
347 adcroft 1.6 _RL deltaTmom
348     _RL deltaTtracer
349     _RL abeps
350     _RL phiMin
351     _RL thetaMin
352     _RL rSphere
353 cnh 1.25 _RL recip_RSphere
354 adcroft 1.6 _RL f0
355 cnh 1.8 _RL freeSurfFac
356 adcroft 1.22 _RL hFacMin
357     _RL hFacMinDz
358 cnh 1.28 _RL hFacMinDp
359 cnh 1.27 _RL hFacMinDr
360 adcroft 1.6 _RL beta
361     _RL viscAh
362     _RL viscAz
363 cnh 1.27 _RL viscAp
364 cnh 1.25 _RL viscAr
365 adcroft 1.6 _RL viscA4
366     _RL diffKhT
367 cnh 1.25 _RL diffKrT
368 adcroft 1.6 _RL diffKzT
369 cnh 1.27 _RL diffKpT
370 adcroft 1.6 _RL diffK4T
371     _RL diffKhS
372 cnh 1.25 _RL diffKrS
373 adcroft 1.6 _RL diffKzS
374 cnh 1.27 _RL diffKpS
375 adcroft 1.6 _RL diffK4S
376     _RL delt
377     _RL tauCD
378     _RL rCD
379     _RL GMmaxslope
380     _RL GMlength
381     _RL GMalpha
382     _RL GMdepth
383     _RL GMkbackground
384 adcroft 1.21 _RL GMmaxval
385 adcroft 1.6 _RL gravity
386 cnh 1.32 _RL recip_gravity
387 cnh 1.8 _RL gBaro
388 adcroft 1.6 _RL rhonil
389 cnh 1.25 _RL recip_rhonil
390     _RL rhoConst
391 cnh 1.26 _RL recip_rhoConst
392     _RL specVol_S(Nr)
393 cnh 1.25 _RL tRef(Nr)
394 cnh 1.26 _RL theta_S(Nr)
395 cnh 1.25 _RL sRef(Nr)
396 cnh 1.11 _RS Fcori(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
397 adcroft 1.6 _RL startTime
398     _RL endTime
399     _RL chkPtFreq
400 cnh 1.7 _RL pChkPtFreq
401 adcroft 1.6 _RL dumpFreq
402 adcroft 1.20 _RL taveFreq
403 cnh 1.9 _RL afFacMom
404     _RL vfFacMom
405     _RL pfFacMom
406     _RL cfFacMom
407     _RL foFacMom
408 cnh 1.15 _RL mTFacMom
409 cnh 1.9 _RL cAdjFreq
410 cnh 1.15 _RL omega
411 cnh 1.18 _RL tauThetaClimRelax
412     _RL lambdaThetaClimRelax
413     _RL tauSaltClimRelax
414     _RL lambdaSaltClimRelax
415 adcroft 1.19 _RL externForcingCycle
416     _RL externForcingPeriod
417 cnh 1.30 _RL horiVertRatio
418 cnh 1.31 _RL recip_horiVertRatio
419 adcroft 1.6
420 cnh 1.1 COMMON /PARM_A/ HeatCapacity_Cp,
421     & Lamba_theta
422 adcroft 1.6 _RL HeatCapacity_Cp
423     _RL Lamba_theta
424 adcroft 1.5
425     C Equation of State (polynomial coeffients)
426     COMMON /PARM_EOS_NL/ eosC,eosSig0,eosRefT,eosRefS
427 cnh 1.25 _RL eosC(9,Nr+1),eosSig0(Nr+1),eosRefT(Nr+1),eosRefS(Nr+1)
428 adcroft 1.5 C Linear equation of state
429     C tAlpha - Linear EOS thermal expansion coefficient ( 1/degree ).
430     C sBeta - Linear EOS haline contraction coefficient.
431 adcroft 1.12 COMMON /PARM_EOS_LIN/ tAlpha,sBeta,eosType
432 adcroft 1.6 _RL tAlpha
433     _RL sBeta
434 adcroft 1.12 character*(6) eosType
435 cnh 1.7

  ViewVC Help
Powered by ViewVC 1.1.22