/[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.88 - (show annotations) (download)
Thu Nov 7 21:51:15 2002 UTC (21 years, 6 months ago) by adcroft
Branch: MAIN
CVS Tags: checkpoint46n_post, checkpoint47a_post, checkpoint47b_post, checkpoint46m_post, checkpoint47
Changes since 1.87: +9 -2 lines
Added new routine quasihydrostaticterms() and flag "quasihydrostatic"
which is false by default and enables QH mode. Exlcusive with nonhydrostatic
flag.

1 C $Header: /u/gcmpack/MITgcm/model/src/ini_parms.F,v 1.87 2002/11/05 18:48:56 adcroft Exp $
2 C $Name: $
3
4 #include "CPP_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: INI_PARMS
8 C !INTERFACE:
9 SUBROUTINE INI_PARMS( myThid )
10 C !DESCRIPTION: \bv
11 C *==========================================================*
12 C | SUBROUTINE INI_PARMS
13 C | o Routine to set model "parameters"
14 C *==========================================================*
15 C | Notes:
16 C | ======
17 C | The present version of this routine is a place-holder.
18 C | A production version needs to handle parameters from an
19 C | external file and possibly reading in some initial field
20 C | values.
21 C *==========================================================*
22 C \ev
23
24 C !USES:
25 IMPLICIT NONE
26 C === Global variables ===
27 #include "SIZE.h"
28 #include "EEPARAMS.h"
29 #include "PARAMS.h"
30 #include "GRID.h"
31 #include "EOS.h"
32
33 C !INPUT/OUTPUT PARAMETERS:
34 C === Routine arguments ===
35 C myThid - Number of this instance of INI_PARMS
36 INTEGER myThid
37
38 C !LOCAL VARIABLES:
39 C === Local variables ===
40 C dxSpacing, dySpacing - Default spacing in X and Y.
41 C Units are that of coordinate system
42 C i.e. cartesian => metres
43 C s. polar => degrees
44 C goptCount - Used to count the nuber of grid options
45 C (only one is allowed! )
46 C msgBuf - Informational/error meesage buffer
47 C errIO - IO error flag
48 C iUnit - Work variable for IO unit number
49 C record - Work variable for IO buffer
50 C K, I, J - Loop counters
51 C xxxDefault - Default value for variable xxx
52 _RL dxSpacing
53 _RL dySpacing
54 CHARACTER*(MAX_LEN_FNAM) delXfile
55 CHARACTER*(MAX_LEN_FNAM) delYfile
56 CHARACTER*(MAX_LEN_MBUF) msgBuf
57 CHARACTER*(MAX_LEN_PREC) record
58 INTEGER goptCount
59 INTEGER K, I, IL, iUnit
60 INTEGER errIO
61 INTEGER IFNBLNK
62 EXTERNAL IFNBLNK
63 INTEGER ILNBLNK
64 EXTERNAL ILNBLNK
65 C Default values for variables which have vertical coordinate system
66 C dependency.
67 _RL viscArDefault
68 _RL diffKrTDefault
69 _RL diffKrSDefault
70 _RL hFacMinDrDefault
71 _RL delRDefault(Nr)
72 _RS rkFacDefault
73 C zCoordInputData :: Variables used to select between different coordinate systems.
74 C pCoordInputData :: The vertical coordinate system in the rest of the model is
75 C rCoordInputData :: written in terms of r. In the model "data" file input data can
76 C coordsSet :: be interms of z, p or r.
77 C :: e.g. delZ or delP or delR for the vertical grid spacing.
78 C :: The following rules apply:
79 C :: All parameters must use the same vertical coordinate system.
80 C :: e.g. delZ and viscAz is legal but
81 C :: delZ and viscAr is an error.
82 C :: Similarly specifyinh delZ and delP is an error.
83 C :: zCoord..., pCoord..., rCoord... are used to flag when z, p or r are
84 C :: used. coordsSet counts how many vertical coordinate systems have been
85 C used to specify variables. coordsSet > 1 is an error.
86 C
87 LOGICAL zCoordInputData
88 LOGICAL pCoordInputData
89 LOGICAL rCoordInputData
90 INTEGER coordsSet
91
92 C Retired main data file parameters. Kept here to trap use of old data files.
93 C zonal_filt_lat - Moved to package "zonal_filt"
94 C nRetired :: Counter used to trap gracefully namelists containing "retired"
95 C :: parameters. These are parameters that are either no-longer used
96 C or that have moved to a different input file and/or namelist.
97 _RL zonal_filt_lat
98 INTEGER nRetired
99 CEOP
100
101 C-- Continuous equation parameters
102 NAMELIST /PARM01/
103 & gravitySign,
104 & gravity, gBaro, rhonil, tAlpha, sBeta, f0, beta, omega,
105 & viscAh, viscAz, viscA4, cosPower, viscAstrain, viscAtension,
106 & diffKhT, diffKzT, diffK4T,
107 & diffKhS, diffKzS, diffK4S,
108 & tRef, sRef, eosType, Integr_GeoPot,
109 & no_slip_sides,no_slip_bottom,
110 & momViscosity, momAdvection, momForcing, useCoriolis,
111 & momPressureForcing, metricTerms, vectorInvariantMomentum,
112 & tempDiffusion, tempAdvection, tempForcing,
113 & saltDiffusion, saltAdvection, saltForcing,
114 & implicSurfPress, implicDiv2DFlow,
115 & implicitFreeSurface, rigidLid, freeSurfFac, hFacMin, hFacMinDz,
116 & exactConserv,uniformLin_PhiSurf,nonlinFreeSurf,hFacInf,hFacSup,
117 & staggerTimeStep,
118 & tempStepping, saltStepping, momStepping, tr1Stepping,
119 & implicitDiffusion, implicitViscosity,
120 & viscAr, diffKrT, diffKrS, hFacMinDr,
121 & viscAp, diffKpT, diffKpS, hFacMinDp,
122 & rhoConst, rhoConstFresh, buoyancyRelation, HeatCapacity_Cp,
123 & writeBinaryPrec, readBinaryPrec, writeStatePrec,
124 & nonHydrostatic, quasiHydrostatic, globalFiles,
125 & allowFreezing, ivdc_kappa,
126 & bottomDragLinear,bottomDragQuadratic,
127 & usePickupBeforeC35, debugMode,
128 & readPickupWithTracer, writePickupWithTracer,
129 & tempAdvScheme, saltAdvScheme, tracerAdvScheme,
130 & multiDimAdvection, useEnergyConservingCoriolis,
131 & useJamartWetPoints, useNHMTerms,
132 & useRealFreshWaterFlux, convertFW2Salt,
133 & temp_EvPrRn, salt_EvPrRn, trac_EvPrRn,
134 & zonal_filt_lat
135
136 C-- Elliptic solver parameters
137 NAMELIST /PARM02/
138 & cg2dMaxIters, cg2dChkResFreq, cg2dTargetResidual,
139 & cg2dTargetResWunit, cg2dpcOffDFac,
140 & cg3dMaxIters, cg3dChkResFreq, cg3dTargetResidual
141
142 C-- Time stepping parammeters
143 NAMELIST /PARM03/
144 & nIter0, nTimeSteps, nEndIter,
145 & deltaT, deltaTmom, deltaTtracer, deltaTfreesurf,
146 & forcing_In_AB, abEps, tauCD, rCD,
147 & startTime, endTime, chkPtFreq,
148 & dumpFreq, taveFreq, tave_lastIter, deltaTClock, diagFreq,
149 & monitorFreq, pChkPtFreq, cAdjFreq,
150 & tauThetaClimRelax, tauSaltClimRelax, tauTr1ClimRelax,
151 & periodicExternalForcing, externForcingPeriod, externForcingCycle
152
153 C-- Gridding parameters
154 NAMELIST /PARM04/
155 & usingCartesianGrid, dxSpacing, dySpacing, delX, delY, delZ,
156 & usingSphericalPolarGrid, phiMin, thetaMin, rSphere,
157 & usingCurvilinearGrid,
158 & delP, delR, rkFac, Ro_SeaLevel, groundAtK1, delRc,
159 & delXfile, delYfile
160
161 C-- Input files
162 NAMELIST /PARM05/
163 & bathyFile, topoFile, hydrogThetaFile, hydrogSaltFile,
164 & zonalWindFile, meridWindFile,
165 & thetaClimFile, saltClimFile,
166 & surfQfile, EmPmRfile, surfQswfile,
167 & uVelInitFile, vVelInitFile, pSurfInitFile,
168 & dQdTFile, ploadFile
169
170 C
171 _BEGIN_MASTER(myThid)
172
173 C Defaults values for input parameters
174 CALL SET_DEFAULTS(
175 O viscArDefault, diffKrTDefault, diffKrSDefault,
176 O hFacMinDrDefault, delRdefault, rkFacDefault,
177 I myThid )
178
179 C-- Initialise "which vertical coordinate system used" flags.
180 zCoordInputData = .FALSE.
181 pCoordInputData = .FALSE.
182 rCoordInputData = .FALSE.
183 usingPCoords = .FALSE.
184 usingZCoords = .FALSE.
185 coordsSet = 0
186
187 C-- Iniialise retired parameters to unlikely value
188 nRetired = 0
189 zonal_filt_lat = UNSET_RL
190
191 C-- Open the parameter file
192 OPEN(UNIT=scrUnit1,STATUS='SCRATCH')
193 OPEN(UNIT=scrUnit2,STATUS='SCRATCH')
194 OPEN(UNIT=modelDataUnit,FILE='data',STATUS='OLD',
195 & IOSTAT=errIO)
196 IF ( errIO .LT. 0 ) THEN
197 WRITE(msgBuf,'(A)')
198 & 'S/R INI_PARMS'
199 CALL PRINT_ERROR( msgBuf , 1)
200 WRITE(msgBuf,'(A)')
201 & 'Unable to open model parameter'
202 CALL PRINT_ERROR( msgBuf , 1)
203 WRITE(msgBuf,'(A)')
204 & 'file "data"'
205 CALL PRINT_ERROR( msgBuf , 1)
206 CALL MODELDATA_EXAMPLE( myThid )
207 STOP 'ABNORMAL END: S/R INI_PARMS'
208 ENDIF
209
210 DO WHILE ( .TRUE. )
211 READ(modelDataUnit,FMT='(A)',END=1001) RECORD
212 IL = MAX(ILNBLNK(RECORD),1)
213 IF ( RECORD(1:1) .NE. commentCharacter )
214 & WRITE(UNIT=scrUnit1,FMT='(A)') RECORD(:IL)
215 WRITE(UNIT=scrUnit2,FMT='(A)') RECORD(:IL)
216 ENDDO
217 1001 CONTINUE
218 CLOSE(modelDataUnit)
219
220 C-- Report contents of model parameter file
221 WRITE(msgBuf,'(A)')
222 &'// ======================================================='
223 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
224 & SQUEEZE_RIGHT , 1)
225 WRITE(msgBuf,'(A)') '// Model parameter file "data"'
226 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
227 & SQUEEZE_RIGHT , 1)
228 WRITE(msgBuf,'(A)')
229 &'// ======================================================='
230 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
231 & SQUEEZE_RIGHT , 1)
232 iUnit = scrUnit2
233 REWIND(iUnit)
234 DO WHILE ( .TRUE. )
235 READ(UNIT=iUnit,FMT='(A)',END=2001) RECORD
236 IL = MAX(ILNBLNK(RECORD),1)
237 WRITE(msgBuf,'(A,A)') '>',RECORD(:IL)
238 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
239 & SQUEEZE_RIGHT , 1)
240 ENDDO
241 2001 CONTINUE
242 CLOSE(iUnit)
243 WRITE(msgBuf,'(A)') ' '
244 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
245 & SQUEEZE_RIGHT , 1)
246
247
248 C-- Read settings from model parameter file "data".
249 iUnit = scrUnit1
250 REWIND(iUnit)
251
252 C-- Set default "physical" parameters
253 viscAz = UNSET_RL
254 viscAr = UNSET_RL
255 viscAp = UNSET_RL
256 diffKzT = UNSET_RL
257 diffKpT = UNSET_RL
258 diffKrT = UNSET_RL
259 diffKzS = UNSET_RL
260 diffKpS = UNSET_RL
261 diffKrS = UNSET_RL
262 gBaro = UNSET_RL
263 rhoConst = UNSET_RL
264 hFacMinDr = UNSET_RL
265 hFacMinDz = UNSET_RL
266 hFacMinDp = UNSET_RL
267 convertFW2Salt = UNSET_RL
268 tAlpha = UNSET_RL
269 sBeta = UNSET_RL
270 READ(UNIT=iUnit,NML=PARM01) !,IOSTAT=errIO)
271 IF ( errIO .LT. 0 ) THEN
272 WRITE(msgBuf,'(A)')
273 & 'S/R INI_PARMS'
274 CALL PRINT_ERROR( msgBuf , 1)
275 WRITE(msgBuf,'(A)')
276 & 'Error reading numerical model '
277 CALL PRINT_ERROR( msgBuf , 1)
278 WRITE(msgBuf,'(A)')
279 & 'parameter file "data"'
280 CALL PRINT_ERROR( msgBuf , 1)
281 WRITE(msgBuf,'(A)')
282 & 'Problem in namelist PARM01'
283 CALL PRINT_ERROR( msgBuf , 1)
284 CALL MODELDATA_EXAMPLE( myThid )
285 STOP 'ABNORMAL END: S/R INI_PARMS'
286 ELSE
287 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM01 : OK'
288 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
289 & SQUEEZE_RIGHT , 1)
290 ENDIF
291 IF ( implicitFreeSurface ) freeSurfFac = 1.D0
292 IF ( rigidLid ) freeSurfFac = 0.D0
293 IF ( gBaro .EQ. UNSET_RL ) gBaro=gravity
294 IF ( rhoConst .EQ. UNSET_RL ) rhoConst=rhoNil
295 C-- Momentum viscosity on/off flag.
296 IF ( momViscosity ) THEN
297 vfFacMom = 1.D0
298 ELSE
299 vfFacMom = 0.D0
300 ENDIF
301 C-- Momentum advection on/off flag.
302 IF ( momAdvection ) THEN
303 afFacMom = 1.D0
304 ELSE
305 afFacMom = 0.D0
306 ENDIF
307 C-- Momentum forcing on/off flag.
308 IF ( momForcing ) THEN
309 foFacMom = 1.D0
310 ELSE
311 foFacMom = 0.D0
312 ENDIF
313 C-- Coriolis term on/off flag.
314 IF ( useCoriolis ) THEN
315 cfFacMom = 1.D0
316 ELSE
317 cfFacMom = 0.D0
318 ENDIF
319 C-- Pressure term on/off flag.
320 IF ( momPressureForcing ) THEN
321 pfFacMom = 1.D0
322 ELSE
323 pfFacMom = 0.D0
324 ENDIF
325 C-- Metric terms on/off flag.
326 IF ( metricTerms ) THEN
327 mTFacMom = 1.D0
328 ELSE
329 mTFacMom = 0.D0
330 ENDIF
331 C-- Non-hydrostatic/quasi-hydrostatic
332 IF (nonHydrostatic.AND.quasiHydrostatic) THEN
333 WRITE(msgBuf,'(A)')
334 & 'Illegal: both nonHydrostatic = quasiHydrostatic = TRUE'
335 CALL PRINT_ERROR( msgBuf , myThid)
336 STOP 'ABNORMAL END: S/R INI_PARMS'
337 ENDIF
338 C-- Advection and Forcing for Temp and salt on/off flags
339 tempAdvection = tempStepping .AND. tempAdvection
340 tempForcing = tempStepping .AND. tempForcing
341 saltAdvection = saltStepping .AND. saltAdvection
342 saltForcing = saltStepping .AND. saltForcing
343 C-- z,p,r coord input switching.
344 IF ( viscAz .NE. UNSET_RL ) zCoordInputData = .TRUE.
345 IF ( viscAp .NE. UNSET_RL ) pCoordInputData = .TRUE.
346 IF ( viscAr .NE. UNSET_RL ) rCoordInputData = .TRUE.
347 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAz
348 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAp
349 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscArDefault
350
351 IF ( diffKzT .NE. UNSET_RL ) zCoordInputData = .TRUE.
352 IF ( diffKpT .NE. UNSET_RL ) pCoordInputData = .TRUE.
353 IF ( diffKrT .NE. UNSET_RL ) rCoordInputData = .TRUE.
354 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKzT
355 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKpT
356 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKrTDefault
357
358 IF ( diffKzS .NE. UNSET_RL ) zCoordInputData = .TRUE.
359 IF ( diffKpS .NE. UNSET_RL ) pCoordInputData = .TRUE.
360 IF ( diffKrS .NE. UNSET_RL ) rCoordInputData = .TRUE.
361 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKzS
362 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKpS
363 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKrSDefault
364
365 IF ( hFacMinDz .NE. UNSET_RL ) zCoordInputData = .TRUE.
366 IF ( hFacMinDp .NE. UNSET_RL ) pCoordInputData = .TRUE.
367 IF ( hFacMinDr .NE. UNSET_RL ) rCoordInputData = .TRUE.
368 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDz
369 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDp
370 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDrDefault
371
372 IF (convertFW2Salt.EQ.UNSET_RL) THEN
373 convertFW2Salt = 35.
374 IF (useRealFreshWaterFlux) convertFW2Salt=-1
375 ENDIF
376
377 IF ( ivdc_kappa .NE. 0. .AND. .NOT. implicitDiffusion ) THEN
378 WRITE(msgBuf,'(A,A)')
379 & 'S/R INI_PARMS: To use ivdc_kappa you must enable implicit',
380 & ' vertical diffusion.'
381 CALL PRINT_ERROR( msgBuf , myThid)
382 STOP 'ABNORMAL END: S/R INI_PARMS'
383 ENDIF
384
385 coordsSet = 0
386 IF ( zCoordInputData ) coordsSet = coordsSet + 1
387 IF ( pCoordInputData ) coordsSet = coordsSet + 1
388 IF ( rCoordInputData ) coordsSet = coordsSet + 1
389 IF ( coordsSet .GT. 1 ) THEN
390 WRITE(msgBuf,'(A)')
391 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
392 CALL PRINT_ERROR( msgBuf , myThid)
393 STOP 'ABNORMAL END: S/R INI_PARMS'
394 ENDIF
395 IF ( rhoConst .LE. 0. ) THEN
396 WRITE(msgBuf,'(A)')
397 & 'S/R INI_PARMS: rhoConst must be greater than 0.'
398 CALL PRINT_ERROR( msgBuf , myThid)
399 STOP 'ABNORMAL END: S/R INI_PARMS'
400 ELSE
401 recip_rhoConst = 1.D0 / rhoConst
402 ENDIF
403 IF ( rhoNil .LE. 0. ) THEN
404 WRITE(msgBuf,'(A)')
405 & 'S/R INI_PARMS: rhoNil must be greater than 0.'
406 CALL PRINT_ERROR( msgBuf , myThid)
407 STOP 'ABNORMAL END: S/R INI_PARMS'
408 ELSE
409 recip_rhoNil = 1.D0 / rhoNil
410 ENDIF
411 IF ( HeatCapacity_Cp .LE. 0. ) THEN
412 WRITE(msgBuf,'(A)')
413 & 'S/R INI_PARMS: HeatCapacity_Cp must be greater than 0.'
414 CALL PRINT_ERROR( msgBuf , myThid)
415 STOP 'ABNORMAL END: S/R INI_PARMS'
416 ELSE
417 recip_Cp = 1.D0 / HeatCapacity_Cp
418 ENDIF
419 IF ( gravity .LE. 0. ) THEN
420 WRITE(msgBuf,'(A)')
421 & 'S/R INI_PARMS: gravity must be greater than 0.'
422 CALL PRINT_ERROR( msgBuf , myThid)
423 STOP 'ABNORMAL END: S/R INI_PARMS'
424 ELSE
425 recip_gravity = 1.D0 / gravity
426 ENDIF
427 C Set globalFiles flag for READ_WRITE_FLD package
428 CALL SET_WRITE_GLOBAL_FLD( globalFiles )
429 C Set globalFiles flag for READ_WRITE_REC package
430 CALL SET_WRITE_GLOBAL_REC( globalFiles )
431 C Set globalFiles flag for READ_WRITE_REC package
432 CALL SET_WRITE_GLOBAL_PICKUP( globalFiles )
433
434 C Check for retired parameters still being used
435 nRetired = 0
436 IF ( zonal_filt_lat .NE. UNSET_RL ) THEN
437 nRetired = nRetired+1
438 WRITE(msgBuf,'(A,A)')
439 & 'S/R INI_PARMS: Paramater "zonal_filt_lat" is',
440 & ' no longer allowed in file "data".'
441 CALL PRINT_ERROR( msgBuf , myThid)
442 WRITE(msgBuf,'(A,A)')
443 & 'S/R INI_PARMS: Paramater "zonal_filt_lat" is',
444 & ' now read from file "data.zonfilt".'
445 CALL PRINT_ERROR( msgBuf , myThid)
446 ENDIF
447
448 C-- Elliptic solver parameters
449 READ(UNIT=iUnit,NML=PARM02) !,IOSTAT=errIO)
450 IF ( errIO .LT. 0 ) THEN
451 WRITE(msgBuf,'(A)')
452 & 'S/R INI_PARMS'
453 CALL PRINT_ERROR( msgBuf , 1)
454 WRITE(msgBuf,'(A)')
455 & 'Error reading numerical model '
456 CALL PRINT_ERROR( msgBuf , 1)
457 WRITE(msgBuf,'(A)')
458 & 'parameter file "data".'
459 CALL PRINT_ERROR( msgBuf , 1)
460 WRITE(msgBuf,'(A)')
461 & 'Problem in namelist PARM02'
462 CALL PRINT_ERROR( msgBuf , 1)
463 CALL MODELDATA_EXAMPLE( myThid )
464 STOP 'ABNORMAL END: S/R INI_PARMS'
465 ELSE
466 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM02 : OK'
467 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
468 & SQUEEZE_RIGHT , 1)
469 ENDIF
470
471 C-- Time stepping parameters
472 rCD = -1.D0
473 READ(UNIT=iUnit,NML=PARM03) !,IOSTAT=errIO)
474 IF ( errIO .LT. 0 ) THEN
475 WRITE(msgBuf,'(A)')
476 & 'S/R INI_PARMS'
477 CALL PRINT_ERROR( msgBuf , 1)
478 WRITE(msgBuf,'(A)')
479 & 'Error reading numerical model '
480 CALL PRINT_ERROR( msgBuf , 1)
481 WRITE(msgBuf,'(A)')
482 & 'parameter file "data"'
483 CALL PRINT_ERROR( msgBuf , 1)
484 WRITE(msgBuf,'(A)')
485 & 'Problem in namelist PARM03'
486 CALL PRINT_ERROR( msgBuf , 1)
487 CALL MODELDATA_EXAMPLE( myThid )
488 STOP 'ABNORMAL END: S/R INI_PARMS'
489 ELSE
490 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM03 : OK'
491 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
492 & SQUEEZE_RIGHT , 1)
493 ENDIF
494 C Process "timestepping" params
495 C o Time step size
496 IF ( deltaT .EQ. 0. ) deltaT = deltaTmom
497 IF ( deltaT .EQ. 0. ) deltaT = deltaTtracer
498 IF ( deltaTmom .EQ. 0. ) deltaTmom = deltaT
499 IF ( deltaTtracer .EQ. 0. ) deltaTtracer = deltaT
500 IF ( deltaTClock .EQ. 0. ) deltaTClock = deltaT
501 C Note that this line should set deltaTreesurf=deltaTmom
502 C but this would change a lot of existing set-ups so we are
503 C obliged to set the default inappropriately.
504 C Be advised that when using asynchronous time stepping
505 C it is better to set deltaTreesurf=deltaTtracer
506 IF ( deltaTfreesurf .EQ. 0. ) deltaTfreesurf = deltaTmom
507 IF ( periodicExternalForcing ) THEN
508 IF ( externForcingCycle*externForcingPeriod .EQ. 0. ) THEN
509 WRITE(msgBuf,'(A)')
510 & 'S/R INI_PARMS: externForcingCycle,externForcingPeriod =0'
511 CALL PRINT_ERROR( msgBuf , 1)
512 STOP 'ABNORMAL END: S/R INI_PARMS'
513 ENDIF
514 IF ( INT(externForcingCycle/externForcingPeriod) .NE.
515 & externForcingCycle/externForcingPeriod ) THEN
516 WRITE(msgBuf,'(A)')
517 & 'S/R INI_PARMS: externForcingCycle <> N*externForcingPeriod'
518 CALL PRINT_ERROR( msgBuf , 1)
519 STOP 'ABNORMAL END: S/R INI_PARMS'
520 ENDIF
521 IF ( externForcingCycle.le.externForcingPeriod ) THEN
522 WRITE(msgBuf,'(A)')
523 & 'S/R INI_PARMS: externForcingCycle < externForcingPeriod'
524 CALL PRINT_ERROR( msgBuf , 1)
525 STOP 'ABNORMAL END: S/R INI_PARMS'
526 ENDIF
527 IF ( externForcingPeriod.lt.deltaTclock ) THEN
528 WRITE(msgBuf,'(A)')
529 & 'S/R INI_PARMS: externForcingPeriod < deltaTclock'
530 CALL PRINT_ERROR( msgBuf , 1)
531 STOP 'ABNORMAL END: S/R INI_PARMS'
532 ENDIF
533 ENDIF
534 C o Convection frequency
535 IF ( cAdjFreq .LT. 0. ) THEN
536 cAdjFreq = deltaTClock
537 ENDIF
538 IF ( ivdc_kappa .NE. 0. .AND. cAdjFreq .NE. 0. ) THEN
539 WRITE(msgBuf,'(A,A)')
540 & 'S/R INI_PARMS: You have enabled both ivdc_kappa and',
541 & ' convective adjustment.'
542 CALL PRINT_ERROR( msgBuf , myThid)
543 STOP 'ABNORMAL END: S/R INI_PARMS'
544 ENDIF
545 C o CD coupling
546 IF ( tauCD .EQ. 0.D0 ) THEN
547 tauCD = deltaTmom
548 ENDIF
549 IF ( rCD .LT. 0. ) THEN
550 rCD = 1. - deltaTMom/tauCD
551 ENDIF
552 C o Temperature climatology relaxation time scale
553 IF ( tauThetaClimRelax .EQ. 0.D0 ) THEN
554 doThetaClimRelax = .FALSE.
555 lambdaThetaClimRelax = 0.D0
556 ELSE
557 doThetaClimRelax = .TRUE.
558 lambdaThetaClimRelax = 1./tauThetaClimRelax
559 ENDIF
560 C o Salinity climatology relaxation time scale
561 IF ( tauSaltClimRelax .EQ. 0.D0 ) THEN
562 doSaltClimRelax = .FALSE.
563 lambdaSaltClimRelax = 0.D0
564 ELSE
565 doSaltClimRelax = .TRUE.
566 lambdaSaltClimRelax = 1./tauSaltClimRelax
567 ENDIF
568 C o Tracer 1 climatology relaxation time scale
569 IF ( tauTr1ClimRelax .EQ. 0.D0 ) THEN
570 doTr1ClimRelax = .FALSE.
571 lambdaTr1ClimRelax = 0.D0
572 ELSE
573 doTr1ClimRelax = .TRUE.
574 lambdaTr1ClimRelax = 1./tauTr1ClimRelax
575 ENDIF
576
577 C o Start time
578 IF ( nIter0 .NE. 0 .AND. startTime .EQ. 0. )
579 & startTime = deltaTClock*float(nIter0)
580 C o nIter0
581 IF ( nIter0 .EQ. 0 .AND. startTime .NE. 0. )
582 & nIter0 = INT( startTime/deltaTClock )
583
584 C o nTimeSteps 1
585 IF ( nTimeSteps .EQ. 0 .AND. nEndIter .NE. 0 )
586 & nTimeSteps = nEndIter-nIter0
587 C o nTimeSteps 2
588 IF ( nTimeSteps .EQ. 0 .AND. endTime .NE. 0. )
589 & nTimeSteps = int(0.5+(endTime-startTime)/deltaTclock)
590 C o nEndIter 1
591 IF ( nEndIter .EQ. 0 .AND. nTimeSteps .NE. 0 )
592 & nEndIter = nIter0+nTimeSteps
593 C o nEndIter 2
594 IF ( nEndIter .EQ. 0 .AND. endTime .NE. 0. )
595 & nEndIter = int(0.5+endTime/deltaTclock)
596 C o End Time 1
597 IF ( endTime .EQ. 0. .AND. nTimeSteps .NE. 0 )
598 & endTime = startTime + deltaTClock*float(nTimeSteps)
599 C o End Time 2
600 IF ( endTime .EQ. 0. .AND. nEndIter .NE. 0 )
601 & endTime = deltaTClock*float(nEndIter)
602
603 C o Consistent?
604 IF ( nEndIter .NE. nIter0+nTimeSteps ) THEN
605 WRITE(msgBuf,'(A)')
606 & 'S/R INI_PARMS: nIter0, nTimeSteps and nEndIter are inconsistent'
607 CALL PRINT_ERROR( msgBuf , 1)
608 WRITE(msgBuf,'(A)')
609 & 'S/R INI_PARMS: Perhaps more than two were set at once'
610 CALL PRINT_ERROR( msgBuf , 1)
611 STOP 'ABNORMAL END: S/R INI_PARMS'
612 ENDIF
613 IF ( nTimeSteps .NE. int(0.5+(endTime-startTime)/deltaTClock) )
614 & THEN
615 WRITE(msgBuf,'(A)')
616 & 'S/R INI_PARMS: both endTime and nTimeSteps have been set'
617 CALL PRINT_ERROR( msgBuf , 1)
618 WRITE(msgBuf,'(A)')
619 & 'S/R INI_PARMS: but are inconsistent'
620 CALL PRINT_ERROR( msgBuf , 1)
621 STOP 'ABNORMAL END: S/R INI_PARMS'
622 ENDIF
623
624 C o Monitor (should also add CPP flag for monitor?)
625 IF (monitorFreq.LT.0.) THEN
626 monitorFreq=0.
627 IF (dumpFreq.NE.0.) monitorFreq=dumpFreq
628 IF (diagFreq.NE.0..AND.diagFreq.LT.monitorFreq)
629 & monitorFreq=diagFreq
630 IF (taveFreq.NE.0..AND.taveFreq.LT.monitorFreq)
631 & monitorFreq=taveFreq
632 IF (chkPtFreq.NE.0..AND.chkPtFreq.LT.monitorFreq)
633 & monitorFreq=chkPtFreq
634 IF (pChkPtFreq.NE.0..AND.pChkPtFreq.LT.monitorFreq)
635 & monitorFreq=pChkPtFreq
636 IF (monitorFreq.EQ.0.) monitorFreq=deltaTclock
637 ENDIF
638
639 C o If taveFreq is finite, then we must make sure the diagnostics
640 C code is being compiled
641 #ifndef ALLOW_TIMEAVE
642 IF (taveFreq.NE.0.) THEN
643 WRITE(msgBuf,'(A)')
644 & 'S/R INI_PARMS: taveFreq <> 0 but you have'
645 CALL PRINT_ERROR( msgBuf , 1)
646 WRITE(msgBuf,'(A)')
647 & 'not compiled the model with the diagnostics routines.'
648 CALL PRINT_ERROR( msgBuf , 1)
649 WRITE(msgBuf,'(A,A)')
650 & 'Re-compile with: #define ALLOW_TIMEAVE',
651 & ' or -DALLOW_TIMEAVE'
652 CALL PRINT_ERROR( msgBuf , 1)
653 STOP 'ABNORMAL END: S/R INI_PARMS'
654 ENDIF
655 #endif
656
657 C-- Grid parameters
658 C In cartesian coords distances are in metres
659 rkFac = UNSET_RS
660 DO K =1,Nr
661 delZ(K) = UNSET_RL
662 delP(K) = UNSET_RL
663 delR(K) = UNSET_RL
664 ENDDO
665 C In spherical polar distances are in degrees
666 recip_rSphere = 1.D0/rSphere
667 dxSpacing = UNSET_RL
668 dySpacing = UNSET_RL
669 delXfile = ' '
670 delYfile = ' '
671 READ(UNIT=iUnit,NML=PARM04,IOSTAT=errIO)
672 IF ( errIO .LT. 0 ) THEN
673 WRITE(msgBuf,'(A)')
674 & 'S/R INI_PARMS'
675 CALL PRINT_ERROR( msgBuf , 1)
676 WRITE(msgBuf,'(A)')
677 & 'Error reading numerical model '
678 CALL PRINT_ERROR( msgBuf , 1)
679 WRITE(msgBuf,'(A)')
680 & 'parameter file "data"'
681 CALL PRINT_ERROR( msgBuf , 1)
682 WRITE(msgBuf,'(A)')
683 & 'Problem in namelist PARM04'
684 CALL PRINT_ERROR( msgBuf , 1)
685 CALL MODELDATA_EXAMPLE( myThid )
686 STOP 'ABNORMAL END: S/R INI_PARMS'
687 ELSE
688 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM04 : OK'
689 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
690 & SQUEEZE_RIGHT , 1)
691 ENDIF
692
693 C X coordinate
694 IF ( delXfile .NE. ' ' ) THEN
695 IF ( delX(1) .NE. UNSET_RL .OR. dxSpacing .NE. UNSET_RL ) THEN
696 WRITE(msgBuf,'(A,A)') 'Too many specifications for delX:',
697 & 'Specify only one of delX, dxSpacing or delXfile'
698 CALL PRINT_ERROR( msgBuf , 1)
699 STOP 'ABNORMAL END: S/R INI_PARMS'
700 ELSE
701 _BEGIN_MASTER( myThid )
702 IF (readBinaryPrec.EQ.precFloat32) THEN
703 OPEN(37,FILE=delXfile,STATUS='OLD',FORM='UNFORMATTED',
704 & ACCESS='DIRECT',RECL=WORDLENGTH*Nx)
705 READ(37,rec=1) delX
706 #ifdef _BYTESWAPIO
707 call MDS_BYTESWAPR4( Nx, delX )
708 #endif
709 CLOSE(37)
710 ELSEIF (readBinaryPrec.EQ.precFloat64) THEN
711 OPEN(37,FILE=delXfile,STATUS='OLD',FORM='UNFORMATTED',
712 & ACCESS='DIRECT',RECL=WORDLENGTH*2*Nx)
713 READ(37,rec=1) delX
714 #ifdef _BYTESWAPIO
715 call MDS_BYTESWAPR8( Nx, delX )
716 #endif
717 CLOSE(37)
718 ENDIF
719 _END_MASTER(myThid)
720 ENDIF
721 ENDIF
722 IF ( dxSpacing .NE. UNSET_RL ) THEN
723 DO i=1,Nx
724 delX(i) = dxSpacing
725 ENDDO
726 ENDIF
727 C Y coordinate
728 IF ( delYfile .NE. ' ' ) THEN
729 IF ( delY(1) .NE. UNSET_RL .OR. dySpacing .NE. UNSET_RL ) THEN
730 WRITE(msgBuf,'(A,A)') 'Too many specifications for delY:',
731 & 'Specify only one of delY, dySpacing or delYfile'
732 CALL PRINT_ERROR( msgBuf , 1)
733 STOP 'ABNORMAL END: S/R INI_PARMS'
734 ELSE
735 _BEGIN_MASTER( myThid )
736 IF (readBinaryPrec.EQ.precFloat32) THEN
737 OPEN(37,FILE=delYfile,STATUS='OLD',FORM='UNFORMATTED',
738 & ACCESS='DIRECT',RECL=WORDLENGTH*Ny)
739 READ(37,rec=1) delY
740 #ifdef _BYTESWAPIO
741 call MDS_BYTESWAPR4( Ny, delY )
742 #endif
743 CLOSE(37)
744 ELSEIF (readBinaryPrec.EQ.precFloat64) THEN
745 OPEN(37,FILE=delYfile,STATUS='OLD',FORM='UNFORMATTED',
746 & ACCESS='DIRECT',RECL=WORDLENGTH*2*Ny)
747 READ(37,rec=1) delY
748 #ifdef _BYTESWAPIO
749 call MDS_BYTESWAPR8( Ny, delY )
750 #endif
751 CLOSE(37)
752 ENDIF
753 _END_MASTER(myThid)
754 ENDIF
755 ENDIF
756 IF ( dySpacing .NE. UNSET_RL ) THEN
757 DO i=1,Ny
758 delY(i) = dySpacing
759 ENDDO
760 ENDIF
761 C
762 IF ( rSphere .NE. 0 ) THEN
763 recip_rSphere = 1.D0/rSphere
764 ELSE
765 recip_rSphere = 0.
766 ENDIF
767 C-- Check for conflicting grid definitions.
768 goptCount = 0
769 IF ( usingCartesianGrid ) goptCount = goptCount+1
770 IF ( usingSphericalPolarGrid ) goptCount = goptCount+1
771 IF ( usingCurvilinearGrid ) goptCount = goptCount+1
772 IF ( goptCount .GT. 1 ) THEN
773 WRITE(msgBuf,'(A)')
774 & 'S/R INI_PARMS: More than one coordinate system requested'
775 CALL PRINT_ERROR( msgBuf , myThid)
776 STOP 'ABNORMAL END: S/R INI_PARMS'
777 ENDIF
778 IF ( goptCount .LT. 1 ) THEN
779 WRITE(msgBuf,'(A)')
780 & 'S/R INI_PARMS: No coordinate system requested'
781 CALL PRINT_ERROR( msgBuf , myThid)
782 STOP 'ABNORMAL END: S/R INI_PARMS'
783 ENDIF
784 C-- Make metric term settings consistent with underlying grid.
785 IF ( usingCartesianGrid ) THEN
786 usingSphericalPolarMterms = .FALSE.
787 metricTerms = .FALSE.
788 useNHMTerms = .FALSE.
789 mTFacMom = 0.
790 useBetaPlaneF = .TRUE.
791 ENDIF
792 IF ( usingSphericalPolarGrid ) THEN
793 useConstantF = .FALSE.
794 useBetaPlaneF = .FALSE.
795 useSphereF = .TRUE.
796 usingSphericalPolarMterms = metricTerms
797 ENDIF
798 IF ( usingCurvilinearGrid ) THEN
799 useSphereF = .TRUE.
800 metricTerms = .FALSE.
801 useNHMTerms = .FALSE.
802 ENDIF
803 C-- set cell Center depth and put Interface at the middle between 2 C
804 setCenterDr = .FALSE.
805 IF (delRc(1).NE.UNSET_RL) setCenterDr=.TRUE.
806 DO K=1,Nr+1
807 IF (delRc(K).EQ.UNSET_RL) setCenterDr = .FALSE.
808 ENDDO
809 C-- p, z, r coord parameters
810 DO K = 1, Nr
811 IF ( delZ(K) .NE. UNSET_RL ) zCoordInputData = .TRUE.
812 IF ( delP(K) .NE. UNSET_RL ) pCoordInputData = .TRUE.
813 IF ( delR(K) .NE. UNSET_RL ) rCoordInputData = .TRUE.
814 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delZ(K)
815 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delP(K)
816 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delRDefault(K)
817 IF (.NOT.setCenterDr .AND. delR(K).EQ.delRDefault(K) ) THEN
818 WRITE(msgBuf,'(A,I4)')
819 & 'S/R INI_PARMS: No value for delZ/delP/delR at K = ',K
820 CALL PRINT_ERROR( msgBuf , 1)
821 STOP 'ABNORMAL END: S/R INI_PARMS'
822 ELSEIF ( setCenterDr .AND. delR(K).NE.delRDefault(K) ) THEN
823 WRITE(msgBuf,'(2A,I4)') 'S/R INI_PARMS:',
824 & ' Cannot specify both delRc and delZ/delP/delR at K=', K
825 CALL PRINT_ERROR( msgBuf , 1)
826 STOP 'ABNORMAL END: S/R INI_PARMS'
827 ENDIF
828 ENDDO
829 C Check for multiple coordinate systems
830 CoordsSet = 0
831 IF ( zCoordInputData ) coordsSet = coordsSet + 1
832 IF ( pCoordInputData ) coordsSet = coordsSet + 1
833 IF ( rCoordInputData ) coordsSet = coordsSet + 1
834 IF ( coordsSet .GT. 1 ) THEN
835 WRITE(msgBuf,'(A)')
836 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
837 CALL PRINT_ERROR( msgBuf , myThid)
838 STOP 'ABNORMAL END: S/R INI_PARMS'
839 ENDIF
840
841 C-- Input files
842 READ(UNIT=iUnit,NML=PARM05) !,IOSTAT=errIO)
843 IF ( errIO .LT. 0 ) THEN
844 WRITE(msgBuf,'(A)')
845 & 'Error reading numerical model '
846 CALL PRINT_ERROR( msgBuf , 1)
847 WRITE(msgBuf,'(A)')
848 & 'parameter file "data"'
849 CALL PRINT_ERROR( msgBuf , 1)
850 WRITE(msgBuf,'(A)')
851 & 'Problem in namelist PARM05'
852 CALL PRINT_ERROR( msgBuf , 1)
853 CALL MODELDATA_EXAMPLE( myThid )
854 STOP 'ABNORMAL END: S/R INI_PARMS'
855 ELSE
856 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM05 : OK'
857 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
858 & SQUEEZE_RIGHT , 1)
859 ENDIF
860 C o If pLoadFile is set, then we should make sure the corresponing
861 C code is being compiled
862 #ifndef ATMOSPHERIC_LOADING
863 IF (pLoadFile.NE.' ') THEN
864 WRITE(msgBuf,'(A)')
865 & 'S/R INI_PARMS: pLoadFile is set but you have not'
866 CALL PRINT_ERROR( msgBuf , 1)
867 WRITE(msgBuf,'(A)')
868 & 'compiled the model with the pressure loading code.'
869 CALL PRINT_ERROR( msgBuf , 1)
870 WRITE(msgBuf,'(A,A)')
871 & 'Re-compile with: #define ATMOSPHERIC_LOADING',
872 & ' or -DATMOSPHERIC_LOADING'
873 CALL PRINT_ERROR( msgBuf , 1)
874 STOP 'ABNORMAL END: S/R INI_PARMS'
875 ENDIF
876 #endif
877
878 C
879 C-- Set factors required for mixing pressure and meters as vertical coordinate.
880 C rkFac is a "sign" parameter which is used where the orientation of the vertical
881 C coordinate (pressure or meters) relative to the vertical index (K) is important.
882 C rkFac = 1 applies when K and the coordinate are in the opposite sense.
883 C rkFac = -1 applies when K and the coordinate are in the same sense.
884 C horiVertRatio is a parameter that maps horizontal units to vertical units.
885 C It is used in certain special cases where lateral and vertical terms are
886 C being combined and a single frame of reference is needed.
887 IF ( zCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
888 rkFac = 1.D0
889 horiVertRatio = 1.D0
890 ENDIF
891 IF ( pCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
892 rkFac = -1.D0
893 horiVertRatio = Gravity * rhoConst
894 ENDIF
895 IF ( rCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
896 rkFac = 1.D0
897 horiVertRatio = 1.D0
898 ENDIF
899 convertEmP2rUnit = 1. _d 0
900 IF (buoyancyRelation.EQ.'ATMOSPHERIC')
901 & horiVertRatio = Gravity * rhoConst
902 IF (buoyancyRelation.EQ.'OCEANICP') THEN
903 horiVertRatio = Gravity * rhoConst
904 convertEmP2rUnit = Gravity * rhoConstFresh
905 ENDIF
906 IF ( rkFac .EQ. UNSET_RS ) rkFac=rkFacDefault
907 recip_rkFac = 1.D0 / rkFac
908 recip_horiVertRatio = 1./horiVertRatio
909 IF ( zCoordInputData ) usingZCoords = .TRUE.
910 IF ( pCoordInputData ) usingPCoords = .TRUE.
911
912 C
913 CLOSE(iUnit)
914
915 C-- Check whether any retired parameters were found.
916 C-- Stop if they were
917 IF ( nRetired .GT. 0 ) THEN
918 WRITE(msgBuf,'(A)')
919 & 'Error reading '
920 CALL PRINT_ERROR( msgBuf , 1)
921 WRITE(msgBuf,'(A)')
922 & 'parameter file "data"'
923 CALL PRINT_ERROR( msgBuf , 1)
924 WRITE(msgBuf,'(A)')
925 & 'some out of date parameters were found in the namelist'
926 CALL PRINT_ERROR( msgBuf , 1)
927 STOP 'ABNORMAL END: S/R INI_PARMS'
928 ENDIF
929
930 _END_MASTER(myThid)
931
932 C-- Everyone else must wait for the parameters to be loaded
933 _BARRIER
934 C
935 RETURN
936 END
937

  ViewVC Help
Powered by ViewVC 1.1.22