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

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

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


Revision 1.35 - (show annotations) (download)
Mon Nov 2 03:34:12 1998 UTC (25 years, 7 months ago) by cnh
Branch: MAIN
CVS Tags: checkpoint17
Changes since 1.34: +28 -34 lines
Changes for TAMC compatability.
Added exp0 a barotropic basin scale box example
Modified exp1 and exp2 to correct SIZE.h for Nr and
variable overlap width support.

1 C $Header: /u/gcmpack/models/MITgcmUV/model/src/ini_parms.F,v 1.34 1998/10/28 03:11:37 cnh Exp $
2
3 #include "CPP_OPTIONS.h"
4
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
18 C === Global variables ===
19 #include "SIZE.h"
20 #include "EEPARAMS.h"
21 #include "PARAMS.h"
22 #include "GRID.h"
23 #include "CG2D.h"
24
25 C === Routine arguments ===
26 C myThid - Number of this instance of INI_PARMS
27 INTEGER myThid
28
29 C === Local variables ===
30 C dxSpacing, dySpacing - Default spacing in X and Y.
31 C Units are that of coordinate system
32 C i.e. cartesian => metres
33 C s. polar => degrees
34 C goptCount - Used to count the nuber of grid options
35 C (only one is allowed! )
36 C msgBuf - Informational/error meesage buffer
37 C errIO - IO error flag
38 C iUnit - Work variable for IO unit number
39 C record - Work variable for IO buffer
40 C K, I, J - Loop counters
41 C xxxDefault - Default value for variable xxx
42 _RL dxSpacing
43 _RL dySpacing
44 CHARACTER*(MAX_LEN_MBUF) msgBuf
45 CHARACTER*(MAX_LEN_PREC) record
46 INTEGER goptCount
47 INTEGER K, I, J, IL, iUnit
48 INTEGER errIO
49 INTEGER IFNBLNK
50 EXTERNAL IFNBLNK
51 INTEGER ILNBLNK
52 EXTERNAL ILNBLNK
53 C Default values for variables which have vertical coordinate system
54 C dependency.
55 _RL viscArDefault
56 _RL diffKrTDefault
57 _RL diffKrSDefault
58 _RL hFacMinDrDefault
59 _RL delRDefault
60 C zCoordInputData - These are used to select between different coordinate systems.
61 C pCoordInputData The vertical coordinate system in the rest of the model is
62 C rCoordInputData written in terms of r. In the model "data" file input data can
63 C coordsSet be interms of z, p or r.
64 C e.g. delZ or delP or delR for the vertical grid spacing.
65 C The following rules apply:
66 C All parameters must use the same vertical coordinate system.
67 C e.g. delZ and viscAz is legal but
68 C delZ and viscAr is an error.
69 C Similarly specifyinh delZ and delP is an error.
70 C zCoord..., pCoord..., rCoord... are used to flag when z, p or r are
71 C used. coordsSet counts how many vertical coordinate systems have been
72 C used to specify variables. coordsSet > 1 is an error.
73 C
74 LOGICAL zCoordInputData
75 LOGICAL pCoordInputData
76 LOGICAL rCoordInputData
77 INTEGER coordsSet
78
79 C-- Continuous equation parameters
80 NAMELIST /PARM01/
81 & gravity, gBaro, rhonil, tAlpha, sBeta, f0, beta,
82 & viscAh, viscAz, viscA4,
83 & diffKhT, diffKzT, diffK4T,
84 & diffKhS, diffKzS, diffK4S,
85 & GMmaxslope,GMlength,GMalpha,GMdepth,GMkbackground,GMmaxval,
86 & tRef, sRef, eosType,
87 & momViscosity, momAdvection, momForcing, useCoriolis,
88 & momPressureForcing, metricTerms,
89 & tempDiffusion, tempAdvection, tempForcing,
90 & saltDiffusion, saltAdvection, saltForcing,
91 & implicitFreeSurface, rigidLid, freeSurfFac, hFacMin, hFacMinDz,
92 & tempStepping, saltStepping, momStepping, implicitDiffusion,
93 & viscAr, diffKrT, diffKrS, hFacMinDr,
94 & viscAp, diffKpT, diffKpS, hFacMinDp,
95 & rhoConst, buoyancyRelation
96
97 C-- Elliptic solver parameters
98 NAMELIST /PARM02/
99 & cg2dMaxIters, cg2dChkResFreq, cg2dTargetResidual, cg2dpcOffDFac,
100 & cg3dMaxIters, cg3dChkResFreq, cg3dTargetResidual
101
102 C-- Time stepping parammeters
103 NAMELIST /PARM03/
104 & nIter0, nTimeSteps, deltaT, deltaTmom, deltaTtracer, abEps,
105 & tauCD, rCD,
106 & startTime, endTime, chkPtFreq, dumpFreq, taveFreq, deltaTClock,
107 & pChkPtFreq, cAdjFreq, tauThetaClimRelax, tauSaltClimRelax,
108 & periodicExternalForcing, externForcingPeriod, externForcingCycle
109
110 C-- Gridding parameters
111 NAMELIST /PARM04/
112 & usingCartesianGrid, delZ, dxSpacing, dySpacing, delX, delY,
113 & usingSphericalPolarGrid, phiMin, thetaMin, rSphere,
114 & l, m, n, delP, delR, rkFac
115
116 C-- Input files
117 NAMELIST /PARM05/
118 & bathyFile, hydrogThetaFile, hydrogSaltFile,
119 & zonalWindFile, meridWindFile, thetaClimFile,
120 & saltClimFile
121
122
123 C
124 _BEGIN_MASTER(myThid)
125
126 C-- Initialise "which vertical coordinate system used" flags.
127 zCoordInputData = .FALSE.
128 pCoordInputData = .FALSE.
129 rCoordInputData = .FALSE.
130 usingPCoords = .FALSE.
131 usingZCoords = .FALSE.
132 coordsSet = 0
133
134 C-- Open the parameter file
135 OPEN(UNIT=scrUnit1,STATUS='SCRATCH')
136 OPEN(UNIT=scrUnit2,STATUS='SCRATCH')
137 OPEN(UNIT=modelDataUnit,FILE='data',STATUS='OLD',
138 & IOSTAT=errIO)
139 IF ( errIO .LT. 0 ) THEN
140 WRITE(msgBuf,'(A)')
141 & 'S/R INI_PARMS'
142 CALL PRINT_ERROR( msgBuf , 1)
143 WRITE(msgBuf,'(A)')
144 & 'Unable to open model parameter'
145 CALL PRINT_ERROR( msgBuf , 1)
146 WRITE(msgBuf,'(A)')
147 & 'file "data"'
148 CALL PRINT_ERROR( msgBuf , 1)
149 CALL MODELDATA_EXAMPLE( myThid )
150 STOP 'ABNORMAL END: S/R INI_PARMS'
151 ENDIF
152
153 DO WHILE ( .TRUE. )
154 READ(modelDataUnit,FMT='(A)',END=1001) RECORD
155 IL = MAX(ILNBLNK(RECORD),1)
156 IF ( RECORD(1:1) .NE. commentCharacter )
157 & WRITE(UNIT=scrUnit1,FMT='(A)') RECORD(:IL)
158 WRITE(UNIT=scrUnit2,FMT='(A)') RECORD(:IL)
159 ENDDO
160 1001 CONTINUE
161 CLOSE(modelDataUnit)
162
163 C-- Report contents of model parameter file
164 WRITE(msgBuf,'(A)')
165 &'// ======================================================='
166 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
167 & SQUEEZE_RIGHT , 1)
168 WRITE(msgBuf,'(A)') '// Model parameter file "data"'
169 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
170 & SQUEEZE_RIGHT , 1)
171 WRITE(msgBuf,'(A)')
172 &'// ======================================================='
173 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
174 & SQUEEZE_RIGHT , 1)
175 iUnit = scrUnit2
176 REWIND(iUnit)
177 DO WHILE ( .TRUE. )
178 READ(UNIT=iUnit,FMT='(A)',END=2001) RECORD
179 IL = MAX(ILNBLNK(RECORD),1)
180 WRITE(msgBuf,'(A,A)') '>',RECORD(:IL)
181 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
182 & SQUEEZE_RIGHT , 1)
183 ENDDO
184 2001 CONTINUE
185 CLOSE(iUnit)
186 WRITE(msgBuf,'(A)') ' '
187 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
188 & SQUEEZE_RIGHT , 1)
189
190
191 C-- Read settings from model parameter file "data".
192 iUnit = scrUnit1
193 REWIND(iUnit)
194
195 C-- Set default "physical" parameters
196 DO K =1,Nr
197 tRef(K) = 30.D0 - FLOAT( K )
198 ENDDO
199 gravity = 9.81D0
200 gBaro = gravity
201 rhoNil = 999.8D0
202 rhoConst = 999.8D0
203 f0 = 1.D-4
204 beta = 1.D-11
205 viscAh = 1.D3
206 diffKhT = 1.D3
207 diffKhS = 1.D3
208 viscArDefault = 1.D-3
209 viscAz = UNSET_RL
210 viscAr = UNSET_RL
211 viscAp = UNSET_RL
212 diffKrTDefault = 1.D-5
213 diffKzT = UNSET_RL
214 diffKpT = UNSET_RL
215 diffKrT = UNSET_RL
216 diffKrSDefault = 1.D-5
217 diffKzS = UNSET_RL
218 diffKpS = UNSET_RL
219 diffKrS = UNSET_RL
220 viscA4 = 0.
221 diffK4T = 0.
222 diffK4S = 0.
223 GMmaxslope = 1.D-2
224 GMlength = 200.D3
225 GMalpha = 0.D0
226 GMdepth = 1000.D0
227 GMkbackground= 0.D0
228 GMmaxval = 2500.D0
229 tAlpha = 2.D-4
230 sBeta = 7.4D-4
231 eosType = 'LINEAR'
232 buoyancyRelation = 'OCEANIC'
233 implicitFreeSurface = .TRUE.
234 rigidLid = .FALSE.
235 freeSurfFac = 1.D0
236 hFacMin = 0.D0
237 hFacMinDrDefault = 0.D0
238 hFacMinDr = UNSET_RL
239 hFacMinDz = UNSET_RL
240 hFacMinDp = UNSET_RL
241 momViscosity = .TRUE.
242 momAdvection = .TRUE.
243 momForcing = .TRUE.
244 useCoriolis = .TRUE.
245 momPressureForcing = .TRUE.
246 momStepping = .TRUE.
247 tempStepping = .TRUE.
248 saltStepping = .TRUE.
249 metricTerms = .TRUE.
250 implicitDiffusion = .FALSE.
251 READ(UNIT=iUnit,NML=PARM01,IOSTAT=errIO)
252 IF ( errIO .LT. 0 ) THEN
253 WRITE(msgBuf,'(A)')
254 & 'S/R INI_PARMS'
255 CALL PRINT_ERROR( msgBuf , 1)
256 WRITE(msgBuf,'(A)')
257 & 'Error reading numerical model '
258 CALL PRINT_ERROR( msgBuf , 1)
259 WRITE(msgBuf,'(A)')
260 & 'parameter file "data"'
261 CALL PRINT_ERROR( msgBuf , 1)
262 WRITE(msgBuf,'(A)')
263 & 'Problem in namelist PARM01'
264 CALL PRINT_ERROR( msgBuf , 1)
265 CALL MODELDATA_EXAMPLE( myThid )
266 STOP 'ABNORMAL END: S/R INI_PARMS'
267 ENDIF
268 IF ( implicitFreeSurface ) freeSurfFac = 1.D0
269 IF ( rigidLid ) freeSurfFac = 0.D0
270 C-- Momentum viscosity on/off flag.
271 IF ( momViscosity ) THEN
272 vfFacMom = 1.D0
273 ELSE
274 vfFacMom = 0.D0
275 ENDIF
276 C-- Momentum advection on/off flag.
277 IF ( momAdvection ) THEN
278 afFacMom = 1.D0
279 ELSE
280 afFacMom = 0.D0
281 ENDIF
282 C-- Momentum forcing on/off flag.
283 IF ( momForcing ) THEN
284 foFacMom = 1.D0
285 ELSE
286 foFacMom = 0.D0
287 ENDIF
288 C-- Coriolis term on/off flag.
289 IF ( useCoriolis ) THEN
290 cfFacMom = 1.D0
291 ELSE
292 cfFacMom = 0.D0
293 ENDIF
294 C-- Pressure term on/off flag.
295 IF ( momPressureForcing ) THEN
296 pfFacMom = 1.D0
297 ELSE
298 pfFacMom = 0.D0
299 ENDIF
300 C-- Metric terms on/off flag.
301 IF ( metricTerms ) THEN
302 mTFacMom = 1.D0
303 ELSE
304 mTFacMom = 1.D0
305 ENDIF
306 C-- z,p,r coord input switching.
307 IF ( viscAz .NE. UNSET_RL ) zCoordInputData = .TRUE.
308 IF ( viscAp .NE. UNSET_RL ) pCoordInputData = .TRUE.
309 IF ( viscAr .NE. UNSET_RL ) rCoordInputData = .TRUE.
310 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAz
311 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAp
312 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscArDefault
313
314 IF ( diffKzT .NE. UNSET_RL ) zCoordInputData = .TRUE.
315 IF ( diffKpT .NE. UNSET_RL ) pCoordInputData = .TRUE.
316 IF ( diffKrT .NE. UNSET_RL ) rCoordInputData = .TRUE.
317 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKzT
318 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKpT
319 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKrTDefault
320
321 IF ( diffKzS .NE. UNSET_RL ) zCoordInputData = .TRUE.
322 IF ( diffKpS .NE. UNSET_RL ) pCoordInputData = .TRUE.
323 IF ( diffKrS .NE. UNSET_RL ) rCoordInputData = .TRUE.
324 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKzS
325 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKpS
326 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKrSDefault
327
328 IF ( hFacMinDz .NE. UNSET_RL ) zCoordInputData = .TRUE.
329 IF ( hFacMinDp .NE. UNSET_RL ) pCoordInputData = .TRUE.
330 IF ( hFacMinDr .NE. UNSET_RL ) rCoordInputData = .TRUE.
331 IF ( hFacMinDz .EQ. UNSET_RL ) hFacMinDr = hFacMinDz
332 IF ( hFacMinDp .EQ. UNSET_RL ) hFacMinDr = hFacMinDp
333 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDrDefault
334
335 IF ( implicitFreeSurface .AND. rigidLid ) THEN
336 WRITE(msgBuf,'(A,A)')
337 & 'S/R INI_PARMS: Cannot select both implicitFreeSurface',
338 & ' and rigidLid.'
339 CALL PRINT_ERROR( msgBuf , myThid)
340 STOP 'ABNORMAL END: S/R INI_PARMS'
341 ENDIF
342 coordsSet = 0
343 IF ( zCoordInputData ) coordsSet = coordsSet + 1
344 IF ( pCoordInputData ) coordsSet = coordsSet + 1
345 IF ( rCoordInputData ) coordsSet = coordsSet + 1
346 IF ( coordsSet .GT. 1 ) THEN
347 WRITE(msgBuf,'(A)')
348 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
349 CALL PRINT_ERROR( msgBuf , myThid)
350 STOP 'ABNORMAL END: S/R INI_PARMS'
351 ENDIF
352 IF ( rhoConst .LE. 0. ) THEN
353 WRITE(msgBuf,'(A)')
354 & 'S/R INI_PARMS: rhoConst must be greater than 0.'
355 CALL PRINT_ERROR( msgBuf , myThid)
356 STOP 'ABNORMAL END: S/R INI_PARMS'
357 ELSE
358 recip_rhoConst = 1.D0 / rhoConst
359 ENDIF
360 IF ( gravity .LE. 0. ) THEN
361 WRITE(msgBuf,'(A)')
362 & 'S/R INI_PARMS: gravity must be greater than 0.'
363 CALL PRINT_ERROR( msgBuf , myThid)
364 STOP 'ABNORMAL END: S/R INI_PARMS'
365 ELSE
366 recip_gravity = 1.D0 / gravity
367 ENDIF
368
369 C-- Elliptic solver parameters
370 cg2dMaxIters = 150
371 cg2dTargetResidual = 1.D-7
372 cg2dChkResFreq = 1
373 cg3dMaxIters = 150
374 cg3dTargetResidual = 1.D-7
375 cg3dChkResFreq = 1
376 cg2dpcOffDFac = 0.51D0
377 READ(UNIT=iUnit,NML=PARM02,IOSTAT=errIO)
378 IF ( errIO .LT. 0 ) THEN
379 WRITE(msgBuf,'(A)')
380 & 'S/R INI_PARMS'
381 CALL PRINT_ERROR( msgBuf , 1)
382 WRITE(msgBuf,'(A)')
383 & 'Error reading numerical model '
384 CALL PRINT_ERROR( msgBuf , 1)
385 WRITE(msgBuf,'(A)')
386 & 'parameter file "data".'
387 CALL PRINT_ERROR( msgBuf , 1)
388 WRITE(msgBuf,'(A)')
389 & 'Problem in namelist PARM02'
390 CALL PRINT_ERROR( msgBuf , 1)
391 CALL MODELDATA_EXAMPLE( myThid )
392 STOP 'ABNORMAL END: S/R INI_PARMS'
393 ENDIF
394
395 C-- Time stepping parameters
396 startTime = 0.
397 nTimeSteps = 0
398 endTime = 0.
399 nIter0 = 0
400 deltaT = 0.
401 deltaTClock = 0.
402 deltaTtracer = 0.
403 deltaTMom = 0.
404 abEps = 0.01
405 pchkPtFreq = 0.
406 chkPtFreq = 3600.*25
407 dumpFreq = 3600.*100
408 taveFreq = 0.
409 writeStatePrec = precFloat64
410 nCheckLev = 1
411 checkPtSuff(1) = 'ckptA'
412 checkPtSuff(2) = 'ckptB'
413 cAdjFreq = -1.D0
414 rCD = -1.D0
415 tauCD = 0.D0
416 tauThetaClimRelax = 0.D0
417 doThetaClimRelax = .FALSE.
418 tauSaltClimRelax = 0.D0
419 doSaltClimRelax = .FALSE.
420 periodicExternalForcing = .FALSE.
421 externForcingPeriod = 0.
422 externForcingCycle = 0.
423 READ(UNIT=iUnit,NML=PARM03,IOSTAT=errIO)
424 IF ( errIO .LT. 0 ) THEN
425 WRITE(msgBuf,'(A)')
426 & 'S/R INI_PARMS'
427 CALL PRINT_ERROR( msgBuf , 1)
428 WRITE(msgBuf,'(A)')
429 & 'Error reading numerical model '
430 CALL PRINT_ERROR( msgBuf , 1)
431 WRITE(msgBuf,'(A)')
432 & 'parameter file "data"'
433 CALL PRINT_ERROR( msgBuf , 1)
434 WRITE(msgBuf,'(A)')
435 & 'Problem in namelist PARM03'
436 CALL PRINT_ERROR( msgBuf , 1)
437 CALL MODELDATA_EXAMPLE( myThid )
438 STOP 'ABNORMAL END: S/R INI_PARMS'
439 ENDIF
440 C Process "timestepping" params
441 C o Time step size
442 IF ( deltaT .EQ. 0. ) deltaT = deltaTmom
443 IF ( deltaT .EQ. 0. ) deltaT = deltaTtracer
444 IF ( deltaTmom .EQ. 0. ) deltaTmom = deltaT
445 IF ( deltaTtracer .EQ. 0. ) deltaTtracer = deltaT
446 IF ( deltaTClock .EQ. 0. ) deltaTClock = deltaT
447 IF ( periodicExternalForcing ) THEN
448 IF ( externForcingCycle*externForcingPeriod .EQ. 0. ) THEN
449 WRITE(msgBuf,'(A)')
450 & 'S/R INI_PARMS: externForcingCycle,externForcingPeriod =0'
451 CALL PRINT_ERROR( msgBuf , 1)
452 STOP 'ABNORMAL END: S/R INI_PARMS'
453 ENDIF
454 IF ( INT(externForcingCycle/externForcingPeriod) .NE.
455 & externForcingCycle/externForcingPeriod ) THEN
456 WRITE(msgBuf,'(A)')
457 & 'S/R INI_PARMS: externForcingCycle <> N*externForcingPeriod'
458 CALL PRINT_ERROR( msgBuf , 1)
459 STOP 'ABNORMAL END: S/R INI_PARMS'
460 ENDIF
461 IF ( externForcingCycle.le.externForcingPeriod ) THEN
462 WRITE(msgBuf,'(A)')
463 & 'S/R INI_PARMS: externForcingCycle < externForcingPeriod'
464 CALL PRINT_ERROR( msgBuf , 1)
465 STOP 'ABNORMAL END: S/R INI_PARMS'
466 ENDIF
467 IF ( externForcingPeriod.lt.deltaTclock ) THEN
468 WRITE(msgBuf,'(A)')
469 & 'S/R INI_PARMS: externForcingPeriod < deltaTclock'
470 CALL PRINT_ERROR( msgBuf , 1)
471 STOP 'ABNORMAL END: S/R INI_PARMS'
472 ENDIF
473 ENDIF
474 C o Convection frequency
475 IF ( cAdjFreq .LT. 0. ) THEN
476 cAdjFreq = deltaTClock
477 ENDIF
478 C o CD coupling
479 IF ( tauCD .EQ. 0.D0 ) THEN
480 tauCD = deltaTmom
481 ENDIF
482 IF ( rCD .LT. 0. ) THEN
483 rCD = 1. - deltaTMom/tauCD
484 ENDIF
485 C o Temperature climatology relaxation time scale
486 IF ( tauThetaClimRelax .EQ. 0.D0 ) THEN
487 doThetaClimRelax = .FALSE.
488 lambdaThetaClimRelax = 0.D0
489 ELSE
490 doThetaClimRelax = .TRUE.
491 lambdaThetaClimRelax = 1./tauThetaClimRelax
492 ENDIF
493 C o Salinity climatology relaxation time scale
494 IF ( tauSaltClimRelax .EQ. 0.D0 ) THEN
495 doSaltClimRelax = .FALSE.
496 lambdaSaltClimRelax = 0.D0
497 ELSE
498 doSaltClimRelax = .TRUE.
499 lambdaSaltClimRelax = 1./tauSaltClimRelax
500 ENDIF
501 C o Time step count
502 IF ( endTime .NE. 0 ) THEN
503 IF ( deltaTClock .NE. 0 ) nTimeSteps =
504 & INT((endTime-startTime)/deltaTClock)
505 ENDIF
506 C o Finish time
507 IF ( endTime .EQ. 0. ) endTime = FLOAT(nTimeSteps)*deltaTClock
508
509 C o If taveFreq is finite, then we must make sure the diagnostics
510 C code is being compiled
511 #ifndef ALLOW_DIAGNOSTICS
512 IF (taveFreq.NE.0.) THEN
513 WRITE(msgBuf,'(A)')
514 & 'S/R INI_PARMS: taveFreq <> 0 but you have'
515 CALL PRINT_ERROR( msgBuf , 1)
516 WRITE(msgBuf,'(A)')
517 & 'not compiled the model with the diagnostics routines.'
518 CALL PRINT_ERROR( msgBuf , 1)
519 WRITE(msgBuf,'(A)')
520 & 'Re-compile with: #define ALLOW_DIAGNOSTICS or -DALLOW_DIAGNOSTICS'
521 CALL PRINT_ERROR( msgBuf , 1)
522 STOP 'ABNORMAL END: S/R INI_PARMS'
523 ENDIF
524 #endif
525
526 C-- Grid parameters
527 C In cartesian coords distances are in metres
528 rkFac = UNSET_I
529 usingCartesianGrid = .TRUE.
530 delRDefault = 1.D2
531 DO K =1,Nr
532 delZ(K) = UNSET_RL
533 delP(K) = UNSET_RL
534 delR(K) = UNSET_RL
535 ENDDO
536 dxSpacing = 20.D0 * 1000.D0
537 dySpacing = 20.D0 * 1000.D0
538 DO i=1,Nx
539 delX(i) = dxSpacing
540 ENDDO
541 DO j=1,Ny
542 delY(j) = dySpacing
543 ENDDO
544 C In spherical polar distances are in degrees
545 usingSphericalPolarGrid = .FALSE.
546 phiMin = -5.0
547 thetaMin = 0.
548 rSphere = 6370. * 1.D3
549 recip_rSphere = 1.D0/rSphere
550 IF ( usingSphericalPolarGrid ) THEN
551 dxSpacing = 1.
552 dySpacing = 1.
553 DO I=1,Nx
554 delX(I) = dxSpacing
555 ENDDO
556 DO J=1,Ny
557 delY(J) = dySpacing
558 ENDDO
559 ENDIF
560
561 READ(UNIT=iUnit,NML=PARM04,IOSTAT=errIO)
562 IF ( errIO .LT. 0 ) THEN
563 WRITE(msgBuf,'(A)')
564 & 'S/R INI_PARMS'
565 CALL PRINT_ERROR( msgBuf , 1)
566 WRITE(msgBuf,'(A)')
567 & 'Error reading numerical model '
568 CALL PRINT_ERROR( msgBuf , 1)
569 WRITE(msgBuf,'(A)')
570 & 'parameter file "data"'
571 CALL PRINT_ERROR( msgBuf , 1)
572 WRITE(msgBuf,'(A)')
573 & 'Problem in namelist PARM04'
574 CALL PRINT_ERROR( msgBuf , 1)
575 CALL MODELDATA_EXAMPLE( myThid )
576 STOP 'ABNORMAL END: S/R INI_PARMS'
577 ENDIF
578 C
579 IF ( rSphere .NE. 0 ) THEN
580 recip_rSphere = 1.D0/rSphere
581 ELSE
582 recip_rSphere = 0.
583 ENDIF
584 C-- Initialize EOS coefficients (3rd order polynomial)
585 IF (eostype.eq.'POLY3') THEN
586 OPEN(37,FILE='POLY3.COEFFS',STATUS='OLD',FORM='FORMATTED')
587 READ(37,*) I
588 IF (I.NE.Nr) THEN
589 WRITE(msgBuf,'(A)')
590 & 'ini_parms: attempt to read POLY3.COEFFS failed'
591 CALL PRINT_ERROR( msgBuf , 1)
592 WRITE(msgBuf,'(A)')
593 & ' because bad # of levels in data'
594 CALL PRINT_ERROR( msgBuf , 1)
595 STOP 'Bad data in POLY3.COEFFS'
596 ENDIF
597 READ(37,*) (eosRefT(K),eosRefS(K),eosSig0(K),K=1,Nr)
598 DO K=1,Nr
599 READ(37,*) (eosC(I,K),I=1,9)
600 ENDDO
601 CLOSE(37)
602 ENDIF
603 C-- Check for conflicting grid definitions.
604 goptCount = 0
605 IF ( usingCartesianGrid ) goptCount = goptCount+1
606 IF ( usingSphericalPolarGrid ) goptCount = goptCount+1
607 IF ( goptCount .NE. 1 ) THEN
608 WRITE(msgBuf,'(A)')
609 & 'S/R INI_PARMS: More than one coordinate system requested'
610 CALL PRINT_ERROR( msgBuf , myThid)
611 STOP 'ABNORMAL END: S/R INI_PARMS'
612 ENDIF
613 C-- Make metric term settings consistent with underlying grid.
614 IF ( usingCartesianGrid ) THEN
615 usingSphericalPolarMterms = .FALSE.
616 metricTerms = .FALSE.
617 mTFacMom = 0
618 useBetaPlaneF = .TRUE.
619 ENDIF
620 IF ( usingSphericalPolarGrid ) THEN
621 useConstantF = .FALSE.
622 useBetaPlaneF = .FALSE.
623 useSphereF = .TRUE.
624 omega = 2.D0 * PI / ( 3600.D0 * 24.D0 )
625 usingSphericalPolarMterms = .TRUE.
626 metricTerms = .TRUE.
627 mTFacMom = 1
628 ENDIF
629 C-- p, z, r coord parameters
630 DO K = 1, Nr
631 IF ( delZ(K) .NE. UNSET_RL ) zCoordInputData = .TRUE.
632 IF ( delP(K) .NE. UNSET_RL ) pCoordInputData = .TRUE.
633 IF ( delR(K) .NE. UNSET_RL ) rCoordInputData = .TRUE.
634 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delZ(K)
635 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delP(K)
636 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delRDefault
637 ENDDO
638 C Check for multiple coordinate systems
639 coordsSet = 0
640 IF ( zCoordInputData ) coordsSet = coordsSet + 1
641 IF ( pCoordInputData ) coordsSet = coordsSet + 1
642 IF ( rCoordInputData ) coordsSet = coordsSet + 1
643 IF ( coordsSet .GT. 1 ) THEN
644 WRITE(msgBuf,'(A)')
645 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
646 CALL PRINT_ERROR( msgBuf , myThid)
647 STOP 'ABNORMAL END: S/R INI_PARMS'
648 ENDIF
649
650 C-- Input files
651 bathyFile = ' '
652 hydrogSaltFile = ' '
653 hydrogThetaFile = ' '
654 zonalWindFile = ' '
655 meridWindFile = ' '
656 thetaClimFile = ' '
657 saltClimFile = ' '
658 READ(UNIT=iUnit,NML=PARM05,IOSTAT=errIO)
659 IF ( errIO .LT. 0 ) THEN
660 WRITE(msgBuf,'(A)')
661 & 'S/R INI_PARMS'
662 CALL PRINT_ERROR( msgBuf , 1)
663 WRITE(msgBuf,'(A)')
664 & 'Error reading numerical model '
665 CALL PRINT_ERROR( msgBuf , 1)
666 WRITE(msgBuf,'(A)')
667 & 'parameter file "data"'
668 CALL PRINT_ERROR( msgBuf , 1)
669 WRITE(msgBuf,'(A)')
670 & 'Problem in namelist PARM05'
671 CALL PRINT_ERROR( msgBuf , 1)
672 CALL MODELDATA_EXAMPLE( myThid )
673 STOP 'ABNORMAL END: S/R INI_PARMS'
674 ENDIF
675
676 C
677 C-- Set factors required for mixing pressure and meters as vertical coordinate.
678 C rkFac is a "sign" parameter which is used where the orientation of the vertical
679 C coordinate (pressure or meters) relative to the vertical index (K) is important.
680 C rkFac = 1 applies when K and the coordinate are in the opposite sense.
681 C rkFac = -1 applies when K and the coordinate are in the same sense.
682 C horiVertRatio is a parameter that maps horizontal units to vertical units.
683 C It is used in certain special cases where lateral and vertical terms are
684 C being combined and a single frame of reference is needed.
685 IF ( zCoordInputData .AND. rkFac .EQ. UNSET_I ) THEN
686 rkFac = 1.D0
687 horiVertRatio = 1.D0
688 ENDIF
689 IF ( pCoordInputData .AND. rkFac .EQ. UNSET_I ) THEN
690 rkFac = -1.D0
691 horiVertRatio = Gravity * rhoConst
692 ENDIF
693 IF ( rCoordInputData .AND. rkFac .EQ. UNSET_I ) THEN
694 rkFac = 1.D0
695 horiVertRatio = 1.D0
696 ENDIF
697 recip_rkFac = 1.D0 / rkFac
698 recip_horiVertRatio = 1./horiVertRatio
699 IF ( zCoordInputData ) usingZCoords = .TRUE.
700 IF ( pCoordInputData ) usingPCoords = .TRUE.
701 C
702 CLOSE(iUnit)
703
704 _END_MASTER(myThid)
705
706 C-- Everyone else must wait for the parameters to be loaded
707 _BARRIER
708 C
709
710 RETURN
711 END
712

  ViewVC Help
Powered by ViewVC 1.1.22