/[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.54 - (show annotations) (download)
Sun Feb 4 14:38:47 2001 UTC (23 years, 3 months ago) by cnh
Branch: MAIN
CVS Tags: checkpoint35
Changes since 1.53: +2 -1 lines
Made sure each .F and .h file had
the CVS keywords Header and Name at its start.
Most had header but very few currently have Name, so
lots of changes!

1 C $Header: /u/gcmpack/models/MITgcmUV/model/src/ini_parms.F,v 1.53 2001/02/02 21:04:48 adcroft Exp $
2 C $Name: $
3
4 #include "CPP_OPTIONS.h"
5
6 SUBROUTINE INI_PARMS( myThid )
7 C /==========================================================\
8 C | SUBROUTINE INI_PARMS |
9 C | o Routine to set model "parameters" |
10 C |==========================================================|
11 C | Notes: |
12 C | ====== |
13 C | The present version of this routine is a place-holder. |
14 C | A production version needs to handle parameters from an |
15 C | external file and possibly reading in some initial field |
16 C | values. |
17 C \==========================================================/
18 IMPLICIT NONE
19
20 C === Global variables ===
21 #include "SIZE.h"
22 #include "EEPARAMS.h"
23 #include "PARAMS.h"
24 #include "GRID.h"
25 #include "CG2D.h"
26
27 C === Routine arguments ===
28 C myThid - Number of this instance of INI_PARMS
29 INTEGER myThid
30
31 C === Local variables ===
32 C dxSpacing, dySpacing - Default spacing in X and Y.
33 C Units are that of coordinate system
34 C i.e. cartesian => metres
35 C s. polar => degrees
36 C goptCount - Used to count the nuber of grid options
37 C (only one is allowed! )
38 C msgBuf - Informational/error meesage buffer
39 C errIO - IO error flag
40 C iUnit - Work variable for IO unit number
41 C record - Work variable for IO buffer
42 C K, I, J - Loop counters
43 C xxxDefault - Default value for variable xxx
44 _RL dxSpacing
45 _RL dySpacing
46 CHARACTER*(MAX_LEN_FNAM) delXfile
47 CHARACTER*(MAX_LEN_FNAM) delYfile
48 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 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 _RL delRDefault(Nr)
64 _RS rkFacDefault
65 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
84 C-- Continuous equation parameters
85 NAMELIST /PARM01/
86 & gravity, gBaro, rhonil, tAlpha, sBeta, f0, beta,
87 & viscAh, viscAz, viscA4, cosPower,
88 & diffKhT, diffKzT, diffK4T,
89 & diffKhS, diffKzS, diffK4S,
90 & tRef, sRef, eosType,
91 & no_slip_sides,no_slip_bottom,
92 & momViscosity, momAdvection, momForcing, useCoriolis,
93 & momPressureForcing, metricTerms,
94 & tempDiffusion, tempAdvection, tempForcing,
95 & saltDiffusion, saltAdvection, saltForcing,
96 & implicitFreeSurface, rigidLid, freeSurfFac, hFacMin, hFacMinDz,
97 & staggerTimeStep,
98 & tempStepping, saltStepping, momStepping,
99 & implicitDiffusion, implicitViscosity,
100 & viscAr, diffKrT, diffKrS, hFacMinDr,
101 & viscAp, diffKpT, diffKpS, hFacMinDp,
102 & rhoConst, buoyancyRelation, HeatCapacity_Cp,
103 & writeBinaryPrec, readBinaryPrec, writeStatePrec,
104 & nonHydrostatic, globalFiles,
105 & allowFreezing, ivdc_kappa,
106 & nShap, zonal_filt_lat, zonal_filt_sinpow, zonal_filt_cospow,
107 & bottomDragLinear,bottomDragQuadratic
108
109 C-- Elliptic solver parameters
110 NAMELIST /PARM02/
111 & cg2dMaxIters, cg2dChkResFreq, cg2dTargetResidual, cg2dpcOffDFac,
112 & cg3dMaxIters, cg3dChkResFreq, cg3dTargetResidual
113
114 C-- Time stepping parammeters
115 NAMELIST /PARM03/
116 & nIter0, nTimeSteps, nEndIter, deltaT, deltaTmom, deltaTtracer,
117 & abEps, tauCD, rCD,
118 & startTime, endTime, chkPtFreq, dumpFreq, taveFreq, deltaTClock,
119 & pChkPtFreq, cAdjFreq, tauThetaClimRelax, tauSaltClimRelax,
120 & periodicExternalForcing, externForcingPeriod, externForcingCycle
121
122 C-- Gridding parameters
123 NAMELIST /PARM04/
124 & usingCartesianGrid, dxSpacing, dySpacing, delX, delY, delZ,
125 & usingSphericalPolarGrid, phiMin, thetaMin, rSphere,
126 & delP, delR, rkFac, Ro_SeaLevel, groundAtK1,
127 & delXfile, delYfile
128
129 C-- Input files
130 NAMELIST /PARM05/
131 & bathyFile, hydrogThetaFile, hydrogSaltFile,
132 & zonalWindFile, meridWindFile,
133 & thetaClimFile, saltClimFile,
134 & surfQfile, EmPmRfile, surfQswfile,
135 & uVelInitFile, vVelInitFile, pSurfInitFile
136
137 C
138 _BEGIN_MASTER(myThid)
139
140 C Defaults values for input parameters
141 CALL SET_DEFAULTS(
142 O viscArDefault, diffKrTDefault, diffKrSDefault,
143 O hFacMinDrDefault, delRdefault, rkFacDefault,
144 I myThid )
145
146 C-- Initialise "which vertical coordinate system used" flags.
147 zCoordInputData = .FALSE.
148 pCoordInputData = .FALSE.
149 rCoordInputData = .FALSE.
150 usingPCoords = .FALSE.
151 usingZCoords = .FALSE.
152 coordsSet = 0
153
154 C-- Open the parameter file
155 OPEN(UNIT=scrUnit1,STATUS='SCRATCH')
156 OPEN(UNIT=scrUnit2,STATUS='SCRATCH')
157 OPEN(UNIT=modelDataUnit,FILE='data',STATUS='OLD',
158 & IOSTAT=errIO)
159 IF ( errIO .LT. 0 ) THEN
160 WRITE(msgBuf,'(A)')
161 & 'S/R INI_PARMS'
162 CALL PRINT_ERROR( msgBuf , 1)
163 WRITE(msgBuf,'(A)')
164 & 'Unable to open model parameter'
165 CALL PRINT_ERROR( msgBuf , 1)
166 WRITE(msgBuf,'(A)')
167 & 'file "data"'
168 CALL PRINT_ERROR( msgBuf , 1)
169 CALL MODELDATA_EXAMPLE( myThid )
170 STOP 'ABNORMAL END: S/R INI_PARMS'
171 ENDIF
172
173 DO WHILE ( .TRUE. )
174 READ(modelDataUnit,FMT='(A)',END=1001) RECORD
175 IL = MAX(ILNBLNK(RECORD),1)
176 IF ( RECORD(1:1) .NE. commentCharacter )
177 & WRITE(UNIT=scrUnit1,FMT='(A)') RECORD(:IL)
178 WRITE(UNIT=scrUnit2,FMT='(A)') RECORD(:IL)
179 ENDDO
180 1001 CONTINUE
181 CLOSE(modelDataUnit)
182
183 C-- Report contents of model parameter file
184 WRITE(msgBuf,'(A)')
185 &'// ======================================================='
186 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
187 & SQUEEZE_RIGHT , 1)
188 WRITE(msgBuf,'(A)') '// Model parameter file "data"'
189 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
190 & SQUEEZE_RIGHT , 1)
191 WRITE(msgBuf,'(A)')
192 &'// ======================================================='
193 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
194 & SQUEEZE_RIGHT , 1)
195 iUnit = scrUnit2
196 REWIND(iUnit)
197 DO WHILE ( .TRUE. )
198 READ(UNIT=iUnit,FMT='(A)',END=2001) RECORD
199 IL = MAX(ILNBLNK(RECORD),1)
200 WRITE(msgBuf,'(A,A)') '>',RECORD(:IL)
201 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
202 & SQUEEZE_RIGHT , 1)
203 ENDDO
204 2001 CONTINUE
205 CLOSE(iUnit)
206 WRITE(msgBuf,'(A)') ' '
207 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
208 & SQUEEZE_RIGHT , 1)
209
210
211 C-- Read settings from model parameter file "data".
212 iUnit = scrUnit1
213 REWIND(iUnit)
214
215 C-- Set default "physical" parameters
216 viscAz = UNSET_RL
217 viscAr = UNSET_RL
218 viscAp = UNSET_RL
219 diffKzT = UNSET_RL
220 diffKpT = UNSET_RL
221 diffKrT = UNSET_RL
222 diffKzS = UNSET_RL
223 diffKpS = UNSET_RL
224 diffKrS = UNSET_RL
225 gBaro = UNSET_RL
226 rhoConst = UNSET_RL
227 hFacMinDr = UNSET_RL
228 hFacMinDz = UNSET_RL
229 hFacMinDp = UNSET_RL
230 READ(UNIT=iUnit,NML=PARM01) !,IOSTAT=errIO)
231 IF ( errIO .LT. 0 ) THEN
232 WRITE(msgBuf,'(A)')
233 & 'S/R INI_PARMS'
234 CALL PRINT_ERROR( msgBuf , 1)
235 WRITE(msgBuf,'(A)')
236 & 'Error reading numerical model '
237 CALL PRINT_ERROR( msgBuf , 1)
238 WRITE(msgBuf,'(A)')
239 & 'parameter file "data"'
240 CALL PRINT_ERROR( msgBuf , 1)
241 WRITE(msgBuf,'(A)')
242 & 'Problem in namelist PARM01'
243 CALL PRINT_ERROR( msgBuf , 1)
244 CALL MODELDATA_EXAMPLE( myThid )
245 STOP 'ABNORMAL END: S/R INI_PARMS'
246 ENDIF
247 IF ( implicitFreeSurface ) freeSurfFac = 1.D0
248 IF ( rigidLid ) freeSurfFac = 0.D0
249 IF ( gBaro .EQ. UNSET_RL ) gBaro=gravity
250 IF ( rhoConst .EQ. UNSET_RL ) rhoConst=rhoNil
251 C-- Momentum viscosity on/off flag.
252 IF ( momViscosity ) THEN
253 vfFacMom = 1.D0
254 ELSE
255 vfFacMom = 0.D0
256 ENDIF
257 C-- Momentum advection on/off flag.
258 IF ( momAdvection ) THEN
259 afFacMom = 1.D0
260 ELSE
261 afFacMom = 0.D0
262 ENDIF
263 C-- Momentum forcing on/off flag.
264 IF ( momForcing ) THEN
265 foFacMom = 1.D0
266 ELSE
267 foFacMom = 0.D0
268 ENDIF
269 C-- Coriolis term on/off flag.
270 IF ( useCoriolis ) THEN
271 cfFacMom = 1.D0
272 ELSE
273 cfFacMom = 0.D0
274 ENDIF
275 C-- Pressure term on/off flag.
276 IF ( momPressureForcing ) THEN
277 pfFacMom = 1.D0
278 ELSE
279 pfFacMom = 0.D0
280 ENDIF
281 C-- Metric terms on/off flag.
282 IF ( metricTerms ) THEN
283 mTFacMom = 1.D0
284 ELSE
285 mTFacMom = 1.D0
286 ENDIF
287 C-- z,p,r coord input switching.
288 IF ( viscAz .NE. UNSET_RL ) zCoordInputData = .TRUE.
289 IF ( viscAp .NE. UNSET_RL ) pCoordInputData = .TRUE.
290 IF ( viscAr .NE. UNSET_RL ) rCoordInputData = .TRUE.
291 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAz
292 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAp
293 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscArDefault
294
295 IF ( diffKzT .NE. UNSET_RL ) zCoordInputData = .TRUE.
296 IF ( diffKpT .NE. UNSET_RL ) pCoordInputData = .TRUE.
297 IF ( diffKrT .NE. UNSET_RL ) rCoordInputData = .TRUE.
298 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKzT
299 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKpT
300 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKrTDefault
301
302 IF ( diffKzS .NE. UNSET_RL ) zCoordInputData = .TRUE.
303 IF ( diffKpS .NE. UNSET_RL ) pCoordInputData = .TRUE.
304 IF ( diffKrS .NE. UNSET_RL ) rCoordInputData = .TRUE.
305 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKzS
306 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKpS
307 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKrSDefault
308
309 IF ( hFacMinDz .NE. UNSET_RL ) zCoordInputData = .TRUE.
310 IF ( hFacMinDp .NE. UNSET_RL ) pCoordInputData = .TRUE.
311 IF ( hFacMinDr .NE. UNSET_RL ) rCoordInputData = .TRUE.
312 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDz
313 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDp
314 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDrDefault
315
316 IF ( ivdc_kappa .NE. 0. .AND. .NOT. implicitDiffusion ) THEN
317 WRITE(msgBuf,'(A,A)')
318 & 'S/R INI_PARMS: To use ivdc_kappa you must enable implicit',
319 & ' vertical diffusion.'
320 CALL PRINT_ERROR( msgBuf , myThid)
321 STOP 'ABNORMAL END: S/R INI_PARMS'
322 ENDIF
323
324 IF ( implicitFreeSurface .AND. rigidLid ) THEN
325 WRITE(msgBuf,'(A,A)')
326 & 'S/R INI_PARMS: Cannot select both implicitFreeSurface',
327 & ' and rigidLid.'
328 CALL PRINT_ERROR( msgBuf , myThid)
329 STOP 'ABNORMAL END: S/R INI_PARMS'
330 ENDIF
331 coordsSet = 0
332 IF ( zCoordInputData ) coordsSet = coordsSet + 1
333 IF ( pCoordInputData ) coordsSet = coordsSet + 1
334 IF ( rCoordInputData ) coordsSet = coordsSet + 1
335 IF ( coordsSet .GT. 1 ) THEN
336 WRITE(msgBuf,'(A)')
337 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
338 CALL PRINT_ERROR( msgBuf , myThid)
339 STOP 'ABNORMAL END: S/R INI_PARMS'
340 ENDIF
341 IF ( rhoConst .LE. 0. ) THEN
342 WRITE(msgBuf,'(A)')
343 & 'S/R INI_PARMS: rhoConst must be greater than 0.'
344 CALL PRINT_ERROR( msgBuf , myThid)
345 STOP 'ABNORMAL END: S/R INI_PARMS'
346 ELSE
347 recip_rhoConst = 1.D0 / rhoConst
348 ENDIF
349 IF ( rhoNil .LE. 0. ) THEN
350 WRITE(msgBuf,'(A)')
351 & 'S/R INI_PARMS: rhoNil must be greater than 0.'
352 CALL PRINT_ERROR( msgBuf , myThid)
353 STOP 'ABNORMAL END: S/R INI_PARMS'
354 ELSE
355 recip_rhoNil = 1.D0 / rhoNil
356 ENDIF
357 IF ( HeatCapacity_Cp .LE. 0. ) THEN
358 WRITE(msgBuf,'(A)')
359 & 'S/R INI_PARMS: HeatCapacity_Cp must be greater than 0.'
360 CALL PRINT_ERROR( msgBuf , myThid)
361 STOP 'ABNORMAL END: S/R INI_PARMS'
362 ELSE
363 recip_Cp = 1.D0 / HeatCapacity_Cp
364 ENDIF
365 IF ( gravity .LE. 0. ) THEN
366 WRITE(msgBuf,'(A)')
367 & 'S/R INI_PARMS: gravity must be greater than 0.'
368 CALL PRINT_ERROR( msgBuf , myThid)
369 STOP 'ABNORMAL END: S/R INI_PARMS'
370 ELSE
371 recip_gravity = 1.D0 / gravity
372 ENDIF
373 C Set globalFiles flag for READ_WRITE_FLD package
374 CALL SET_WRITE_GLOBAL_FLD( globalFiles )
375 C Set globalFiles flag for READ_WRITE_REC package
376 CALL SET_WRITE_GLOBAL_REC( globalFiles )
377 C Set globalFiles flag for READ_WRITE_REC package
378 CALL SET_WRITE_GLOBAL_PICKUP( globalFiles )
379
380 C-- Elliptic solver parameters
381 READ(UNIT=iUnit,NML=PARM02) !,IOSTAT=errIO)
382 IF ( errIO .LT. 0 ) THEN
383 WRITE(msgBuf,'(A)')
384 & 'S/R INI_PARMS'
385 CALL PRINT_ERROR( msgBuf , 1)
386 WRITE(msgBuf,'(A)')
387 & 'Error reading numerical model '
388 CALL PRINT_ERROR( msgBuf , 1)
389 WRITE(msgBuf,'(A)')
390 & 'parameter file "data".'
391 CALL PRINT_ERROR( msgBuf , 1)
392 WRITE(msgBuf,'(A)')
393 & 'Problem in namelist PARM02'
394 CALL PRINT_ERROR( msgBuf , 1)
395 CALL MODELDATA_EXAMPLE( myThid )
396 STOP 'ABNORMAL END: S/R INI_PARMS'
397 ENDIF
398
399 C-- Time stepping parameters
400 rCD = -1.D0
401 READ(UNIT=iUnit,NML=PARM03) !,IOSTAT=errIO)
402 IF ( errIO .LT. 0 ) THEN
403 WRITE(msgBuf,'(A)')
404 & 'S/R INI_PARMS'
405 CALL PRINT_ERROR( msgBuf , 1)
406 WRITE(msgBuf,'(A)')
407 & 'Error reading numerical model '
408 CALL PRINT_ERROR( msgBuf , 1)
409 WRITE(msgBuf,'(A)')
410 & 'parameter file "data"'
411 CALL PRINT_ERROR( msgBuf , 1)
412 WRITE(msgBuf,'(A)')
413 & 'Problem in namelist PARM03'
414 CALL PRINT_ERROR( msgBuf , 1)
415 CALL MODELDATA_EXAMPLE( myThid )
416 STOP 'ABNORMAL END: S/R INI_PARMS'
417 ENDIF
418 C Process "timestepping" params
419 C o Time step size
420 IF ( deltaT .EQ. 0. ) deltaT = deltaTmom
421 IF ( deltaT .EQ. 0. ) deltaT = deltaTtracer
422 IF ( deltaTmom .EQ. 0. ) deltaTmom = deltaT
423 IF ( deltaTtracer .EQ. 0. ) deltaTtracer = deltaT
424 IF ( deltaTClock .EQ. 0. ) deltaTClock = deltaT
425 IF ( periodicExternalForcing ) THEN
426 IF ( externForcingCycle*externForcingPeriod .EQ. 0. ) THEN
427 WRITE(msgBuf,'(A)')
428 & 'S/R INI_PARMS: externForcingCycle,externForcingPeriod =0'
429 CALL PRINT_ERROR( msgBuf , 1)
430 STOP 'ABNORMAL END: S/R INI_PARMS'
431 ENDIF
432 IF ( INT(externForcingCycle/externForcingPeriod) .NE.
433 & externForcingCycle/externForcingPeriod ) THEN
434 WRITE(msgBuf,'(A)')
435 & 'S/R INI_PARMS: externForcingCycle <> N*externForcingPeriod'
436 CALL PRINT_ERROR( msgBuf , 1)
437 STOP 'ABNORMAL END: S/R INI_PARMS'
438 ENDIF
439 IF ( externForcingCycle.le.externForcingPeriod ) THEN
440 WRITE(msgBuf,'(A)')
441 & 'S/R INI_PARMS: externForcingCycle < externForcingPeriod'
442 CALL PRINT_ERROR( msgBuf , 1)
443 STOP 'ABNORMAL END: S/R INI_PARMS'
444 ENDIF
445 IF ( externForcingPeriod.lt.deltaTclock ) THEN
446 WRITE(msgBuf,'(A)')
447 & 'S/R INI_PARMS: externForcingPeriod < deltaTclock'
448 CALL PRINT_ERROR( msgBuf , 1)
449 STOP 'ABNORMAL END: S/R INI_PARMS'
450 ENDIF
451 ENDIF
452 C o Convection frequency
453 IF ( cAdjFreq .LT. 0. ) THEN
454 cAdjFreq = deltaTClock
455 ENDIF
456 IF ( ivdc_kappa .NE. 0. .AND. cAdjFreq .NE. 0. ) THEN
457 WRITE(msgBuf,'(A,A)')
458 & 'S/R INI_PARMS: You have enabled both ivdc_kappa and',
459 & ' convective adjustment.'
460 CALL PRINT_ERROR( msgBuf , myThid)
461 STOP 'ABNORMAL END: S/R INI_PARMS'
462 ENDIF
463 C o CD coupling
464 IF ( tauCD .EQ. 0.D0 ) THEN
465 tauCD = deltaTmom
466 ENDIF
467 IF ( rCD .LT. 0. ) THEN
468 rCD = 1. - deltaTMom/tauCD
469 ENDIF
470 C o Temperature climatology relaxation time scale
471 IF ( tauThetaClimRelax .EQ. 0.D0 ) THEN
472 doThetaClimRelax = .FALSE.
473 lambdaThetaClimRelax = 0.D0
474 ELSE
475 doThetaClimRelax = .TRUE.
476 lambdaThetaClimRelax = 1./tauThetaClimRelax
477 ENDIF
478 C o Salinity climatology relaxation time scale
479 IF ( tauSaltClimRelax .EQ. 0.D0 ) THEN
480 doSaltClimRelax = .FALSE.
481 lambdaSaltClimRelax = 0.D0
482 ELSE
483 doSaltClimRelax = .TRUE.
484 lambdaSaltClimRelax = 1./tauSaltClimRelax
485 ENDIF
486
487 C o Start time
488 IF ( nIter0 .NE. 0 .AND. startTime .EQ. 0. )
489 & startTime = deltaTClock*float(nIter0)
490 C o nIter0
491 IF ( nIter0 .EQ. 0 .AND. startTime .NE. 0. )
492 & nIter0 = INT( startTime/deltaTClock )
493
494 C o nTimeSteps 1
495 IF ( nTimeSteps .EQ. 0 .AND. nEndIter .NE. 0 )
496 & nTimeSteps = nEndIter-nIter0
497 C o nTimeSteps 2
498 IF ( nTimeSteps .EQ. 0 .AND. endTime .NE. 0. )
499 & nTimeSteps = int(0.5+(endTime-startTime)/deltaTclock)
500 C o nEndIter 1
501 IF ( nEndIter .EQ. 0 .AND. nTimeSteps .NE. 0 )
502 & nEndIter = nIter0+nTimeSteps
503 C o nEndIter 2
504 IF ( nEndIter .EQ. 0 .AND. endTime .NE. 0. )
505 & nEndIter = int(0.5+endTime/deltaTclock)
506 C o End Time 1
507 IF ( endTime .EQ. 0. .AND. nTimeSteps .NE. 0 )
508 & endTime = startTime + deltaTClock*float(nTimeSteps)
509 C o End Time 2
510 IF ( endTime .EQ. 0. .AND. nEndIter .NE. 0 )
511 & endTime = deltaTClock*float(nEndIter)
512
513 C o Consistent?
514 IF ( nEndIter .NE. nIter0+nTimeSteps ) THEN
515 WRITE(msgBuf,'(A)')
516 & 'S/R INI_PARMS: nIter0, nTimeSteps and nEndIter are inconsistent'
517 CALL PRINT_ERROR( msgBuf , 1)
518 WRITE(msgBuf,'(A)')
519 & 'S/R INI_PARMS: Perhaps more than two were set at once'
520 CALL PRINT_ERROR( msgBuf , 1)
521 STOP 'ABNORMAL END: S/R INI_PARMS'
522 ENDIF
523 IF ( nTimeSteps .NE. int(0.5+(endTime-startTime)/deltaTClock) )
524 & THEN
525 WRITE(msgBuf,'(A)')
526 & 'S/R INI_PARMS: both endTime and nTimeSteps have been set'
527 CALL PRINT_ERROR( msgBuf , 1)
528 WRITE(msgBuf,'(A)')
529 & 'S/R INI_PARMS: but are inconsistent'
530 CALL PRINT_ERROR( msgBuf , 1)
531 STOP 'ABNORMAL END: S/R INI_PARMS'
532 ENDIF
533
534 C o If taveFreq is finite, then we must make sure the diagnostics
535 C code is being compiled
536 #ifndef INCLUDE_DIAGNOSTICS_INTERFACE_CODE
537 IF (taveFreq.NE.0.) THEN
538 WRITE(msgBuf,'(A)')
539 & 'S/R INI_PARMS: taveFreq <> 0 but you have'
540 CALL PRINT_ERROR( msgBuf , 1)
541 WRITE(msgBuf,'(A)')
542 & 'not compiled the model with the diagnostics routines.'
543 CALL PRINT_ERROR( msgBuf , 1)
544 WRITE(msgBuf,'(A,A)')
545 & 'Re-compile with: #define INCLUDE_DIAGNOSTICS_INTERFACE_CODE',
546 & ' or -DINCLUDE_DIAGNOSTICS_INTERFACE_CODE'
547 CALL PRINT_ERROR( msgBuf , 1)
548 STOP 'ABNORMAL END: S/R INI_PARMS'
549 ENDIF
550 #endif
551
552 C-- Grid parameters
553 C In cartesian coords distances are in metres
554 rkFac = UNSET_RS
555 DO K =1,Nr
556 delZ(K) = UNSET_RL
557 delP(K) = UNSET_RL
558 delR(K) = UNSET_RL
559 ENDDO
560 C In spherical polar distances are in degrees
561 recip_rSphere = 1.D0/rSphere
562 dxSpacing = UNSET_RL
563 dySpacing = UNSET_RL
564 delXfile = ' '
565 delYfile = ' '
566 READ(UNIT=iUnit,NML=PARM04) !,IOSTAT=errIO)
567 IF ( errIO .LT. 0 ) THEN
568 WRITE(msgBuf,'(A)')
569 & 'S/R INI_PARMS'
570 CALL PRINT_ERROR( msgBuf , 1)
571 WRITE(msgBuf,'(A)')
572 & 'Error reading numerical model '
573 CALL PRINT_ERROR( msgBuf , 1)
574 WRITE(msgBuf,'(A)')
575 & 'parameter file "data"'
576 CALL PRINT_ERROR( msgBuf , 1)
577 WRITE(msgBuf,'(A)')
578 & 'Problem in namelist PARM04'
579 CALL PRINT_ERROR( msgBuf , 1)
580 CALL MODELDATA_EXAMPLE( myThid )
581 STOP 'ABNORMAL END: S/R INI_PARMS'
582 ENDIF
583
584 C X coordinate
585 IF ( delXfile .NE. ' ' ) THEN
586 IF ( delX(1) .NE. UNSET_RL .OR. dxSpacing .NE. UNSET_RL ) THEN
587 WRITE(msgBuf,'(A,A)') 'Too many specifications for delX:',
588 & 'Specify only one of delX, dxSpacing or delXfile'
589 CALL PRINT_ERROR( msgBuf , 1)
590 STOP 'ABNORMAL END: S/R INI_PARMS'
591 ELSE
592 _BEGIN_MASTER( myThid )
593 IF (readBinaryPrec.EQ.precFloat32) THEN
594 OPEN(37,FILE=delXfile,STATUS='OLD',FORM='UNFORMATTED',
595 & ACCESS='DIRECT',RECL=WORDLENGTH*Nx)
596 READ(37,rec=1) delX
597 #ifdef _BYTESWAPIO
598 call MDS_BYTESWAPR4( Nx, delX )
599 #endif
600 CLOSE(37)
601 ELSEIF (readBinaryPrec.EQ.precFloat64) THEN
602 OPEN(37,FILE=delXfile,STATUS='OLD',FORM='UNFORMATTED',
603 & ACCESS='DIRECT',RECL=WORDLENGTH*2*Nx)
604 READ(37,rec=1) delX
605 #ifdef _BYTESWAPIO
606 call MDS_BYTESWAPR8( Nx, delX )
607 #endif
608 CLOSE(37)
609 ENDIF
610 _END_MASTER(myThid)
611 ENDIF
612 ENDIF
613 IF ( dxSpacing .NE. UNSET_RL ) THEN
614 DO i=1,Nx
615 delX(i) = dxSpacing
616 ENDDO
617 ENDIF
618 C Y coordinate
619 IF ( delYfile .NE. ' ' ) THEN
620 IF ( delY(1) .NE. UNSET_RL .OR. dySpacing .NE. UNSET_RL ) THEN
621 WRITE(msgBuf,'(A,A)') 'Too many specifications for delY:',
622 & 'Specify only one of delY, dySpacing or delYfile'
623 CALL PRINT_ERROR( msgBuf , 1)
624 STOP 'ABNORMAL END: S/R INI_PARMS'
625 ELSE
626 _BEGIN_MASTER( myThid )
627 IF (readBinaryPrec.EQ.precFloat32) THEN
628 OPEN(37,FILE=delYfile,STATUS='OLD',FORM='UNFORMATTED',
629 & ACCESS='DIRECT',RECL=WORDLENGTH*Ny)
630 READ(37,rec=1) delY
631 #ifdef _BYTESWAPIO
632 call MDS_BYTESWAPR4( Ny, delY )
633 #endif
634 CLOSE(37)
635 ELSEIF (readBinaryPrec.EQ.precFloat64) THEN
636 OPEN(37,FILE=delYfile,STATUS='OLD',FORM='UNFORMATTED',
637 & ACCESS='DIRECT',RECL=WORDLENGTH*2*Ny)
638 READ(37,rec=1) delY
639 #ifdef _BYTESWAPIO
640 call MDS_BYTESWAPR8( Ny, delY )
641 #endif
642 CLOSE(37)
643 ENDIF
644 _END_MASTER(myThid)
645 ENDIF
646 ENDIF
647 IF ( dySpacing .NE. UNSET_RL ) THEN
648 DO i=1,Ny
649 delY(i) = dySpacing
650 ENDDO
651 ENDIF
652 C
653 IF ( rSphere .NE. 0 ) THEN
654 recip_rSphere = 1.D0/rSphere
655 ELSE
656 recip_rSphere = 0.
657 ENDIF
658 C-- Initialize EOS coefficients (3rd order polynomial)
659 IF (eostype.eq.'POLY3') THEN
660 OPEN(37,FILE='POLY3.COEFFS',STATUS='OLD',FORM='FORMATTED')
661 READ(37,*) I
662 IF (I.NE.Nr) THEN
663 WRITE(msgBuf,'(A)')
664 & 'ini_parms: attempt to read POLY3.COEFFS failed'
665 CALL PRINT_ERROR( msgBuf , 1)
666 WRITE(msgBuf,'(A)')
667 & ' because bad # of levels in data'
668 CALL PRINT_ERROR( msgBuf , 1)
669 STOP 'Bad data in POLY3.COEFFS'
670 ENDIF
671 READ(37,*) (eosRefT(K),eosRefS(K),eosSig0(K),K=1,Nr)
672 DO K=1,Nr
673 READ(37,*) (eosC(I,K),I=1,9)
674 ENDDO
675 CLOSE(37)
676 ENDIF
677 C-- Check for conflicting grid definitions.
678 goptCount = 0
679 IF ( usingCartesianGrid ) goptCount = goptCount+1
680 IF ( usingSphericalPolarGrid ) goptCount = goptCount+1
681 IF ( goptCount .NE. 1 ) THEN
682 WRITE(msgBuf,'(A)')
683 & 'S/R INI_PARMS: More than one coordinate system requested'
684 CALL PRINT_ERROR( msgBuf , myThid)
685 STOP 'ABNORMAL END: S/R INI_PARMS'
686 ENDIF
687 C-- Make metric term settings consistent with underlying grid.
688 IF ( usingCartesianGrid ) THEN
689 usingSphericalPolarMterms = .FALSE.
690 metricTerms = .FALSE.
691 mTFacMom = 0
692 useBetaPlaneF = .TRUE.
693 ENDIF
694 IF ( usingSphericalPolarGrid ) THEN
695 useConstantF = .FALSE.
696 useBetaPlaneF = .FALSE.
697 useSphereF = .TRUE.
698 omega = 2.D0 * PI / ( 3600.D0 * 24.D0 )
699 usingSphericalPolarMterms = .TRUE.
700 metricTerms = .TRUE.
701 mTFacMom = 1
702 ENDIF
703 C-- p, z, r coord parameters
704 DO K = 1, Nr
705 IF ( delZ(K) .NE. UNSET_RL ) zCoordInputData = .TRUE.
706 IF ( delP(K) .NE. UNSET_RL ) pCoordInputData = .TRUE.
707 IF ( delR(K) .NE. UNSET_RL ) rCoordInputData = .TRUE.
708 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delZ(K)
709 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delP(K)
710 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delRDefault(K)
711 IF ( delR(K) .EQ. 0. ) THEN
712 WRITE(msgBuf,'(A,I4)')
713 & 'S/R INI_PARMS: No value for delZ/delP/delR at K = ',K
714 CALL PRINT_ERROR( msgBuf , 1)
715 STOP 'ABNORMAL END: S/R INI_PARMS'
716 ENDIF
717 ENDDO
718 C Check for multiple coordinate systems
719 CoordsSet = 0
720 IF ( zCoordInputData ) coordsSet = coordsSet + 1
721 IF ( pCoordInputData ) coordsSet = coordsSet + 1
722 IF ( rCoordInputData ) coordsSet = coordsSet + 1
723 IF ( coordsSet .GT. 1 ) THEN
724 WRITE(msgBuf,'(A)')
725 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
726 CALL PRINT_ERROR( msgBuf , myThid)
727 STOP 'ABNORMAL END: S/R INI_PARMS'
728 ENDIF
729
730 C-- Input files
731 READ(UNIT=iUnit,NML=PARM05) !,IOSTAT=errIO)
732 IF ( errIO .LT. 0 ) THEN
733 WRITE(msgBuf,'(A)')
734 & 'Error reading numerical model '
735 CALL PRINT_ERROR( msgBuf , 1)
736 WRITE(msgBuf,'(A)')
737 & 'parameter file "data"'
738 CALL PRINT_ERROR( msgBuf , 1)
739 WRITE(msgBuf,'(A)')
740 & 'Problem in namelist PARM05'
741 CALL PRINT_ERROR( msgBuf , 1)
742 CALL MODELDATA_EXAMPLE( myThid )
743 STOP 'ABNORMAL END: S/R INI_PARMS'
744 ENDIF
745
746 C
747 C-- Set factors required for mixing pressure and meters as vertical coordinate.
748 C rkFac is a "sign" parameter which is used where the orientation of the vertical
749 C coordinate (pressure or meters) relative to the vertical index (K) is important.
750 C rkFac = 1 applies when K and the coordinate are in the opposite sense.
751 C rkFac = -1 applies when K and the coordinate are in the same sense.
752 C horiVertRatio is a parameter that maps horizontal units to vertical units.
753 C It is used in certain special cases where lateral and vertical terms are
754 C being combined and a single frame of reference is needed.
755 IF ( zCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
756 rkFac = 1.D0
757 horiVertRatio = 1.D0
758 ENDIF
759 IF ( pCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
760 rkFac = -1.D0
761 horiVertRatio = Gravity * rhoConst
762 ENDIF
763 IF ( rCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
764 rkFac = 1.D0
765 horiVertRatio = 1.D0
766 ENDIF
767 IF (buoyancyRelation.EQ.'ATMOSPHERIC')
768 & horiVertRatio = Gravity * rhoConst
769 IF ( rkFac .EQ. UNSET_RS ) rkFac=rkFacDefault
770 recip_rkFac = 1.D0 / rkFac
771 recip_horiVertRatio = 1./horiVertRatio
772 IF ( zCoordInputData ) usingZCoords = .TRUE.
773 IF ( pCoordInputData ) usingPCoords = .TRUE.
774
775 C
776 CLOSE(iUnit)
777
778 _END_MASTER(myThid)
779
780 C-- Everyone else must wait for the parameters to be loaded
781 _BARRIER
782 C
783
784 RETURN
785 END
786

  ViewVC Help
Powered by ViewVC 1.1.22