/[MITgcm]/MITgcm/model/src/ini_parms.F
ViewVC logotype

Annotation of /MITgcm/model/src/ini_parms.F

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


Revision 1.44 - (hide annotations) (download)
Tue Jun 29 18:33:26 1999 UTC (24 years, 10 months ago) by adcroft
Branch: MAIN
Changes since 1.43: +2 -2 lines
Added COS(latitude)^cosPower dependence to viscosity terms.
New parameter "cosPower" defaults to 0.0 which removes latitudinal
dependence.

1 adcroft 1.44 C $Header: /u/gcmpack/models/MITgcmUV/model/src/ini_parms.F,v 1.43 1999/05/18 18:01:13 adcroft Exp $
2 cnh 1.1
3 adcroft 1.22 #include "CPP_OPTIONS.h"
4 cnh 1.1
5     SUBROUTINE INI_PARMS( myThid )
6     C /==========================================================\
7     C | SUBROUTINE INI_PARMS |
8     C | o Routine to set model "parameters" |
9     C |==========================================================|
10     C | Notes: |
11     C | ====== |
12     C | The present version of this routine is a place-holder. |
13     C | A production version needs to handle parameters from an |
14     C | external file and possibly reading in some initial field |
15     C | values. |
16     C \==========================================================/
17 adcroft 1.38 IMPLICIT NONE
18 cnh 1.1
19     C === Global variables ===
20     #include "SIZE.h"
21     #include "EEPARAMS.h"
22     #include "PARAMS.h"
23 cnh 1.28 #include "GRID.h"
24 cnh 1.1 #include "CG2D.h"
25 adcroft 1.43 #ifdef ALLOW_KPP
26     #include "KPPMIX.h"
27     #endif
28 cnh 1.1
29     C === Routine arguments ===
30     C myThid - Number of this instance of INI_PARMS
31     INTEGER myThid
32    
33     C === Local variables ===
34     C dxSpacing, dySpacing - Default spacing in X and Y.
35     C Units are that of coordinate system
36     C i.e. cartesian => metres
37     C s. polar => degrees
38     C goptCount - Used to count the nuber of grid options
39     C (only one is allowed! )
40     C msgBuf - Informational/error meesage buffer
41     C errIO - IO error flag
42     C iUnit - Work variable for IO unit number
43     C record - Work variable for IO buffer
44     C K, I, J - Loop counters
45 cnh 1.28 C xxxDefault - Default value for variable xxx
46     _RL dxSpacing
47     _RL dySpacing
48 cnh 1.1 CHARACTER*(MAX_LEN_MBUF) msgBuf
49     CHARACTER*(MAX_LEN_PREC) record
50     INTEGER goptCount
51     INTEGER K, I, J, IL, iUnit
52     INTEGER errIO
53     INTEGER IFNBLNK
54     EXTERNAL IFNBLNK
55     INTEGER ILNBLNK
56     EXTERNAL ILNBLNK
57 cnh 1.28 C Default values for variables which have vertical coordinate system
58     C dependency.
59     _RL viscArDefault
60     _RL diffKrTDefault
61     _RL diffKrSDefault
62     _RL hFacMinDrDefault
63 adcroft 1.39 _RL delRDefault(Nr)
64 adcroft 1.41 _RS rkFacDefault
65 cnh 1.28 C zCoordInputData - These are used to select between different coordinate systems.
66     C pCoordInputData The vertical coordinate system in the rest of the model is
67     C rCoordInputData written in terms of r. In the model "data" file input data can
68     C coordsSet be interms of z, p or r.
69     C e.g. delZ or delP or delR for the vertical grid spacing.
70     C The following rules apply:
71     C All parameters must use the same vertical coordinate system.
72     C e.g. delZ and viscAz is legal but
73     C delZ and viscAr is an error.
74     C Similarly specifyinh delZ and delP is an error.
75     C zCoord..., pCoord..., rCoord... are used to flag when z, p or r are
76     C used. coordsSet counts how many vertical coordinate systems have been
77     C used to specify variables. coordsSet > 1 is an error.
78     C
79     LOGICAL zCoordInputData
80     LOGICAL pCoordInputData
81     LOGICAL rCoordInputData
82     INTEGER coordsSet
83 cnh 1.1
84     C-- Continuous equation parameters
85     NAMELIST /PARM01/
86 cnh 1.8 & gravity, gBaro, rhonil, tAlpha, sBeta, f0, beta,
87 adcroft 1.44 & viscAh, viscAz, viscA4, cosPower,
88 cnh 1.27 & diffKhT, diffKzT, diffK4T,
89 cnh 1.1 & diffKhS, diffKzS, diffK4S,
90 adcroft 1.23 & GMmaxslope,GMlength,GMalpha,GMdepth,GMkbackground,GMmaxval,
91 adcroft 1.11 & tRef, sRef, eosType,
92 adcroft 1.39 & no_slip_sides,no_slip_bottom,
93 cnh 1.1 & momViscosity, momAdvection, momForcing, useCoriolis,
94 cnh 1.14 & momPressureForcing, metricTerms,
95 cnh 1.1 & tempDiffusion, tempAdvection, tempForcing,
96 cnh 1.8 & saltDiffusion, saltAdvection, saltForcing,
97 adcroft 1.24 & implicitFreeSurface, rigidLid, freeSurfFac, hFacMin, hFacMinDz,
98 cnh 1.27 & tempStepping, saltStepping, momStepping, implicitDiffusion,
99     & viscAr, diffKrT, diffKrS, hFacMinDr,
100 cnh 1.29 & viscAp, diffKpT, diffKpS, hFacMinDp,
101 adcroft 1.41 & rhoConst, buoyancyRelation, HeatCapacity_Cp,
102 adcroft 1.40 & writeBinaryPrec, readBinaryPrec, writeStatePrec,
103 adcroft 1.42 & openBoundaries, nonHydrostatic, globalFiles,
104     & allowFreezing
105 cnh 1.1
106     C-- Elliptic solver parameters
107     NAMELIST /PARM02/
108 cnh 1.34 & cg2dMaxIters, cg2dChkResFreq, cg2dTargetResidual, cg2dpcOffDFac,
109     & cg3dMaxIters, cg3dChkResFreq, cg3dTargetResidual
110 cnh 1.1
111     C-- Time stepping parammeters
112     NAMELIST /PARM03/
113 cnh 1.34 & nIter0, nTimeSteps, deltaT, deltaTmom, deltaTtracer, abEps,
114     & tauCD, rCD,
115 adcroft 1.20 & startTime, endTime, chkPtFreq, dumpFreq, taveFreq, deltaTClock,
116     & pChkPtFreq, cAdjFreq, tauThetaClimRelax, tauSaltClimRelax,
117 adcroft 1.19 & periodicExternalForcing, externForcingPeriod, externForcingCycle
118 cnh 1.1
119     C-- Gridding parameters
120     NAMELIST /PARM04/
121     & usingCartesianGrid, delZ, dxSpacing, dySpacing, delX, delY,
122     & usingSphericalPolarGrid, phiMin, thetaMin, rSphere,
123 cnh 1.28 & l, m, n, delP, delR, rkFac
124 cnh 1.1
125 cnh 1.15 C-- Input files
126     NAMELIST /PARM05/
127     & bathyFile, hydrogThetaFile, hydrogSaltFile,
128 adcroft 1.41 & zonalWindFile, meridWindFile,
129     & thetaClimFile, saltClimFile,
130     & surfQfile, EmPmRfile
131 cnh 1.15
132 adcroft 1.37 C-- Open Boundaries
133     NAMELIST /PARM06/
134     & OB_Jnorth, OB_Jsouth, OB_Ieast, OB_Iwest
135 cnh 1.28
136 adcroft 1.43 #ifdef ALLOW_KPP
137     C-- KPP vertical mixing parameters
138     NAMELIST /PARM07/
139     & usingKPPmixing, KPPmixingMaps, KPPwriteState,
140     & minKPPviscAz, maxKPPviscAz, minKPPghat, maxKPPghat,
141     & minKPPdiffKzT, maxKPPdiffKzT, minKPPdiffKzS, maxKPPdiffKzS,
142     & epsln, epsilon, vonk, conc1, conam, concm, conc2, zetam,
143     & conas, concs, conc3, zetas,
144     & Ricr, cekman, cmonob, concv, hbf, Vtc,
145     & zmin, zmax, umin, umax,
146     & num_v_smooth_Ri, num_v_smooth_BV,
147     & num_z_smooth_sh, num_m_smooth_sh,
148     & Riinfty, BVSQcon, difm0, difs0, dift0,
149     & difmiw, difsiw, diftiw, difmcon, difscon, diftcon,
150     & cstar, cg
151     #endif
152    
153 cnh 1.1 C
154     _BEGIN_MASTER(myThid)
155    
156 adcroft 1.39 C Defaults values for input parameters
157     CALL SET_DEFAULTS(
158     O viscArDefault, diffKrTDefault, diffKrSDefault,
159 adcroft 1.41 O hFacMinDrDefault, delRdefault, rkFacDefault,
160 adcroft 1.39 I myThid )
161    
162 cnh 1.28 C-- Initialise "which vertical coordinate system used" flags.
163     zCoordInputData = .FALSE.
164     pCoordInputData = .FALSE.
165     rCoordInputData = .FALSE.
166 cnh 1.29 usingPCoords = .FALSE.
167     usingZCoords = .FALSE.
168 cnh 1.28 coordsSet = 0
169    
170 cnh 1.1 C-- Open the parameter file
171     OPEN(UNIT=scrUnit1,STATUS='SCRATCH')
172     OPEN(UNIT=scrUnit2,STATUS='SCRATCH')
173 cnh 1.34 OPEN(UNIT=modelDataUnit,FILE='data',STATUS='OLD',
174 cnh 1.35 & IOSTAT=errIO)
175     IF ( errIO .LT. 0 ) THEN
176 cnh 1.1 WRITE(msgBuf,'(A)')
177     & 'S/R INI_PARMS'
178     CALL PRINT_ERROR( msgBuf , 1)
179     WRITE(msgBuf,'(A)')
180     & 'Unable to open model parameter'
181     CALL PRINT_ERROR( msgBuf , 1)
182     WRITE(msgBuf,'(A)')
183     & 'file "data"'
184     CALL PRINT_ERROR( msgBuf , 1)
185     CALL MODELDATA_EXAMPLE( myThid )
186     STOP 'ABNORMAL END: S/R INI_PARMS'
187 cnh 1.35 ENDIF
188 cnh 1.1
189 cnh 1.35 DO WHILE ( .TRUE. )
190     READ(modelDataUnit,FMT='(A)',END=1001) RECORD
191     IL = MAX(ILNBLNK(RECORD),1)
192     IF ( RECORD(1:1) .NE. commentCharacter )
193     & WRITE(UNIT=scrUnit1,FMT='(A)') RECORD(:IL)
194     WRITE(UNIT=scrUnit2,FMT='(A)') RECORD(:IL)
195     ENDDO
196 cnh 1.1 1001 CONTINUE
197     CLOSE(modelDataUnit)
198    
199     C-- Report contents of model parameter file
200     WRITE(msgBuf,'(A)')
201     &'// ======================================================='
202 cnh 1.34 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
203     & SQUEEZE_RIGHT , 1)
204 cnh 1.1 WRITE(msgBuf,'(A)') '// Model parameter file "data"'
205 cnh 1.34 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
206     & SQUEEZE_RIGHT , 1)
207 cnh 1.1 WRITE(msgBuf,'(A)')
208     &'// ======================================================='
209     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
210     & SQUEEZE_RIGHT , 1)
211     iUnit = scrUnit2
212     REWIND(iUnit)
213 cnh 1.35 DO WHILE ( .TRUE. )
214 cnh 1.1 READ(UNIT=iUnit,FMT='(A)',END=2001) RECORD
215     IL = MAX(ILNBLNK(RECORD),1)
216     WRITE(msgBuf,'(A,A)') '>',RECORD(:IL)
217 cnh 1.34 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
218     & SQUEEZE_RIGHT , 1)
219 cnh 1.35 ENDDO
220 cnh 1.1 2001 CONTINUE
221     CLOSE(iUnit)
222     WRITE(msgBuf,'(A)') ' '
223     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
224     & SQUEEZE_RIGHT , 1)
225    
226    
227     C-- Read settings from model parameter file "data".
228     iUnit = scrUnit1
229     REWIND(iUnit)
230    
231     C-- Set default "physical" parameters
232 cnh 1.28 viscAz = UNSET_RL
233     viscAr = UNSET_RL
234     viscAp = UNSET_RL
235     diffKzT = UNSET_RL
236     diffKpT = UNSET_RL
237     diffKrT = UNSET_RL
238     diffKzS = UNSET_RL
239     diffKpS = UNSET_RL
240     diffKrS = UNSET_RL
241 adcroft 1.39 gBaro = UNSET_RL
242     rhoConst = UNSET_RL
243 cnh 1.28 hFacMinDr = UNSET_RL
244     hFacMinDz = UNSET_RL
245     hFacMinDp = UNSET_RL
246 adcroft 1.41 READ(UNIT=iUnit,NML=PARM01) !,IOSTAT=errIO)
247 cnh 1.35 IF ( errIO .LT. 0 ) THEN
248 cnh 1.1 WRITE(msgBuf,'(A)')
249     & 'S/R INI_PARMS'
250     CALL PRINT_ERROR( msgBuf , 1)
251     WRITE(msgBuf,'(A)')
252     & 'Error reading numerical model '
253     CALL PRINT_ERROR( msgBuf , 1)
254     WRITE(msgBuf,'(A)')
255     & 'parameter file "data"'
256     CALL PRINT_ERROR( msgBuf , 1)
257     WRITE(msgBuf,'(A)')
258     & 'Problem in namelist PARM01'
259     CALL PRINT_ERROR( msgBuf , 1)
260     CALL MODELDATA_EXAMPLE( myThid )
261     STOP 'ABNORMAL END: S/R INI_PARMS'
262 cnh 1.35 ENDIF
263 cnh 1.28 IF ( implicitFreeSurface ) freeSurfFac = 1.D0
264     IF ( rigidLid ) freeSurfFac = 0.D0
265 adcroft 1.39 IF ( gBaro .EQ. UNSET_RL ) gBaro=gravity
266     IF ( rhoConst .EQ. UNSET_RL ) rhoConst=rhoNil
267 cnh 1.28 C-- Momentum viscosity on/off flag.
268 cnh 1.9 IF ( momViscosity ) THEN
269 cnh 1.28 vfFacMom = 1.D0
270 cnh 1.9 ELSE
271 cnh 1.28 vfFacMom = 0.D0
272 cnh 1.9 ENDIF
273 cnh 1.28 C-- Momentum advection on/off flag.
274 cnh 1.9 IF ( momAdvection ) THEN
275 cnh 1.28 afFacMom = 1.D0
276 cnh 1.9 ELSE
277 cnh 1.28 afFacMom = 0.D0
278 cnh 1.9 ENDIF
279 cnh 1.28 C-- Momentum forcing on/off flag.
280 cnh 1.9 IF ( momForcing ) THEN
281 cnh 1.28 foFacMom = 1.D0
282 cnh 1.9 ELSE
283 cnh 1.28 foFacMom = 0.D0
284 cnh 1.9 ENDIF
285 cnh 1.28 C-- Coriolis term on/off flag.
286 cnh 1.9 IF ( useCoriolis ) THEN
287 cnh 1.28 cfFacMom = 1.D0
288 cnh 1.9 ELSE
289 cnh 1.28 cfFacMom = 0.D0
290 cnh 1.9 ENDIF
291 cnh 1.28 C-- Pressure term on/off flag.
292 cnh 1.9 IF ( momPressureForcing ) THEN
293 cnh 1.28 pfFacMom = 1.D0
294 cnh 1.9 ELSE
295 cnh 1.28 pfFacMom = 0.D0
296 cnh 1.9 ENDIF
297 cnh 1.28 C-- Metric terms on/off flag.
298 cnh 1.14 IF ( metricTerms ) THEN
299 cnh 1.28 mTFacMom = 1.D0
300 cnh 1.14 ELSE
301 cnh 1.28 mTFacMom = 1.D0
302 cnh 1.14 ENDIF
303 cnh 1.28 C-- z,p,r coord input switching.
304     IF ( viscAz .NE. UNSET_RL ) zCoordInputData = .TRUE.
305     IF ( viscAp .NE. UNSET_RL ) pCoordInputData = .TRUE.
306     IF ( viscAr .NE. UNSET_RL ) rCoordInputData = .TRUE.
307     IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAz
308     IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAp
309     IF ( viscAr .EQ. UNSET_RL ) viscAr = viscArDefault
310    
311     IF ( diffKzT .NE. UNSET_RL ) zCoordInputData = .TRUE.
312     IF ( diffKpT .NE. UNSET_RL ) pCoordInputData = .TRUE.
313     IF ( diffKrT .NE. UNSET_RL ) rCoordInputData = .TRUE.
314     IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKzT
315     IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKpT
316     IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKrTDefault
317    
318     IF ( diffKzS .NE. UNSET_RL ) zCoordInputData = .TRUE.
319     IF ( diffKpS .NE. UNSET_RL ) pCoordInputData = .TRUE.
320     IF ( diffKrS .NE. UNSET_RL ) rCoordInputData = .TRUE.
321     IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKzS
322     IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKpS
323     IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKrSDefault
324    
325     IF ( hFacMinDz .NE. UNSET_RL ) zCoordInputData = .TRUE.
326     IF ( hFacMinDp .NE. UNSET_RL ) pCoordInputData = .TRUE.
327     IF ( hFacMinDr .NE. UNSET_RL ) rCoordInputData = .TRUE.
328     IF ( hFacMinDz .EQ. UNSET_RL ) hFacMinDr = hFacMinDz
329     IF ( hFacMinDp .EQ. UNSET_RL ) hFacMinDr = hFacMinDp
330     IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDrDefault
331 cnh 1.8
332     IF ( implicitFreeSurface .AND. rigidLid ) THEN
333 cnh 1.34 WRITE(msgBuf,'(A,A)')
334     & 'S/R INI_PARMS: Cannot select both implicitFreeSurface',
335     & ' and rigidLid.'
336 cnh 1.28 CALL PRINT_ERROR( msgBuf , myThid)
337     STOP 'ABNORMAL END: S/R INI_PARMS'
338     ENDIF
339     coordsSet = 0
340     IF ( zCoordInputData ) coordsSet = coordsSet + 1
341     IF ( pCoordInputData ) coordsSet = coordsSet + 1
342     IF ( rCoordInputData ) coordsSet = coordsSet + 1
343     IF ( coordsSet .GT. 1 ) THEN
344     WRITE(msgBuf,'(A)')
345     & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
346 cnh 1.8 CALL PRINT_ERROR( msgBuf , myThid)
347     STOP 'ABNORMAL END: S/R INI_PARMS'
348     ENDIF
349 cnh 1.28 IF ( rhoConst .LE. 0. ) THEN
350     WRITE(msgBuf,'(A)')
351     & 'S/R INI_PARMS: rhoConst must be greater than 0.'
352     CALL PRINT_ERROR( msgBuf , myThid)
353     STOP 'ABNORMAL END: S/R INI_PARMS'
354     ELSE
355     recip_rhoConst = 1.D0 / rhoConst
356 adcroft 1.38 ENDIF
357     IF ( rhoNil .LE. 0. ) THEN
358     WRITE(msgBuf,'(A)')
359     & 'S/R INI_PARMS: rhoNil must be greater than 0.'
360     CALL PRINT_ERROR( msgBuf , myThid)
361     STOP 'ABNORMAL END: S/R INI_PARMS'
362     ELSE
363     recip_rhoNil = 1.D0 / rhoNil
364 cnh 1.33 ENDIF
365 adcroft 1.39 IF ( HeatCapacity_Cp .LE. 0. ) THEN
366     WRITE(msgBuf,'(A)')
367     & 'S/R INI_PARMS: HeatCapacity_Cp must be greater than 0.'
368     CALL PRINT_ERROR( msgBuf , myThid)
369     STOP 'ABNORMAL END: S/R INI_PARMS'
370     ELSE
371     recip_Cp = 1.D0 / HeatCapacity_Cp
372     ENDIF
373 cnh 1.33 IF ( gravity .LE. 0. ) THEN
374     WRITE(msgBuf,'(A)')
375     & 'S/R INI_PARMS: gravity must be greater than 0.'
376     CALL PRINT_ERROR( msgBuf , myThid)
377     STOP 'ABNORMAL END: S/R INI_PARMS'
378     ELSE
379     recip_gravity = 1.D0 / gravity
380 cnh 1.28 ENDIF
381 adcroft 1.42 C Set globalFiles flag for READ_WRITE_FLD package
382     CALL SET_WRITE_GLOBAL_FLD( globalFiles )
383     C Set globalFiles flag for READ_WRITE_REC package
384     CALL SET_WRITE_GLOBAL_REC( globalFiles )
385     C Set globalFiles flag for READ_WRITE_REC package
386     CALL SET_WRITE_GLOBAL_PICKUP( globalFiles )
387 cnh 1.1
388     C-- Elliptic solver parameters
389 adcroft 1.41 READ(UNIT=iUnit,NML=PARM02) !,IOSTAT=errIO)
390 cnh 1.35 IF ( errIO .LT. 0 ) THEN
391 cnh 1.1 WRITE(msgBuf,'(A)')
392     & 'S/R INI_PARMS'
393     CALL PRINT_ERROR( msgBuf , 1)
394     WRITE(msgBuf,'(A)')
395     & 'Error reading numerical model '
396     CALL PRINT_ERROR( msgBuf , 1)
397     WRITE(msgBuf,'(A)')
398     & 'parameter file "data".'
399     CALL PRINT_ERROR( msgBuf , 1)
400     WRITE(msgBuf,'(A)')
401     & 'Problem in namelist PARM02'
402     CALL PRINT_ERROR( msgBuf , 1)
403     CALL MODELDATA_EXAMPLE( myThid )
404     STOP 'ABNORMAL END: S/R INI_PARMS'
405 cnh 1.35 ENDIF
406 cnh 1.1
407     C-- Time stepping parameters
408 cnh 1.28 rCD = -1.D0
409 adcroft 1.41 READ(UNIT=iUnit,NML=PARM03) !,IOSTAT=errIO)
410 cnh 1.35 IF ( errIO .LT. 0 ) THEN
411 cnh 1.1 WRITE(msgBuf,'(A)')
412     & 'S/R INI_PARMS'
413     CALL PRINT_ERROR( msgBuf , 1)
414     WRITE(msgBuf,'(A)')
415     & 'Error reading numerical model '
416     CALL PRINT_ERROR( msgBuf , 1)
417     WRITE(msgBuf,'(A)')
418     & 'parameter file "data"'
419     CALL PRINT_ERROR( msgBuf , 1)
420     WRITE(msgBuf,'(A)')
421     & 'Problem in namelist PARM03'
422     CALL PRINT_ERROR( msgBuf , 1)
423     CALL MODELDATA_EXAMPLE( myThid )
424     STOP 'ABNORMAL END: S/R INI_PARMS'
425 cnh 1.35 ENDIF
426 cnh 1.4 C Process "timestepping" params
427     C o Time step size
428     IF ( deltaT .EQ. 0. ) deltaT = deltaTmom
429     IF ( deltaT .EQ. 0. ) deltaT = deltaTtracer
430     IF ( deltaTmom .EQ. 0. ) deltaTmom = deltaT
431     IF ( deltaTtracer .EQ. 0. ) deltaTtracer = deltaT
432 cnh 1.7 IF ( deltaTClock .EQ. 0. ) deltaTClock = deltaT
433 adcroft 1.19 IF ( periodicExternalForcing ) THEN
434     IF ( externForcingCycle*externForcingPeriod .EQ. 0. ) THEN
435     WRITE(msgBuf,'(A)')
436     & 'S/R INI_PARMS: externForcingCycle,externForcingPeriod =0'
437     CALL PRINT_ERROR( msgBuf , 1)
438     STOP 'ABNORMAL END: S/R INI_PARMS'
439     ENDIF
440     IF ( INT(externForcingCycle/externForcingPeriod) .NE.
441     & externForcingCycle/externForcingPeriod ) THEN
442     WRITE(msgBuf,'(A)')
443     & 'S/R INI_PARMS: externForcingCycle <> N*externForcingPeriod'
444     CALL PRINT_ERROR( msgBuf , 1)
445     STOP 'ABNORMAL END: S/R INI_PARMS'
446     ENDIF
447     IF ( externForcingCycle.le.externForcingPeriod ) THEN
448     WRITE(msgBuf,'(A)')
449     & 'S/R INI_PARMS: externForcingCycle < externForcingPeriod'
450     CALL PRINT_ERROR( msgBuf , 1)
451     STOP 'ABNORMAL END: S/R INI_PARMS'
452     ENDIF
453     IF ( externForcingPeriod.lt.deltaTclock ) THEN
454     WRITE(msgBuf,'(A)')
455     & 'S/R INI_PARMS: externForcingPeriod < deltaTclock'
456     CALL PRINT_ERROR( msgBuf , 1)
457     STOP 'ABNORMAL END: S/R INI_PARMS'
458     ENDIF
459     ENDIF
460 cnh 1.9 C o Convection frequency
461     IF ( cAdjFreq .LT. 0. ) THEN
462     cAdjFreq = deltaTClock
463     ENDIF
464 cnh 1.14 C o CD coupling
465 cnh 1.28 IF ( tauCD .EQ. 0.D0 ) THEN
466 cnh 1.14 tauCD = deltaTmom
467     ENDIF
468     IF ( rCD .LT. 0. ) THEN
469     rCD = 1. - deltaTMom/tauCD
470     ENDIF
471 cnh 1.18 C o Temperature climatology relaxation time scale
472 cnh 1.28 IF ( tauThetaClimRelax .EQ. 0.D0 ) THEN
473 cnh 1.18 doThetaClimRelax = .FALSE.
474 cnh 1.28 lambdaThetaClimRelax = 0.D0
475 cnh 1.18 ELSE
476     doThetaClimRelax = .TRUE.
477     lambdaThetaClimRelax = 1./tauThetaClimRelax
478     ENDIF
479     C o Salinity climatology relaxation time scale
480 cnh 1.28 IF ( tauSaltClimRelax .EQ. 0.D0 ) THEN
481 cnh 1.18 doSaltClimRelax = .FALSE.
482 cnh 1.28 lambdaSaltClimRelax = 0.D0
483 cnh 1.18 ELSE
484     doSaltClimRelax = .TRUE.
485     lambdaSaltClimRelax = 1./tauSaltClimRelax
486     ENDIF
487 adcroft 1.41
488     C o Start time
489     IF ( nIter0 .NE. 0 .AND. startTime .EQ. 0. )
490     & startTime = deltaTClock*float(nIter0)
491     C o nIter0
492     IF ( nIter0 .EQ. 0 .AND. startTime .NE. 0. )
493     & nIter0 = INT( startTime/deltaTClock )
494     C o End Time
495     IF ( nTimeSteps .NE. 0 .AND. endTime .EQ. 0. )
496     & endTime = startTime + deltaTClock*float(nTimeSteps)
497     C o nTimeSteps
498     IF ( nTimeSteps .EQ. 0 .AND. endTime .NE. 0. )
499     & nTimeSteps = INT((endTime-startTime)/deltaTClock)
500     C o Consistent?
501     IF ( endTime-startTime .NE. deltaTClock*float(nTimeSteps) ) THEN
502     WRITE(msgBuf,'(A)')
503     & 'S/R INI_PARMS: both endTime and nTimeSteps have been set'
504     CALL PRINT_ERROR( msgBuf , 1)
505     WRITE(msgBuf,'(A)')
506     & 'S/R INI_PARMS: but are inconsistent'
507     CALL PRINT_ERROR( msgBuf , 1)
508     STOP 'ABNORMAL END: S/R INI_PARMS'
509 cnh 1.4 ENDIF
510 adcroft 1.21
511     C o If taveFreq is finite, then we must make sure the diagnostics
512     C code is being compiled
513 cnh 1.36 #ifndef INCLUDE_DIAGNOSTICS_INTERFACE_CODE
514 adcroft 1.21 IF (taveFreq.NE.0.) THEN
515     WRITE(msgBuf,'(A)')
516     & 'S/R INI_PARMS: taveFreq <> 0 but you have'
517     CALL PRINT_ERROR( msgBuf , 1)
518     WRITE(msgBuf,'(A)')
519     & 'not compiled the model with the diagnostics routines.'
520     CALL PRINT_ERROR( msgBuf , 1)
521 cnh 1.36 WRITE(msgBuf,'(A,A)')
522     & 'Re-compile with: #define INCLUDE_DIAGNOSTICS_INTERFACE_CODE',
523     & ' or -DINCLUDE_DIAGNOSTICS_INTERFACE_CODE'
524 adcroft 1.21 CALL PRINT_ERROR( msgBuf , 1)
525     STOP 'ABNORMAL END: S/R INI_PARMS'
526     ENDIF
527     #endif
528 cnh 1.1
529     C-- Grid parameters
530     C In cartesian coords distances are in metres
531 adcroft 1.41 rkFac = UNSET_RS
532 cnh 1.26 DO K =1,Nr
533 cnh 1.28 delZ(K) = UNSET_RL
534     delP(K) = UNSET_RL
535     delR(K) = UNSET_RL
536 cnh 1.1 ENDDO
537     C In spherical polar distances are in degrees
538 cnh 1.28 recip_rSphere = 1.D0/rSphere
539 adcroft 1.39 dxSpacing = UNSET_RL
540     dySpacing = UNSET_RL
541 adcroft 1.41 READ(UNIT=iUnit,NML=PARM04) !,IOSTAT=errIO)
542 cnh 1.35 IF ( errIO .LT. 0 ) THEN
543 cnh 1.1 WRITE(msgBuf,'(A)')
544     & 'S/R INI_PARMS'
545     CALL PRINT_ERROR( msgBuf , 1)
546     WRITE(msgBuf,'(A)')
547     & 'Error reading numerical model '
548     CALL PRINT_ERROR( msgBuf , 1)
549     WRITE(msgBuf,'(A)')
550     & 'parameter file "data"'
551     CALL PRINT_ERROR( msgBuf , 1)
552     WRITE(msgBuf,'(A)')
553     & 'Problem in namelist PARM04'
554     CALL PRINT_ERROR( msgBuf , 1)
555     CALL MODELDATA_EXAMPLE( myThid )
556     STOP 'ABNORMAL END: S/R INI_PARMS'
557 cnh 1.35 ENDIF
558 cnh 1.28 C
559 adcroft 1.39 IF ( dxSpacing .NE. UNSET_RL ) THEN
560     DO i=1,Nx
561     delX(i) = dxSpacing
562     ENDDO
563     ENDIF
564     IF ( dySpacing .NE. UNSET_RL ) THEN
565     DO j=1,Ny
566     delY(j) = dySpacing
567     ENDDO
568     ENDIF
569 cnh 1.14 IF ( rSphere .NE. 0 ) THEN
570 cnh 1.28 recip_rSphere = 1.D0/rSphere
571 cnh 1.14 ELSE
572 cnh 1.26 recip_rSphere = 0.
573 cnh 1.14 ENDIF
574 cnh 1.28 C-- Initialize EOS coefficients (3rd order polynomial)
575 adcroft 1.11 IF (eostype.eq.'POLY3') THEN
576     OPEN(37,FILE='POLY3.COEFFS',STATUS='OLD',FORM='FORMATTED')
577     READ(37,*) I
578 cnh 1.26 IF (I.NE.Nr) THEN
579 cnh 1.28 WRITE(msgBuf,'(A)')
580     & 'ini_parms: attempt to read POLY3.COEFFS failed'
581     CALL PRINT_ERROR( msgBuf , 1)
582     WRITE(msgBuf,'(A)')
583     & ' because bad # of levels in data'
584     CALL PRINT_ERROR( msgBuf , 1)
585 adcroft 1.11 STOP 'Bad data in POLY3.COEFFS'
586     ENDIF
587 cnh 1.26 READ(37,*) (eosRefT(K),eosRefS(K),eosSig0(K),K=1,Nr)
588     DO K=1,Nr
589 adcroft 1.11 READ(37,*) (eosC(I,K),I=1,9)
590     ENDDO
591     CLOSE(37)
592     ENDIF
593 cnh 1.28 C-- Check for conflicting grid definitions.
594 cnh 1.1 goptCount = 0
595     IF ( usingCartesianGrid ) goptCount = goptCount+1
596     IF ( usingSphericalPolarGrid ) goptCount = goptCount+1
597     IF ( goptCount .NE. 1 ) THEN
598     WRITE(msgBuf,'(A)')
599     & 'S/R INI_PARMS: More than one coordinate system requested'
600     CALL PRINT_ERROR( msgBuf , myThid)
601     STOP 'ABNORMAL END: S/R INI_PARMS'
602 cnh 1.14 ENDIF
603 cnh 1.28 C-- Make metric term settings consistent with underlying grid.
604 cnh 1.14 IF ( usingCartesianGrid ) THEN
605     usingSphericalPolarMterms = .FALSE.
606     metricTerms = .FALSE.
607     mTFacMom = 0
608 cnh 1.18 useBetaPlaneF = .TRUE.
609 cnh 1.14 ENDIF
610     IF ( usingSphericalPolarGrid ) THEN
611     useConstantF = .FALSE.
612     useBetaPlaneF = .FALSE.
613     useSphereF = .TRUE.
614 cnh 1.28 omega = 2.D0 * PI / ( 3600.D0 * 24.D0 )
615 cnh 1.14 usingSphericalPolarMterms = .TRUE.
616     metricTerms = .TRUE.
617     mTFacMom = 1
618 cnh 1.1 ENDIF
619 cnh 1.28 C-- p, z, r coord parameters
620     DO K = 1, Nr
621     IF ( delZ(K) .NE. UNSET_RL ) zCoordInputData = .TRUE.
622     IF ( delP(K) .NE. UNSET_RL ) pCoordInputData = .TRUE.
623     IF ( delR(K) .NE. UNSET_RL ) rCoordInputData = .TRUE.
624     IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delZ(K)
625     IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delP(K)
626 adcroft 1.39 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delRDefault(K)
627 adcroft 1.41 IF ( delR(K) .EQ. 0. ) THEN
628     WRITE(msgBuf,'(A,I4)')
629     & 'S/R INI_PARMS: No value for delZ/delP/delR at K = ',K
630     CALL PRINT_ERROR( msgBuf , 1)
631     STOP 'ABNORMAL END: S/R INI_PARMS'
632     ENDIF
633 cnh 1.28 ENDDO
634     C Check for multiple coordinate systems
635 adcroft 1.39 CoordsSet = 0
636 cnh 1.28 IF ( zCoordInputData ) coordsSet = coordsSet + 1
637     IF ( pCoordInputData ) coordsSet = coordsSet + 1
638     IF ( rCoordInputData ) coordsSet = coordsSet + 1
639     IF ( coordsSet .GT. 1 ) THEN
640     WRITE(msgBuf,'(A)')
641     & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
642     CALL PRINT_ERROR( msgBuf , myThid)
643     STOP 'ABNORMAL END: S/R INI_PARMS'
644     ENDIF
645 cnh 1.15
646     C-- Input files
647 adcroft 1.41 READ(UNIT=iUnit,NML=PARM05) !,IOSTAT=errIO)
648 cnh 1.35 IF ( errIO .LT. 0 ) THEN
649 cnh 1.15 WRITE(msgBuf,'(A)')
650     & 'S/R INI_PARMS'
651     CALL PRINT_ERROR( msgBuf , 1)
652     WRITE(msgBuf,'(A)')
653     & 'Error reading numerical model '
654     CALL PRINT_ERROR( msgBuf , 1)
655     WRITE(msgBuf,'(A)')
656     & 'parameter file "data"'
657     CALL PRINT_ERROR( msgBuf , 1)
658     WRITE(msgBuf,'(A)')
659     & 'Problem in namelist PARM05'
660     CALL PRINT_ERROR( msgBuf , 1)
661     CALL MODELDATA_EXAMPLE( myThid )
662     STOP 'ABNORMAL END: S/R INI_PARMS'
663 cnh 1.35 ENDIF
664 cnh 1.25
665 cnh 1.28 C
666 cnh 1.30 C-- Set factors required for mixing pressure and meters as vertical coordinate.
667     C rkFac is a "sign" parameter which is used where the orientation of the vertical
668     C coordinate (pressure or meters) relative to the vertical index (K) is important.
669     C rkFac = 1 applies when K and the coordinate are in the opposite sense.
670     C rkFac = -1 applies when K and the coordinate are in the same sense.
671     C horiVertRatio is a parameter that maps horizontal units to vertical units.
672     C It is used in certain special cases where lateral and vertical terms are
673     C being combined and a single frame of reference is needed.
674 adcroft 1.41 IF ( zCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
675 cnh 1.30 rkFac = 1.D0
676     horiVertRatio = 1.D0
677     ENDIF
678 adcroft 1.41 IF ( pCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
679 cnh 1.30 rkFac = -1.D0
680     horiVertRatio = Gravity * rhoConst
681     ENDIF
682 adcroft 1.41 IF ( rCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
683 cnh 1.30 rkFac = 1.D0
684     horiVertRatio = 1.D0
685     ENDIF
686 adcroft 1.41 IF ( rkFac .EQ. UNSET_RS ) rkFac=rkFacDefault
687 cnh 1.28 recip_rkFac = 1.D0 / rkFac
688 cnh 1.32 recip_horiVertRatio = 1./horiVertRatio
689 cnh 1.29 IF ( zCoordInputData ) usingZCoords = .TRUE.
690     IF ( pCoordInputData ) usingPCoords = .TRUE.
691 adcroft 1.37
692     C-- OBCS
693     IF (openBoundaries) THEN
694     READ(UNIT=iUnit,NML=PARM06)
695     DO J=1,Ny
696     if (OB_Ieast(J).lt.0) OB_Ieast(J)=OB_Ieast(J)+Nx+1
697     ENDDO
698     DO I=1,Nx
699     if (OB_Jnorth(I).lt.0) OB_Jnorth(I)=OB_Jnorth(I)+Ny+1
700     ENDDO
701     write(0,*) 'OB Jn =',OB_Jnorth
702     write(0,*) 'OB Js =',OB_Jsouth
703     write(0,*) 'OB Ie =',OB_Ieast
704     write(0,*) 'OB Iw =',OB_Iwest
705     ENDIF
706 adcroft 1.43
707     #ifdef ALLOW_KPP
708     C-- KPP vertical mixing parameters
709     usingKPPmixing = .FALSE.
710     KPPmixingMaps = .FALSE.
711     KPPwriteState = .FALSE.
712     minKPPghat = 0.
713     maxKPPghat = 1.e10
714     minKPPviscAz = 0.
715     DO K =1,Nr
716     maxKPPviscAz(K) = 100.
717     ENDDO
718     minKPPdiffKzT = 0.
719     maxKPPdiffKzT = 100.
720     minKPPdiffKzS = 0.
721     maxKPPdiffKzS = 100.
722    
723     c-----------------------------------------------------------------------
724     c define some non-dimensional constants and
725     c the vertical mixing coefficients in m-k-s units
726     c-----------------------------------------------------------------------
727    
728     epsln = 1.e-20 !
729     epsilon = 0.1 ! nondimensional extent of surface layer
730     vonk = 0.40 ! von Karmans constant
731     conc1 = 5.0 ! scalar coefficients
732     conam = 1.257 ! .
733     concm = 8.380 ! .
734     conc2 = 16.0 ! .
735     zetam = -0.2 ! .
736     conas = -28.86 ! .
737     concs = 98.96 ! .
738     conc3 = 16.0 ! .
739     zetas = -1.0 ! .
740    
741     c parameters for subroutine "bldepth"
742    
743     Ricr = 0.30 ! critical bulk Richardson Number
744     cekman = 0.7 ! coefficient for ekman depth
745     cmonob = 1.0 ! coefficient for Monin-Obukhov depth
746     concv = 1.8 ! ratio of interior buoyancy frequency to
747     ! buoyancy frequency at entrainment depth
748     hbf = 1.0 ! fraction of bounadry layer depth to
749     ! which absorbed solar radiation
750     ! contributes to surface buoyancy forcing
751     Vtc = concv * sqrt(0.2/concs/epsilon) / vonk**2 / Ricr
752     ! non-dimensional coefficient for velocity
753     ! scale of turbulant velocity shear
754    
755     c parameters and common arrays for subroutines "kmixinit" and "wscale"
756    
757     zmin = -4.e-7 ! minimum limit for zehat in table (m3/s3)
758     zmax = 0.0 ! maximum limit for zehat in table (m3/s3)
759     umin = 0.0 ! minimum limit for ustar in table (m/s)
760     umax = .04 ! maximum limit for ustar in table (m/s)
761    
762     c parameters for subroutine "Ri_iwmix"
763    
764     num_v_smooth_Ri = 1 ! number of times Ri is vertically smoothed
765     num_v_smooth_BV = 1 ! number of times BV is vertically smoothed
766     num_z_smooth_sh = 0 ! number of times shear is zonally smoothed
767     num_m_smooth_sh = 0 ! number of times shear is meridionally smoothed
768     Riinfty = 0.7 ! local Richardson Number limit for shear instability
769     BVSQcon = -0.2e-4 ! Brunt-Vaisala squared (1/s^2)
770     difm0 = 0.005 ! max visc due to shear instability (m^2/s)
771     difs0 = 0.005 ! max tracer diffusivity .. (m^2/s)
772     dift0 = 0.005 ! max heat diffusivity .. (m^2/s)
773     difmiw = 0.001 ! viscosity from background internal wave (m^2/s)
774     difsiw = 0.00003 ! tracer diffusivity .. (m^2/s)
775     diftiw = 0.00003 ! heat diffusivity .. (m^2/s)
776     difmcon = 0.1 ! viscosity due to convective instability (m^2/s)
777     difscon = 0.1 ! tracer diffusivity .. (m^2/s)
778     diftcon = 0.1 ! heat diffusivity .. (m^2/s)
779    
780     c parameters for subroutine "blmix"
781    
782     cstar = 10. ! proportionality coefficient for nonlocal transport
783    
784     ! non-dimensional coefficient for counter-gradient term
785     cg = cstar * vonk * (concs * vonk * epsilon)**(1./3.)
786    
787     READ(UNIT=iUnit,NML=PARM07,IOSTAT=errIO,err=13)
788     IF ( errIO .GE. 0 ) GOTO 14
789     13 CONTINUE
790     WRITE(msgBuf,'(A)')
791     & 'S/R INI_PARMS'
792     CALL PRINT_ERROR( msgBuf , 1)
793     WRITE(msgBuf,'(A)')
794     & 'Error reading numerical model '
795     CALL PRINT_ERROR( msgBuf , 1)
796     WRITE(msgBuf,'(A)')
797     & 'parameter file "data"'
798     CALL PRINT_ERROR( msgBuf , 1)
799     WRITE(msgBuf,'(A)')
800     & 'Problem in namelist PARM07'
801     CALL PRINT_ERROR( msgBuf , 1)
802     CALL MODELDATA_EXAMPLE( myThid )
803     STOP 'ABNORMAL END: S/R INI_PARMS'
804     14 CONTINUE
805     #endif
806 adcroft 1.37
807 cnh 1.25 C
808     CLOSE(iUnit)
809 cnh 1.1
810     _END_MASTER(myThid)
811    
812     C-- Everyone else must wait for the parameters to be loaded
813     _BARRIER
814     C
815    
816     RETURN
817     END
818    

  ViewVC Help
Powered by ViewVC 1.1.22