/[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.85 - (show annotations) (download)
Mon Oct 7 16:20:39 2002 UTC (21 years, 7 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint46l_post, checkpoint46l_pre, checkpoint46j_post, checkpoint46k_post
Changes since 1.84: +8 -5 lines
* update timeave pkg for wVel diagnostic ; put convertEmP2rUnit in PARAMS.h

1 C $Header: /u/gcmpack/MITgcm/model/src/ini_parms.F,v 1.84 2002/09/25 19:36:50 mlosch 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, globalFiles,
125 & allowFreezing, ivdc_kappa,
126 & bottomDragLinear,bottomDragQuadratic,
127 & usePickupBeforeC35, debugMode,
128 & readPickupWithTracer, writePickupWithTracer,
129 & tempAdvScheme, saltAdvScheme, tracerAdvScheme,
130 & multiDimAdvection, useEnergyConservingCoriolis,
131 & useJamartWetPoints,
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, gravitySign,
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-- Advection and Forcing for Temp and salt on/off flags
332 tempAdvection = tempStepping .AND. tempAdvection
333 tempForcing = tempStepping .AND. tempForcing
334 saltAdvection = saltStepping .AND. saltAdvection
335 saltForcing = saltStepping .AND. saltForcing
336 C-- z,p,r coord input switching.
337 IF ( viscAz .NE. UNSET_RL ) zCoordInputData = .TRUE.
338 IF ( viscAp .NE. UNSET_RL ) pCoordInputData = .TRUE.
339 IF ( viscAr .NE. UNSET_RL ) rCoordInputData = .TRUE.
340 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAz
341 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAp
342 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscArDefault
343
344 IF ( diffKzT .NE. UNSET_RL ) zCoordInputData = .TRUE.
345 IF ( diffKpT .NE. UNSET_RL ) pCoordInputData = .TRUE.
346 IF ( diffKrT .NE. UNSET_RL ) rCoordInputData = .TRUE.
347 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKzT
348 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKpT
349 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKrTDefault
350
351 IF ( diffKzS .NE. UNSET_RL ) zCoordInputData = .TRUE.
352 IF ( diffKpS .NE. UNSET_RL ) pCoordInputData = .TRUE.
353 IF ( diffKrS .NE. UNSET_RL ) rCoordInputData = .TRUE.
354 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKzS
355 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKpS
356 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKrSDefault
357
358 IF ( hFacMinDz .NE. UNSET_RL ) zCoordInputData = .TRUE.
359 IF ( hFacMinDp .NE. UNSET_RL ) pCoordInputData = .TRUE.
360 IF ( hFacMinDr .NE. UNSET_RL ) rCoordInputData = .TRUE.
361 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDz
362 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDp
363 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDrDefault
364
365 IF (convertFW2Salt.EQ.UNSET_RL) THEN
366 convertFW2Salt = 35.
367 IF (useRealFreshWaterFlux) convertFW2Salt=-1
368 ENDIF
369
370 IF ( ivdc_kappa .NE. 0. .AND. .NOT. implicitDiffusion ) THEN
371 WRITE(msgBuf,'(A,A)')
372 & 'S/R INI_PARMS: To use ivdc_kappa you must enable implicit',
373 & ' vertical diffusion.'
374 CALL PRINT_ERROR( msgBuf , myThid)
375 STOP 'ABNORMAL END: S/R INI_PARMS'
376 ENDIF
377
378 coordsSet = 0
379 IF ( zCoordInputData ) coordsSet = coordsSet + 1
380 IF ( pCoordInputData ) coordsSet = coordsSet + 1
381 IF ( rCoordInputData ) coordsSet = coordsSet + 1
382 IF ( coordsSet .GT. 1 ) THEN
383 WRITE(msgBuf,'(A)')
384 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
385 CALL PRINT_ERROR( msgBuf , myThid)
386 STOP 'ABNORMAL END: S/R INI_PARMS'
387 ENDIF
388 IF ( rhoConst .LE. 0. ) THEN
389 WRITE(msgBuf,'(A)')
390 & 'S/R INI_PARMS: rhoConst must be greater than 0.'
391 CALL PRINT_ERROR( msgBuf , myThid)
392 STOP 'ABNORMAL END: S/R INI_PARMS'
393 ELSE
394 recip_rhoConst = 1.D0 / rhoConst
395 ENDIF
396 IF ( rhoNil .LE. 0. ) THEN
397 WRITE(msgBuf,'(A)')
398 & 'S/R INI_PARMS: rhoNil must be greater than 0.'
399 CALL PRINT_ERROR( msgBuf , myThid)
400 STOP 'ABNORMAL END: S/R INI_PARMS'
401 ELSE
402 recip_rhoNil = 1.D0 / rhoNil
403 ENDIF
404 IF ( HeatCapacity_Cp .LE. 0. ) THEN
405 WRITE(msgBuf,'(A)')
406 & 'S/R INI_PARMS: HeatCapacity_Cp must be greater than 0.'
407 CALL PRINT_ERROR( msgBuf , myThid)
408 STOP 'ABNORMAL END: S/R INI_PARMS'
409 ELSE
410 recip_Cp = 1.D0 / HeatCapacity_Cp
411 ENDIF
412 IF ( gravity .LE. 0. ) THEN
413 WRITE(msgBuf,'(A)')
414 & 'S/R INI_PARMS: gravity must be greater than 0.'
415 CALL PRINT_ERROR( msgBuf , myThid)
416 STOP 'ABNORMAL END: S/R INI_PARMS'
417 ELSE
418 recip_gravity = 1.D0 / gravity
419 ENDIF
420 C Set globalFiles flag for READ_WRITE_FLD package
421 CALL SET_WRITE_GLOBAL_FLD( globalFiles )
422 C Set globalFiles flag for READ_WRITE_REC package
423 CALL SET_WRITE_GLOBAL_REC( globalFiles )
424 C Set globalFiles flag for READ_WRITE_REC package
425 CALL SET_WRITE_GLOBAL_PICKUP( globalFiles )
426
427 C Check for retired parameters still being used
428 nRetired = 0
429 IF ( zonal_filt_lat .NE. UNSET_RL ) THEN
430 nRetired = nRetired+1
431 WRITE(msgBuf,'(A,A)')
432 & 'S/R INI_PARMS: Paramater "zonal_filt_lat" is',
433 & ' no longer allowed in file "data".'
434 CALL PRINT_ERROR( msgBuf , myThid)
435 WRITE(msgBuf,'(A,A)')
436 & 'S/R INI_PARMS: Paramater "zonal_filt_lat" is',
437 & ' now read from file "data.zonfilt".'
438 CALL PRINT_ERROR( msgBuf , myThid)
439 ENDIF
440
441 C-- Elliptic solver parameters
442 READ(UNIT=iUnit,NML=PARM02) !,IOSTAT=errIO)
443 IF ( errIO .LT. 0 ) THEN
444 WRITE(msgBuf,'(A)')
445 & 'S/R INI_PARMS'
446 CALL PRINT_ERROR( msgBuf , 1)
447 WRITE(msgBuf,'(A)')
448 & 'Error reading numerical model '
449 CALL PRINT_ERROR( msgBuf , 1)
450 WRITE(msgBuf,'(A)')
451 & 'parameter file "data".'
452 CALL PRINT_ERROR( msgBuf , 1)
453 WRITE(msgBuf,'(A)')
454 & 'Problem in namelist PARM02'
455 CALL PRINT_ERROR( msgBuf , 1)
456 CALL MODELDATA_EXAMPLE( myThid )
457 STOP 'ABNORMAL END: S/R INI_PARMS'
458 ELSE
459 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM02 : OK'
460 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
461 & SQUEEZE_RIGHT , 1)
462 ENDIF
463
464 C-- Time stepping parameters
465 rCD = -1.D0
466 READ(UNIT=iUnit,NML=PARM03) !,IOSTAT=errIO)
467 IF ( errIO .LT. 0 ) THEN
468 WRITE(msgBuf,'(A)')
469 & 'S/R INI_PARMS'
470 CALL PRINT_ERROR( msgBuf , 1)
471 WRITE(msgBuf,'(A)')
472 & 'Error reading numerical model '
473 CALL PRINT_ERROR( msgBuf , 1)
474 WRITE(msgBuf,'(A)')
475 & 'parameter file "data"'
476 CALL PRINT_ERROR( msgBuf , 1)
477 WRITE(msgBuf,'(A)')
478 & 'Problem in namelist PARM03'
479 CALL PRINT_ERROR( msgBuf , 1)
480 CALL MODELDATA_EXAMPLE( myThid )
481 STOP 'ABNORMAL END: S/R INI_PARMS'
482 ELSE
483 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM03 : OK'
484 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
485 & SQUEEZE_RIGHT , 1)
486 ENDIF
487 C Process "timestepping" params
488 C o Time step size
489 IF ( deltaT .EQ. 0. ) deltaT = deltaTmom
490 IF ( deltaT .EQ. 0. ) deltaT = deltaTtracer
491 IF ( deltaTmom .EQ. 0. ) deltaTmom = deltaT
492 IF ( deltaTtracer .EQ. 0. ) deltaTtracer = deltaT
493 IF ( deltaTClock .EQ. 0. ) deltaTClock = deltaT
494 C Note that this line should set deltaTreesurf=deltaTmom
495 C but this would change a lot of existing set-ups so we are
496 C obliged to set the default inappropriately.
497 C Be advised that when using asynchronous time stepping
498 C it is better to set deltaTreesurf=deltaTtracer
499 IF ( deltaTfreesurf .EQ. 0. ) deltaTfreesurf = deltaTmom
500 IF ( periodicExternalForcing ) THEN
501 IF ( externForcingCycle*externForcingPeriod .EQ. 0. ) THEN
502 WRITE(msgBuf,'(A)')
503 & 'S/R INI_PARMS: externForcingCycle,externForcingPeriod =0'
504 CALL PRINT_ERROR( msgBuf , 1)
505 STOP 'ABNORMAL END: S/R INI_PARMS'
506 ENDIF
507 IF ( INT(externForcingCycle/externForcingPeriod) .NE.
508 & externForcingCycle/externForcingPeriod ) THEN
509 WRITE(msgBuf,'(A)')
510 & 'S/R INI_PARMS: externForcingCycle <> N*externForcingPeriod'
511 CALL PRINT_ERROR( msgBuf , 1)
512 STOP 'ABNORMAL END: S/R INI_PARMS'
513 ENDIF
514 IF ( externForcingCycle.le.externForcingPeriod ) THEN
515 WRITE(msgBuf,'(A)')
516 & 'S/R INI_PARMS: externForcingCycle < externForcingPeriod'
517 CALL PRINT_ERROR( msgBuf , 1)
518 STOP 'ABNORMAL END: S/R INI_PARMS'
519 ENDIF
520 IF ( externForcingPeriod.lt.deltaTclock ) THEN
521 WRITE(msgBuf,'(A)')
522 & 'S/R INI_PARMS: externForcingPeriod < deltaTclock'
523 CALL PRINT_ERROR( msgBuf , 1)
524 STOP 'ABNORMAL END: S/R INI_PARMS'
525 ENDIF
526 ENDIF
527 C o Convection frequency
528 IF ( cAdjFreq .LT. 0. ) THEN
529 cAdjFreq = deltaTClock
530 ENDIF
531 IF ( ivdc_kappa .NE. 0. .AND. cAdjFreq .NE. 0. ) THEN
532 WRITE(msgBuf,'(A,A)')
533 & 'S/R INI_PARMS: You have enabled both ivdc_kappa and',
534 & ' convective adjustment.'
535 CALL PRINT_ERROR( msgBuf , myThid)
536 STOP 'ABNORMAL END: S/R INI_PARMS'
537 ENDIF
538 C o CD coupling
539 IF ( tauCD .EQ. 0.D0 ) THEN
540 tauCD = deltaTmom
541 ENDIF
542 IF ( rCD .LT. 0. ) THEN
543 rCD = 1. - deltaTMom/tauCD
544 ENDIF
545 C o Temperature climatology relaxation time scale
546 IF ( tauThetaClimRelax .EQ. 0.D0 ) THEN
547 doThetaClimRelax = .FALSE.
548 lambdaThetaClimRelax = 0.D0
549 ELSE
550 doThetaClimRelax = .TRUE.
551 lambdaThetaClimRelax = 1./tauThetaClimRelax
552 ENDIF
553 C o Salinity climatology relaxation time scale
554 IF ( tauSaltClimRelax .EQ. 0.D0 ) THEN
555 doSaltClimRelax = .FALSE.
556 lambdaSaltClimRelax = 0.D0
557 ELSE
558 doSaltClimRelax = .TRUE.
559 lambdaSaltClimRelax = 1./tauSaltClimRelax
560 ENDIF
561 C o Tracer 1 climatology relaxation time scale
562 IF ( tauTr1ClimRelax .EQ. 0.D0 ) THEN
563 doTr1ClimRelax = .FALSE.
564 lambdaTr1ClimRelax = 0.D0
565 ELSE
566 doTr1ClimRelax = .TRUE.
567 lambdaTr1ClimRelax = 1./tauTr1ClimRelax
568 ENDIF
569
570 C o Start time
571 IF ( nIter0 .NE. 0 .AND. startTime .EQ. 0. )
572 & startTime = deltaTClock*float(nIter0)
573 C o nIter0
574 IF ( nIter0 .EQ. 0 .AND. startTime .NE. 0. )
575 & nIter0 = INT( startTime/deltaTClock )
576
577 C o nTimeSteps 1
578 IF ( nTimeSteps .EQ. 0 .AND. nEndIter .NE. 0 )
579 & nTimeSteps = nEndIter-nIter0
580 C o nTimeSteps 2
581 IF ( nTimeSteps .EQ. 0 .AND. endTime .NE. 0. )
582 & nTimeSteps = int(0.5+(endTime-startTime)/deltaTclock)
583 C o nEndIter 1
584 IF ( nEndIter .EQ. 0 .AND. nTimeSteps .NE. 0 )
585 & nEndIter = nIter0+nTimeSteps
586 C o nEndIter 2
587 IF ( nEndIter .EQ. 0 .AND. endTime .NE. 0. )
588 & nEndIter = int(0.5+endTime/deltaTclock)
589 C o End Time 1
590 IF ( endTime .EQ. 0. .AND. nTimeSteps .NE. 0 )
591 & endTime = startTime + deltaTClock*float(nTimeSteps)
592 C o End Time 2
593 IF ( endTime .EQ. 0. .AND. nEndIter .NE. 0 )
594 & endTime = deltaTClock*float(nEndIter)
595
596 C o Consistent?
597 IF ( nEndIter .NE. nIter0+nTimeSteps ) THEN
598 WRITE(msgBuf,'(A)')
599 & 'S/R INI_PARMS: nIter0, nTimeSteps and nEndIter are inconsistent'
600 CALL PRINT_ERROR( msgBuf , 1)
601 WRITE(msgBuf,'(A)')
602 & 'S/R INI_PARMS: Perhaps more than two were set at once'
603 CALL PRINT_ERROR( msgBuf , 1)
604 STOP 'ABNORMAL END: S/R INI_PARMS'
605 ENDIF
606 IF ( nTimeSteps .NE. int(0.5+(endTime-startTime)/deltaTClock) )
607 & THEN
608 WRITE(msgBuf,'(A)')
609 & 'S/R INI_PARMS: both endTime and nTimeSteps have been set'
610 CALL PRINT_ERROR( msgBuf , 1)
611 WRITE(msgBuf,'(A)')
612 & 'S/R INI_PARMS: but are inconsistent'
613 CALL PRINT_ERROR( msgBuf , 1)
614 STOP 'ABNORMAL END: S/R INI_PARMS'
615 ENDIF
616
617 C o Monitor (should also add CPP flag for monitor?)
618 IF (monitorFreq.LT.0.) THEN
619 monitorFreq=0.
620 IF (dumpFreq.NE.0.) monitorFreq=dumpFreq
621 IF (diagFreq.NE.0..AND.diagFreq.LT.monitorFreq)
622 & monitorFreq=diagFreq
623 IF (taveFreq.NE.0..AND.taveFreq.LT.monitorFreq)
624 & monitorFreq=taveFreq
625 IF (chkPtFreq.NE.0..AND.chkPtFreq.LT.monitorFreq)
626 & monitorFreq=chkPtFreq
627 IF (pChkPtFreq.NE.0..AND.pChkPtFreq.LT.monitorFreq)
628 & monitorFreq=pChkPtFreq
629 IF (monitorFreq.EQ.0.) monitorFreq=deltaTclock
630 ENDIF
631
632 C o If taveFreq is finite, then we must make sure the diagnostics
633 C code is being compiled
634 #ifndef ALLOW_TIMEAVE
635 IF (taveFreq.NE.0.) THEN
636 WRITE(msgBuf,'(A)')
637 & 'S/R INI_PARMS: taveFreq <> 0 but you have'
638 CALL PRINT_ERROR( msgBuf , 1)
639 WRITE(msgBuf,'(A)')
640 & 'not compiled the model with the diagnostics routines.'
641 CALL PRINT_ERROR( msgBuf , 1)
642 WRITE(msgBuf,'(A,A)')
643 & 'Re-compile with: #define ALLOW_TIMEAVE',
644 & ' or -DALLOW_TIMEAVE'
645 CALL PRINT_ERROR( msgBuf , 1)
646 STOP 'ABNORMAL END: S/R INI_PARMS'
647 ENDIF
648 #endif
649
650 C-- Grid parameters
651 C In cartesian coords distances are in metres
652 rkFac = UNSET_RS
653 DO K =1,Nr
654 delZ(K) = UNSET_RL
655 delP(K) = UNSET_RL
656 delR(K) = UNSET_RL
657 ENDDO
658 C In spherical polar distances are in degrees
659 recip_rSphere = 1.D0/rSphere
660 dxSpacing = UNSET_RL
661 dySpacing = UNSET_RL
662 delXfile = ' '
663 delYfile = ' '
664 READ(UNIT=iUnit,NML=PARM04,IOSTAT=errIO)
665 IF ( errIO .LT. 0 ) THEN
666 WRITE(msgBuf,'(A)')
667 & 'S/R INI_PARMS'
668 CALL PRINT_ERROR( msgBuf , 1)
669 WRITE(msgBuf,'(A)')
670 & 'Error reading numerical model '
671 CALL PRINT_ERROR( msgBuf , 1)
672 WRITE(msgBuf,'(A)')
673 & 'parameter file "data"'
674 CALL PRINT_ERROR( msgBuf , 1)
675 WRITE(msgBuf,'(A)')
676 & 'Problem in namelist PARM04'
677 CALL PRINT_ERROR( msgBuf , 1)
678 CALL MODELDATA_EXAMPLE( myThid )
679 STOP 'ABNORMAL END: S/R INI_PARMS'
680 ELSE
681 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM04 : OK'
682 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
683 & SQUEEZE_RIGHT , 1)
684 ENDIF
685
686 C X coordinate
687 IF ( delXfile .NE. ' ' ) THEN
688 IF ( delX(1) .NE. UNSET_RL .OR. dxSpacing .NE. UNSET_RL ) THEN
689 WRITE(msgBuf,'(A,A)') 'Too many specifications for delX:',
690 & 'Specify only one of delX, dxSpacing or delXfile'
691 CALL PRINT_ERROR( msgBuf , 1)
692 STOP 'ABNORMAL END: S/R INI_PARMS'
693 ELSE
694 _BEGIN_MASTER( myThid )
695 IF (readBinaryPrec.EQ.precFloat32) THEN
696 OPEN(37,FILE=delXfile,STATUS='OLD',FORM='UNFORMATTED',
697 & ACCESS='DIRECT',RECL=WORDLENGTH*Nx)
698 READ(37,rec=1) delX
699 #ifdef _BYTESWAPIO
700 call MDS_BYTESWAPR4( Nx, delX )
701 #endif
702 CLOSE(37)
703 ELSEIF (readBinaryPrec.EQ.precFloat64) THEN
704 OPEN(37,FILE=delXfile,STATUS='OLD',FORM='UNFORMATTED',
705 & ACCESS='DIRECT',RECL=WORDLENGTH*2*Nx)
706 READ(37,rec=1) delX
707 #ifdef _BYTESWAPIO
708 call MDS_BYTESWAPR8( Nx, delX )
709 #endif
710 CLOSE(37)
711 ENDIF
712 _END_MASTER(myThid)
713 ENDIF
714 ENDIF
715 IF ( dxSpacing .NE. UNSET_RL ) THEN
716 DO i=1,Nx
717 delX(i) = dxSpacing
718 ENDDO
719 ENDIF
720 C Y coordinate
721 IF ( delYfile .NE. ' ' ) THEN
722 IF ( delY(1) .NE. UNSET_RL .OR. dySpacing .NE. UNSET_RL ) THEN
723 WRITE(msgBuf,'(A,A)') 'Too many specifications for delY:',
724 & 'Specify only one of delY, dySpacing or delYfile'
725 CALL PRINT_ERROR( msgBuf , 1)
726 STOP 'ABNORMAL END: S/R INI_PARMS'
727 ELSE
728 _BEGIN_MASTER( myThid )
729 IF (readBinaryPrec.EQ.precFloat32) THEN
730 OPEN(37,FILE=delYfile,STATUS='OLD',FORM='UNFORMATTED',
731 & ACCESS='DIRECT',RECL=WORDLENGTH*Ny)
732 READ(37,rec=1) delY
733 #ifdef _BYTESWAPIO
734 call MDS_BYTESWAPR4( Ny, delY )
735 #endif
736 CLOSE(37)
737 ELSEIF (readBinaryPrec.EQ.precFloat64) THEN
738 OPEN(37,FILE=delYfile,STATUS='OLD',FORM='UNFORMATTED',
739 & ACCESS='DIRECT',RECL=WORDLENGTH*2*Ny)
740 READ(37,rec=1) delY
741 #ifdef _BYTESWAPIO
742 call MDS_BYTESWAPR8( Ny, delY )
743 #endif
744 CLOSE(37)
745 ENDIF
746 _END_MASTER(myThid)
747 ENDIF
748 ENDIF
749 IF ( dySpacing .NE. UNSET_RL ) THEN
750 DO i=1,Ny
751 delY(i) = dySpacing
752 ENDDO
753 ENDIF
754 C
755 IF ( rSphere .NE. 0 ) THEN
756 recip_rSphere = 1.D0/rSphere
757 ELSE
758 recip_rSphere = 0.
759 ENDIF
760 C-- Check for conflicting grid definitions.
761 goptCount = 0
762 IF ( usingCartesianGrid ) goptCount = goptCount+1
763 IF ( usingSphericalPolarGrid ) goptCount = goptCount+1
764 IF ( usingCurvilinearGrid ) goptCount = goptCount+1
765 IF ( goptCount .GT. 1 ) THEN
766 WRITE(msgBuf,'(A)')
767 & 'S/R INI_PARMS: More than one coordinate system requested'
768 CALL PRINT_ERROR( msgBuf , myThid)
769 STOP 'ABNORMAL END: S/R INI_PARMS'
770 ENDIF
771 IF ( goptCount .LT. 1 ) THEN
772 WRITE(msgBuf,'(A)')
773 & 'S/R INI_PARMS: No coordinate system requested'
774 CALL PRINT_ERROR( msgBuf , myThid)
775 STOP 'ABNORMAL END: S/R INI_PARMS'
776 ENDIF
777 C-- Make metric term settings consistent with underlying grid.
778 IF ( usingCartesianGrid ) THEN
779 usingSphericalPolarMterms = .FALSE.
780 metricTerms = .FALSE.
781 mTFacMom = 0.
782 useBetaPlaneF = .TRUE.
783 ENDIF
784 IF ( usingSphericalPolarGrid ) THEN
785 useConstantF = .FALSE.
786 useBetaPlaneF = .FALSE.
787 useSphereF = .TRUE.
788 usingSphericalPolarMterms = metricTerms
789 ENDIF
790 IF ( usingCurvilinearGrid ) THEN
791 useSphereF = .TRUE.
792 ENDIF
793 C-- set cell Center depth and put Interface at the middle between 2 C
794 setCenterDr = .FALSE.
795 IF (delRc(1).NE.UNSET_RL) setCenterDr=.TRUE.
796 DO K=1,Nr+1
797 IF (delRc(K).EQ.UNSET_RL) setCenterDr = .FALSE.
798 ENDDO
799 C-- p, z, r coord parameters
800 DO K = 1, Nr
801 IF ( delZ(K) .NE. UNSET_RL ) zCoordInputData = .TRUE.
802 IF ( delP(K) .NE. UNSET_RL ) pCoordInputData = .TRUE.
803 IF ( delR(K) .NE. UNSET_RL ) rCoordInputData = .TRUE.
804 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delZ(K)
805 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delP(K)
806 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delRDefault(K)
807 IF (.NOT.setCenterDr .AND. delR(K).EQ.delRDefault(K) ) THEN
808 WRITE(msgBuf,'(A,I4)')
809 & 'S/R INI_PARMS: No value for delZ/delP/delR at K = ',K
810 CALL PRINT_ERROR( msgBuf , 1)
811 STOP 'ABNORMAL END: S/R INI_PARMS'
812 ELSEIF ( setCenterDr .AND. delR(K).NE.delRDefault(K) ) THEN
813 WRITE(msgBuf,'(2A,I4)') 'S/R INI_PARMS:',
814 & ' Cannot specify both delRc and delZ/delP/delR at K=', K
815 CALL PRINT_ERROR( msgBuf , 1)
816 STOP 'ABNORMAL END: S/R INI_PARMS'
817 ENDIF
818 ENDDO
819 C Check for multiple coordinate systems
820 CoordsSet = 0
821 IF ( zCoordInputData ) coordsSet = coordsSet + 1
822 IF ( pCoordInputData ) coordsSet = coordsSet + 1
823 IF ( rCoordInputData ) coordsSet = coordsSet + 1
824 IF ( coordsSet .GT. 1 ) THEN
825 WRITE(msgBuf,'(A)')
826 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
827 CALL PRINT_ERROR( msgBuf , myThid)
828 STOP 'ABNORMAL END: S/R INI_PARMS'
829 ENDIF
830
831 C-- Input files
832 READ(UNIT=iUnit,NML=PARM05) !,IOSTAT=errIO)
833 IF ( errIO .LT. 0 ) THEN
834 WRITE(msgBuf,'(A)')
835 & 'Error reading numerical model '
836 CALL PRINT_ERROR( msgBuf , 1)
837 WRITE(msgBuf,'(A)')
838 & 'parameter file "data"'
839 CALL PRINT_ERROR( msgBuf , 1)
840 WRITE(msgBuf,'(A)')
841 & 'Problem in namelist PARM05'
842 CALL PRINT_ERROR( msgBuf , 1)
843 CALL MODELDATA_EXAMPLE( myThid )
844 STOP 'ABNORMAL END: S/R INI_PARMS'
845 ELSE
846 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM05 : OK'
847 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
848 & SQUEEZE_RIGHT , 1)
849 ENDIF
850 C o If pLoadFile is set, then we should make sure the corresponing
851 C code is being compiled
852 #ifndef ATMOSPHERIC_LOADING
853 IF (pLoadFile.NE.' ') THEN
854 WRITE(msgBuf,'(A)')
855 & 'S/R INI_PARMS: pLoadFile is set but you have not'
856 CALL PRINT_ERROR( msgBuf , 1)
857 WRITE(msgBuf,'(A)')
858 & 'compiled the model with the pressure loading code.'
859 CALL PRINT_ERROR( msgBuf , 1)
860 WRITE(msgBuf,'(A,A)')
861 & 'Re-compile with: #define ATMOSPHERIC_LOADING',
862 & ' or -DATMOSPHERIC_LOADING'
863 CALL PRINT_ERROR( msgBuf , 1)
864 STOP 'ABNORMAL END: S/R INI_PARMS'
865 ENDIF
866 #endif
867
868 C
869 C-- Set factors required for mixing pressure and meters as vertical coordinate.
870 C rkFac is a "sign" parameter which is used where the orientation of the vertical
871 C coordinate (pressure or meters) relative to the vertical index (K) is important.
872 C rkFac = 1 applies when K and the coordinate are in the opposite sense.
873 C rkFac = -1 applies when K and the coordinate are in the same sense.
874 C horiVertRatio is a parameter that maps horizontal units to vertical units.
875 C It is used in certain special cases where lateral and vertical terms are
876 C being combined and a single frame of reference is needed.
877 IF ( zCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
878 rkFac = 1.D0
879 horiVertRatio = 1.D0
880 ENDIF
881 IF ( pCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
882 rkFac = -1.D0
883 horiVertRatio = Gravity * rhoConst
884 ENDIF
885 IF ( rCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
886 rkFac = 1.D0
887 horiVertRatio = 1.D0
888 ENDIF
889 convertEmP2rUnit = 1. _d 0
890 IF (buoyancyRelation.EQ.'ATMOSPHERIC')
891 & horiVertRatio = Gravity * rhoConst
892 IF (buoyancyRelation.EQ.'OCEANICP') THEN
893 horiVertRatio = Gravity * rhoConst
894 convertEmP2rUnit = Gravity * rhoConstFresh
895 ENDIF
896 IF ( rkFac .EQ. UNSET_RS ) rkFac=rkFacDefault
897 recip_rkFac = 1.D0 / rkFac
898 recip_horiVertRatio = 1./horiVertRatio
899 IF ( zCoordInputData ) usingZCoords = .TRUE.
900 IF ( pCoordInputData ) usingPCoords = .TRUE.
901
902 C
903 CLOSE(iUnit)
904
905 C-- Check whether any retired parameters were found.
906 C-- Stop if they were
907 IF ( nRetired .GT. 0 ) THEN
908 WRITE(msgBuf,'(A)')
909 & 'Error reading '
910 CALL PRINT_ERROR( msgBuf , 1)
911 WRITE(msgBuf,'(A)')
912 & 'parameter file "data"'
913 CALL PRINT_ERROR( msgBuf , 1)
914 WRITE(msgBuf,'(A)')
915 & 'some out of date parameters were found in the namelist'
916 CALL PRINT_ERROR( msgBuf , 1)
917 STOP 'ABNORMAL END: S/R INI_PARMS'
918 ENDIF
919
920 _END_MASTER(myThid)
921
922 C-- Everyone else must wait for the parameters to be loaded
923 _BARRIER
924 C
925 RETURN
926 END
927

  ViewVC Help
Powered by ViewVC 1.1.22