/[MITgcm]/MITgcm_contrib/osse/codemod/ini_parms.F
ViewVC logotype

Contents of /MITgcm_contrib/osse/codemod/ini_parms.F

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


Revision 1.1 - (show annotations) (download)
Tue Jun 22 19:44:40 2004 UTC (21 years, 1 month ago) by afe
Branch: MAIN
attempts to get tank config working with current checkpoint

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

  ViewVC Help
Powered by ViewVC 1.1.22