/[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.174 - (show annotations) (download)
Tue Dec 13 16:41:47 2005 UTC (18 years, 5 months ago) by jmc
Branch: MAIN
Changes since 1.173: +2 -2 lines
add new parameter (implicitIntGravWave) even if the code is not yet ready;
(hard to maintain a different version of theses files during testing phase)

1 C $Header: /u/gcmpack/MITgcm/model/src/ini_parms.F,v 1.173 2005/10/24 22:40:28 jmc Exp $
2 C $Name: $
3
4 #include "PACKAGES_CONFIG.h"
5 #include "CPP_OPTIONS.h"
6
7 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
8 CBOP
9 C !ROUTINE: INI_PARMS
10 C !INTERFACE:
11 SUBROUTINE INI_PARMS( myThid )
12
13 C !DESCRIPTION:
14 C Routine to set model "parameters".
15
16 C !USES:
17 IMPLICIT NONE
18 #include "SIZE.h"
19 #include "EEPARAMS.h"
20 #include "PARAMS.h"
21 #include "GRID.h"
22 #include "EOS.h"
23
24 C !INPUT/OUTPUT PARAMETERS:
25 C myThid - Number of this instance of INI_PARMS
26 INTEGER myThid
27
28 C !LOCAL VARIABLES:
29 C dxSpacing, dySpacing - Default spacing in X and Y.
30 C Units are that of coordinate system
31 C i.e. cartesian => metres
32 C s. polar => degrees
33 C deltaTtracer :: Timestep for tracer equations ( s )
34 C goptCount - Used to count the nuber of grid options
35 C (only one is allowed! )
36 C msgBuf - Informational/error meesage buffer
37 C errIO - IO error flag
38 C iUnit - Work variable for IO unit number
39 C record - Work variable for IO buffer
40 C K, I, J - Loop counters
41 C xxxDefault - Default value for variable xxx
42 _RL dxSpacing
43 _RL dySpacing
44 _RL deltaTtracer
45 CHARACTER*(MAX_LEN_MBUF) msgBuf
46 CHARACTER*(MAX_LEN_PREC) record
47 INTEGER goptCount
48 INTEGER K, i, j, IL, iUnit
49 INTEGER errIO
50 INTEGER IFNBLNK
51 EXTERNAL IFNBLNK
52 INTEGER ILNBLNK
53 EXTERNAL ILNBLNK
54 C Default values for variables which have vertical coordinate system
55 C dependency.
56 _RL viscArDefault
57 _RL diffKrTDefault
58 _RL diffKrSDefault
59 _RL hFacMinDrDefault
60 _RL delRDefault(Nr)
61 _RS rkFacDefault
62 _RL tracerDefault
63 C zCoordInputData :: Variables used to select between different coordinate systems.
64 C pCoordInputData :: The vertical coordinate system in the rest of the model is
65 C rCoordInputData :: written in terms of r. In the model "data" file input data can
66 C coordsSet :: be interms of z, p or r.
67 C :: e.g. delZ or delP or delR for the vertical grid spacing.
68 C :: The following rules apply:
69 C :: All parameters must use the same vertical coordinate system.
70 C :: e.g. delZ and viscAz is legal but
71 C :: delZ and viscAr is an error.
72 C :: Similarly specifyinh delZ and delP is an error.
73 C :: zCoord..., pCoord..., rCoord... are used to flag when z, p or r are
74 C :: used. coordsSet counts how many vertical coordinate systems have been
75 C used to specify variables. coordsSet > 1 is an error.
76 C
77
78 LOGICAL zCoordInputData
79 LOGICAL pCoordInputData
80 LOGICAL rCoordInputData
81 INTEGER coordsSet
82 LOGICAL diffKrSet
83
84 C Variables which have vertical coordinate system dependency.
85 C delZ :: Vertical grid spacing ( m ).
86 C delP :: Vertical grid spacing ( Pa ).
87 C viscAz :: Eddy viscosity coeff. for mixing of
88 C momentum vertically ( m^2/s )
89 C viscAp :: Eddy viscosity coeff. for mixing of
90 C momentum vertically ( Pa^2/s )
91 C diffKzT :: Laplacian diffusion coeff. for mixing of
92 C heat vertically ( m^2/s )
93 C diffKpT :: Laplacian diffusion coeff. for mixing of
94 C heat vertically ( Pa^2/s )
95 C diffKzS :: Laplacian diffusion coeff. for mixing of
96 C salt vertically ( m^2/s )
97 C diffKpS :: Laplacian diffusion coeff. for mixing of
98 C salt vertically ( Pa^2/s )
99 _RL delZ(Nr)
100 _RL delP(Nr)
101 _RL viscAz
102 _RL viscAp
103 _RL diffKzT
104 _RL diffKpT
105 _RL diffKrT
106 _RL diffKzS
107 _RL diffKpS
108 _RL diffKrS
109
110 C Retired main data file parameters. Kept here to trap use of old data files.
111 C tracerAdvScheme :: tracer advection scheme (old passive tracer code)
112 C trac_EvPrRn :: tracer conc. in Rain & Evap (old passive tracer code)
113 C saltDiffusion :: diffusion of salinity on/off (flag not used)
114 C tempDiffusion :: diffusion of temperature on/off (flag not used)
115 C zonal_filt_lat :: Moved to package "zonal_filt"
116 C groundAtK1 :: put the surface(k=1) at the ground (replaced by usingPCoords)
117 C rkFac :: removed from namelist ; replaced by -rkSign
118 C viscAstrain :: replaced by standard viscosity coeff & useStrainTensionVisc
119 C viscAtension :: replaced by standard viscosity coeff & useStrainTensionVisc
120 C nRetired :: Counter used to trap gracefully namelists containing "retired"
121 C :: parameters. These are parameters that are either no-longer used
122 C or that have moved to a different input file and/or namelist.
123 C useAnisotropicViscAgridMax :: Changed to be default behavior. Can
124 C use old method by setting useAreaViscLength=.true.
125 LOGICAL tempDiffusion, saltDiffusion
126 INTEGER tracerAdvScheme
127 _RL trac_EvPrRn
128 _RL zonal_filt_lat
129 _RL rkFac
130 _RL viscAstrain, viscAtension
131 LOGICAL groundAtK1, useAnisotropicViscAGridMax
132 INTEGER nRetired
133
134 C-- Continuous equation parameters
135 NAMELIST /PARM01/
136 & gravitySign, nh_Am2,
137 & gravity, gBaro, rhonil, tAlpha, sBeta,
138 & f0, beta, omega, rotationPeriod,
139 & viscAh, viscAhW, viscAhMax,
140 & viscAhGrid, viscAhGridMax, viscAhGridMin,
141 & viscC2leith, viscC4leith,
142 & useFullLeith, useAnisotropicViscAGridMax, useStrainTensionVisc,
143 & useAreaViscLength,
144 & viscC2leithD, viscC4leithD, viscC2Smag, viscC4Smag,
145 & viscAhD, viscAhZ, viscA4D, viscA4Z,
146 & viscA4, viscA4W,
147 & viscA4Max, viscA4Grid, viscA4GridMax, viscA4GridMin,
148 & viscA4ReMax, viscAhReMax,
149 & viscAz, cosPower, viscAstrain, viscAtension,
150 & diffKhT, diffKzT, diffK4T,
151 & diffKhS, diffKzS, diffK4S,
152 & tRef, sRef, eosType, integr_GeoPot, selectFindRoSurf,
153 & atm_Cp, atm_Rd, atm_Rq,
154 & no_slip_sides, sideDragFactor,
155 & no_slip_bottom, bottomDragLinear, bottomDragQuadratic,
156 & momViscosity, momAdvection, momForcing, useCoriolis,
157 & useConstantF,
158 & momPressureForcing, metricTerms, vectorInvariantMomentum,
159 & tempDiffusion, tempAdvection, tempForcing,
160 & saltDiffusion, saltAdvection, saltForcing,
161 & implicSurfPress, implicDiv2DFlow,
162 & implicitFreeSurface, rigidLid, freeSurfFac, hFacMin, hFacMinDz,
163 & exactConserv,uniformLin_PhiSurf,nonlinFreeSurf,hFacInf,hFacSup,
164 & select_rStar,
165 & implicitIntGravWave, staggerTimeStep,
166 & tempStepping, saltStepping, momStepping,
167 & implicitDiffusion, implicitViscosity,
168 & tempImplVertAdv, saltImplVertAdv, momImplVertAdv,
169 & viscAr, diffKrT, diffKrS, diffKrNrT, diffKrNrS, hFacMinDr,
170 & viscAp, diffKpT, diffKpS, hFacMinDp,
171 & diffKrBL79surf, diffKrBL79deep, diffKrBL79scl, diffKrBL79Ho,
172 & rhoConst, rhoConstFresh, buoyancyRelation, HeatCapacity_Cp,
173 & writeBinaryPrec, readBinaryPrec, writeStatePrec,
174 & nonHydrostatic, quasiHydrostatic, globalFiles, useSingleCpuIO,
175 & allowFreezing, useOldFreezing, ivdc_kappa,
176 & usePickupBeforeC35, usePickupBeforeC54, debugMode, debugLevel,
177 & tempAdvScheme, tempVertAdvScheme,
178 & saltAdvScheme, saltVertAdvScheme, tracerAdvScheme,
179 & multiDimAdvection, useEnergyConservingCoriolis,
180 & useCDscheme, useJamartWetPoints, useJamartMomAdv, useNHMTerms,
181 & SadournyCoriolis, upwindVorticity, highOrderVorticity,
182 & useAbsVorticity, upwindShear, selectKEscheme,
183 & useRealFreshWaterFlux, convertFW2Salt,
184 & temp_EvPrRn, salt_EvPrRn, trac_EvPrRn,
185 & zonal_filt_lat,
186 & inAdExact
187
188 C-- Elliptic solver parameters
189 NAMELIST /PARM02/
190 & cg2dMaxIters, cg2dChkResFreq, cg2dTargetResidual,
191 & cg2dTargetResWunit, cg2dpcOffDFac, cg2dPreCondFreq,
192 & cg3dMaxIters, cg3dChkResFreq, cg3dTargetResidual
193
194 C-- Time stepping parammeters
195 NAMELIST /PARM03/
196 & nIter0, nTimeSteps, nEndIter, pickupSuff,
197 & deltaT, deltaTClock, deltaTmom,
198 & deltaTtracer, dTtracerLev, deltaTfreesurf,
199 & forcing_In_AB, doAB_onGtGs,
200 & abEps, alph_AB, beta_AB, startFromPickupAB2,
201 & tauCD, rCD,
202 & baseTime, startTime, endTime, chkPtFreq,
203 & dumpFreq, dumpInitAndLast, adjDumpFreq, taveFreq, tave_lastIter,
204 & diagFreq, monitorFreq, adjMonitorFreq, pChkPtFreq, cAdjFreq,
205 & outputTypesInclusive,
206 & tauThetaClimRelax, tauSaltClimRelax, latBandClimRelax,
207 & tauTr1ClimRelax,
208 & periodicExternalForcing, externForcingPeriod, externForcingCycle,
209 & calendarDumps
210
211 C-- Gridding parameters
212 NAMELIST /PARM04/
213 & usingCartesianGrid, dxSpacing, dySpacing, delX, delY, delZ,
214 & usingSphericalPolarGrid, phiMin, thetaMin, rSphere,
215 & usingCurvilinearGrid, usingCylindricalGrid,
216 & delP, delR, rkFac, Ro_SeaLevel, groundAtK1, delRc,
217 & delXFile, delYFile, horizGridFile
218
219 C-- Input files
220 NAMELIST /PARM05/
221 & bathyFile, topoFile,
222 & hydrogThetaFile, hydrogSaltFile,
223 & zonalWindFile, meridWindFile,
224 & thetaClimFile, saltClimFile,
225 & surfQfile, surfQnetFile, surfQswFile, EmPmRfile, saltFluxFile,
226 & lambdaThetaFile, lambdaSaltFile,
227 & uVelInitFile, vVelInitFile, pSurfInitFile,
228 & dQdTFile, ploadFile,tCylIn,tCylOut,
229 & eddyTauxFile, eddyTauyFile,
230 & mdsioLocalDir,
231 & the_run_name
232 CEOP
233
234 C
235 _BEGIN_MASTER(myThid)
236
237 C Defaults values for input parameters
238 CALL SET_DEFAULTS(
239 O viscArDefault, diffKrTDefault, diffKrSDefault,
240 O hFacMinDrDefault, delRdefault, rkFacDefault,
241 I myThid )
242
243 C-- Initialise "which vertical coordinate system used" flags.
244 zCoordInputData = .FALSE.
245 pCoordInputData = .FALSE.
246 rCoordInputData = .FALSE.
247 coordsSet = 0
248
249 C-- Initialise retired parameters to unlikely value
250 nRetired = 0
251 tempDiffusion = .TRUE.
252 saltDiffusion = .TRUE.
253 tracerAdvScheme = UNSET_I
254 trac_EvPrRn = UNSET_RL
255 zonal_filt_lat = UNSET_RL
256 gravitySign = UNSET_RL
257 rkFac = UNSET_RL
258 groundAtK1 = .FALSE.
259 viscAstrain = UNSET_RL
260 viscAtension = UNSET_RL
261 useAnisotropicViscAgridMax=.TRUE.
262
263 C-- Open the parameter file
264 #ifdef TARGET_BGL
265 OPEN(UNIT=scrUnit1,FILE='scratch1',STATUS='UNKNOWN')
266 OPEN(UNIT=scrUnit2,FILE='scratch2',STATUS='UNKNOWN')
267 #else
268 OPEN(UNIT=scrUnit1,STATUS='SCRATCH')
269 OPEN(UNIT=scrUnit2,STATUS='SCRATCH')
270 #endif
271 OPEN(UNIT=modelDataUnit,FILE='data',STATUS='OLD',
272 & IOSTAT=errIO)
273 IF ( errIO .LT. 0 ) THEN
274 WRITE(msgBuf,'(A)')
275 & 'S/R INI_PARMS'
276 CALL PRINT_ERROR( msgBuf , 1)
277 WRITE(msgBuf,'(A)')
278 & 'Unable to open model parameter'
279 CALL PRINT_ERROR( msgBuf , 1)
280 WRITE(msgBuf,'(A)')
281 & 'file "data"'
282 CALL PRINT_ERROR( msgBuf , 1)
283 CALL MODELDATA_EXAMPLE( myThid )
284 STOP 'ABNORMAL END: S/R INI_PARMS'
285 ENDIF
286
287 DO WHILE ( .TRUE. )
288 READ(modelDataUnit,FMT='(A)',END=1001) RECORD
289 IL = MAX(ILNBLNK(RECORD),1)
290 IF ( RECORD(1:1) .NE. commentCharacter ) THEN
291 CALL NML_SET_TERMINATOR( RECORD )
292 WRITE(UNIT=scrUnit1,FMT='(A)') RECORD(:IL)
293 ENDIF
294 WRITE(UNIT=scrUnit2,FMT='(A)') RECORD(:IL)
295 ENDDO
296 1001 CONTINUE
297 CLOSE(modelDataUnit)
298
299 C-- Report contents of model parameter file
300 WRITE(msgBuf,'(A)')
301 &'// ======================================================='
302 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
303 & SQUEEZE_RIGHT , 1)
304 WRITE(msgBuf,'(A)') '// Model parameter file "data"'
305 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
306 & SQUEEZE_RIGHT , 1)
307 WRITE(msgBuf,'(A)')
308 &'// ======================================================='
309 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
310 & SQUEEZE_RIGHT , 1)
311 iUnit = scrUnit2
312 REWIND(iUnit)
313 DO WHILE ( .TRUE. )
314 READ(UNIT=iUnit,FMT='(A)',END=2001) RECORD
315 IL = MAX(ILNBLNK(RECORD),1)
316 WRITE(msgBuf,'(A,A)') '>',RECORD(:IL)
317 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
318 & SQUEEZE_RIGHT , 1)
319 ENDDO
320 2001 CONTINUE
321 CLOSE(iUnit)
322 WRITE(msgBuf,'(A)') ' '
323 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
324 & SQUEEZE_RIGHT , 1)
325
326
327 C-- Read settings from model parameter file "data".
328 iUnit = scrUnit1
329 REWIND(iUnit)
330
331 C-- Set default "physical" parameters
332 viscAhW = UNSET_RL
333 viscA4W = UNSET_RL
334 viscAhD = UNSET_RL
335 viscAhZ = UNSET_RL
336 viscA4D = UNSET_RL
337 viscA4Z = UNSET_RL
338 viscAz = UNSET_RL
339 viscAr = UNSET_RL
340 viscAp = UNSET_RL
341 diffKzT = UNSET_RL
342 diffKpT = UNSET_RL
343 diffKrT = UNSET_RL
344 diffKzS = UNSET_RL
345 diffKpS = UNSET_RL
346 diffKrS = UNSET_RL
347 DO k=1,Nr
348 diffKrNrT(k) = UNSET_RL
349 diffKrNrS(k) = UNSET_RL
350 tRef(k) = UNSET_RL
351 sRef(k) = UNSET_RL
352 ENDDO
353 gBaro = UNSET_RL
354 rhoConst = UNSET_RL
355 omega = UNSET_RL
356 hFacMinDr = UNSET_RL
357 hFacMinDz = UNSET_RL
358 hFacMinDp = UNSET_RL
359 rhoConstFresh = UNSET_RL
360 convertFW2Salt = UNSET_RL
361 tAlpha = UNSET_RL
362 sBeta = UNSET_RL
363 tempVertAdvScheme = 0
364 saltVertAdvScheme = 0
365 C-- z,p,r coord input switching.
366 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; starts to read PARM01'
367 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
368 & SQUEEZE_RIGHT , 1)
369 READ(UNIT=iUnit,NML=PARM01) !,IOSTAT=errIO)
370 IF ( errIO .LT. 0 ) THEN
371 WRITE(msgBuf,'(A)')
372 & 'S/R INI_PARMS'
373 CALL PRINT_ERROR( msgBuf , 1)
374 WRITE(msgBuf,'(A)')
375 & 'Error reading numerical model '
376 CALL PRINT_ERROR( msgBuf , 1)
377 WRITE(msgBuf,'(A)')
378 & 'parameter file "data"'
379 CALL PRINT_ERROR( msgBuf , 1)
380 WRITE(msgBuf,'(A)')
381 & 'Problem in namelist PARM01'
382 CALL PRINT_ERROR( msgBuf , 1)
383 CALL MODELDATA_EXAMPLE( myThid )
384 STOP 'ABNORMAL END: S/R INI_PARMS'
385 ELSE
386 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM01 : OK'
387 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
388 & SQUEEZE_RIGHT , 1)
389 ENDIF
390
391 C- set the type of vertical coordinate and type of fluid
392 C according to buoyancyRelation
393 usingPCoords = .FALSE.
394 usingZCoords = .FALSE.
395 fluidIsAir = .FALSE.
396 fluidIsWater = .FALSE.
397 IF ( buoyancyRelation.EQ.'ATMOSPHERIC' ) THEN
398 usingPCoords = .TRUE.
399 fluidIsAir = .TRUE.
400 ELSEIF ( buoyancyRelation.EQ.'OCEANICP') THEN
401 usingPCoords = .TRUE.
402 fluidIsWater = .TRUE.
403 ELSEIF ( buoyancyRelation.EQ.'OCEANIC' ) THEN
404 usingZCoords = .TRUE.
405 fluidIsWater = .TRUE.
406 ELSE
407 WRITE(msgBuf,'(2A)') 'S/R INI_PARMS:',
408 & ' Bad value of buoyancyRelation '
409 CALL PRINT_ERROR( msgBuf , myThid)
410 STOP 'ABNORMAL END: S/R INI_PARMS'
411 ENDIF
412
413 IF ( .NOT.rigidLid .AND.
414 & .NOT.implicitFreeSurface ) THEN
415 C- No barotropic solver selected => use implicitFreeSurface as default
416 WRITE(msgBuf,'(A)')
417 & 'S/R INI_PARMS: No request for barotropic solver'
418 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
419 & SQUEEZE_RIGHT , myThid)
420 WRITE(msgBuf,'(A)')
421 & 'S/R INI_PARMS: => Use implicitFreeSurface as default'
422 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
423 & SQUEEZE_RIGHT , myThid)
424 implicitFreeSurface = .TRUE.
425 ENDIF
426 IF ( implicitFreeSurface ) freeSurfFac = 1.D0
427 IF ( rigidLid ) freeSurfFac = 0.D0
428 IF ( gBaro .EQ. UNSET_RL ) gBaro=gravity
429 IF ( rhoConst .EQ. UNSET_RL ) rhoConst=rhoNil
430 IF ( rhoConstFresh .EQ. UNSET_RL ) rhoConstFresh=rhoConst
431 IF ( omega .EQ. UNSET_RL ) THEN
432 omega = 0. _d 0
433 IF ( rotationPeriod .NE. 0. _d 0 )
434 & omega = 2.D0 * PI / rotationPeriod
435 ELSEIF ( omega .EQ. 0. _d 0 ) THEN
436 rotationPeriod = 0. _d 0
437 ELSE
438 rotationPeriod = 2.D0 * PI / omega
439 ENDIF
440 IF (atm_Rd .EQ. UNSET_RL) THEN
441 atm_Rd = atm_Cp * atm_kappa
442 ELSE
443 atm_kappa = atm_Rd / atm_Cp
444 ENDIF
445 C- set default vertical profile for temperature: tRef
446 tracerDefault = 20.
447 IF ( fluidIsAir ) tracerDefault = 300.
448 DO k=1,Nr
449 IF (tRef(k).EQ.UNSET_RL) tRef(k) = tracerDefault
450 tracerDefault = tRef(k)
451 ENDDO
452 C- set default vertical profile for salinity/water-vapour: sRef
453 tracerDefault = 30.
454 IF ( fluidIsAir ) tracerDefault = 0.
455 DO k=1,Nr
456 IF (sRef(k).EQ.UNSET_RL) sRef(k) = tracerDefault
457 tracerDefault = sRef(k)
458 ENDDO
459 C-- On/Off flags for each terms of the momentum equation
460 nonHydrostatic = momStepping .AND. nonHydrostatic
461 quasiHydrostatic = momStepping .AND. quasiHydrostatic
462 momAdvection = momStepping .AND. momAdvection
463 momViscosity = momStepping .AND. momViscosity
464 momForcing = momStepping .AND. momForcing
465 useCoriolis = momStepping .AND. useCoriolis
466 useCDscheme = momStepping .AND. useCDscheme
467 momPressureForcing= momStepping .AND. momPressureForcing
468 momImplVertAdv = momAdvection .AND. momImplVertAdv
469 implicitViscosity= momViscosity .AND. implicitViscosity
470 C-- Momentum viscosity on/off flag.
471 IF ( momViscosity ) THEN
472 vfFacMom = 1.D0
473 ELSE
474 vfFacMom = 0.D0
475 ENDIF
476 C-- Momentum advection on/off flag.
477 IF ( momAdvection ) THEN
478 afFacMom = 1.D0
479 ELSE
480 afFacMom = 0.D0
481 ENDIF
482 C-- Momentum forcing on/off flag.
483 IF ( momForcing ) THEN
484 foFacMom = 1.D0
485 ELSE
486 foFacMom = 0.D0
487 ENDIF
488 C-- Coriolis term on/off flag.
489 IF ( useCoriolis ) THEN
490 cfFacMom = 1.D0
491 ELSE
492 cfFacMom = 0.D0
493 ENDIF
494 C-- Pressure term on/off flag.
495 IF ( momPressureForcing ) THEN
496 pfFacMom = 1.D0
497 ELSE
498 pfFacMom = 0.D0
499 ENDIF
500 C-- Metric terms on/off flag.
501 IF ( metricTerms ) THEN
502 mTFacMom = 1.D0
503 ELSE
504 mTFacMom = 0.D0
505 ENDIF
506 C-- Non-hydrostatic/quasi-hydrostatic
507 IF (nonHydrostatic.AND.quasiHydrostatic) THEN
508 WRITE(msgBuf,'(A)')
509 & 'Illegal: both nonHydrostatic = quasiHydrostatic = TRUE'
510 CALL PRINT_ERROR( msgBuf , myThid)
511 STOP 'ABNORMAL END: S/R INI_PARMS'
512 ENDIF
513 C-- Advection and Forcing for Temp and salt on/off flags
514 tempAdvection = tempStepping .AND. tempAdvection
515 tempForcing = tempStepping .AND. tempForcing
516 saltAdvection = saltStepping .AND. saltAdvection
517 saltForcing = saltStepping .AND. saltForcing
518 tempImplVertAdv = tempAdvection .AND. tempImplVertAdv
519 saltImplVertAdv = saltAdvection .AND. saltImplVertAdv
520 IF (tempVertAdvScheme.EQ.0) tempVertAdvScheme = tempAdvScheme
521 IF (saltVertAdvScheme.EQ.0) saltVertAdvScheme = saltAdvScheme
522 C-- horizontal viscosity for vertical moments
523 IF ( viscAhW .EQ. UNSET_RL ) viscAhW = viscAh
524 IF ( viscA4W .EQ. UNSET_RL ) viscA4W = viscA4
525 C-- horizontal viscosity (acting on Divergence or Vorticity)
526 IF ( viscAhD .EQ. UNSET_RL ) viscAhD = viscAh
527 IF ( viscAhZ .EQ. UNSET_RL ) viscAhZ = viscAh
528 IF ( viscA4D .EQ. UNSET_RL ) viscA4D = viscA4
529 IF ( viscA4Z .EQ. UNSET_RL ) viscA4Z = viscA4
530 C-- z,p,r coord input switching.
531 IF ( viscAz .NE. UNSET_RL ) zCoordInputData = .TRUE.
532 IF ( viscAp .NE. UNSET_RL ) pCoordInputData = .TRUE.
533 IF ( viscAr .NE. UNSET_RL ) rCoordInputData = .TRUE.
534 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAz
535 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscAp
536 IF ( viscAr .EQ. UNSET_RL ) viscAr = viscArDefault
537
538 IF ( diffKzT .NE. UNSET_RL ) zCoordInputData = .TRUE.
539 IF ( diffKpT .NE. UNSET_RL ) pCoordInputData = .TRUE.
540 IF ( diffKrT .NE. UNSET_RL ) rCoordInputData = .TRUE.
541 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKzT
542 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKpT
543 IF ( diffKrT .EQ. UNSET_RL ) diffKrT = diffKrTDefault
544 diffKrSet = .TRUE.
545 DO k=1,Nr
546 IF ( diffKrNrT(k).EQ. UNSET_RL ) diffKrSet = .FALSE.
547 ENDDO
548 IF ( .NOT.diffKrSet ) THEN
549 DO k=1,Nr
550 diffKrNrT(k) = diffKrT
551 ENDDO
552 ELSEIF ( diffKrT.NE.diffKrTDefault ) THEN
553 WRITE(msgBuf,'(2A)') 'S/R INI_PARMS: Cannot set both ',
554 & 'diffKrNrT and diffKrT (or Kp,Kz) in input file data'
555 CALL PRINT_ERROR( msgBuf , myThid)
556 STOP 'ABNORMAL END: S/R INI_PARMS'
557 ENDIF
558
559 IF ( diffKzS .NE. UNSET_RL ) zCoordInputData = .TRUE.
560 IF ( diffKpS .NE. UNSET_RL ) pCoordInputData = .TRUE.
561 IF ( diffKrS .NE. UNSET_RL ) rCoordInputData = .TRUE.
562 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKzS
563 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKpS
564 IF ( diffKrS .EQ. UNSET_RL ) diffKrS = diffKrSDefault
565 diffKrSet = .TRUE.
566 DO k=1,Nr
567 IF ( diffKrNrS(k).EQ. UNSET_RL ) diffKrSet = .FALSE.
568 ENDDO
569 IF ( .NOT.diffKrSet ) THEN
570 DO k=1,Nr
571 diffKrNrS(k) = diffKrS
572 ENDDO
573 ELSEIF ( diffKrS.NE.diffKrSDefault ) THEN
574 WRITE(msgBuf,'(2A)') 'S/R INI_PARMS: Cannot set both ',
575 & 'diffKrNrS and diffKrS (or Kp,Kz) in input file data'
576 CALL PRINT_ERROR( msgBuf , myThid)
577 STOP 'ABNORMAL END: S/R INI_PARMS'
578 ENDIF
579
580 IF ( hFacMinDz .NE. UNSET_RL ) zCoordInputData = .TRUE.
581 IF ( hFacMinDp .NE. UNSET_RL ) pCoordInputData = .TRUE.
582 IF ( hFacMinDr .NE. UNSET_RL ) rCoordInputData = .TRUE.
583 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDz
584 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDp
585 IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr = hFacMinDrDefault
586
587 IF (convertFW2Salt.EQ.UNSET_RL) THEN
588 convertFW2Salt = 35.
589 IF (useRealFreshWaterFlux) convertFW2Salt=-1
590 ENDIF
591
592 IF ( ivdc_kappa .NE. 0. .AND. .NOT. implicitDiffusion ) THEN
593 WRITE(msgBuf,'(A,A)')
594 & 'S/R INI_PARMS: To use ivdc_kappa you must enable implicit',
595 & ' vertical diffusion.'
596 CALL PRINT_ERROR( msgBuf , myThid)
597 STOP 'ABNORMAL END: S/R INI_PARMS'
598 ENDIF
599
600 coordsSet = 0
601 IF ( zCoordInputData ) coordsSet = coordsSet + 1
602 IF ( pCoordInputData ) coordsSet = coordsSet + 1
603 IF ( rCoordInputData ) coordsSet = coordsSet + 1
604 IF ( coordsSet .GT. 1 ) THEN
605 WRITE(msgBuf,'(A)')
606 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
607 CALL PRINT_ERROR( msgBuf , myThid)
608 STOP 'ABNORMAL END: S/R INI_PARMS'
609 ENDIF
610 IF ( rhoConst .LE. 0. ) THEN
611 WRITE(msgBuf,'(A)')
612 & 'S/R INI_PARMS: rhoConst must be greater than 0.'
613 CALL PRINT_ERROR( msgBuf , myThid)
614 STOP 'ABNORMAL END: S/R INI_PARMS'
615 ELSE
616 recip_rhoConst = 1.D0 / rhoConst
617 ENDIF
618 IF ( rhoNil .LE. 0. ) THEN
619 WRITE(msgBuf,'(A)')
620 & 'S/R INI_PARMS: rhoNil must be greater than 0.'
621 CALL PRINT_ERROR( msgBuf , myThid)
622 STOP 'ABNORMAL END: S/R INI_PARMS'
623 ELSE
624 recip_rhoNil = 1.D0 / rhoNil
625 ENDIF
626 IF ( HeatCapacity_Cp .LE. 0. ) THEN
627 WRITE(msgBuf,'(A)')
628 & 'S/R INI_PARMS: HeatCapacity_Cp must be greater than 0.'
629 CALL PRINT_ERROR( msgBuf , myThid)
630 STOP 'ABNORMAL END: S/R INI_PARMS'
631 ELSE
632 recip_Cp = 1.D0 / HeatCapacity_Cp
633 ENDIF
634 IF ( gravity .LE. 0. ) THEN
635 WRITE(msgBuf,'(A)')
636 & 'S/R INI_PARMS: gravity must be greater than 0.'
637 CALL PRINT_ERROR( msgBuf , myThid)
638 STOP 'ABNORMAL END: S/R INI_PARMS'
639 ELSE
640 recip_gravity = 1.D0 / gravity
641 ENDIF
642 C This flags are now passed to RW and MDSIO packages in ini_model_io.F
643 C Set globalFiles flag for READ_WRITE_FLD package
644 c CALL SET_WRITE_GLOBAL_FLD( globalFiles )
645 C Set globalFiles flag for READ_WRITE_REC package
646 c CALL SET_WRITE_GLOBAL_REC( globalFiles )
647 C Set globalFiles flag for READ_WRITE_REC package
648 c CALL SET_WRITE_GLOBAL_PICKUP( globalFiles )
649
650 C Check for retired parameters still being used
651 nRetired = 0
652 IF ( zonal_filt_lat .NE. UNSET_RL ) THEN
653 nRetired = nRetired+1
654 WRITE(msgBuf,'(A,A)')
655 & 'S/R INI_PARMS: Paramater "zonal_filt_lat" is',
656 & ' no longer allowed in file "data".'
657 CALL PRINT_ERROR( msgBuf , myThid)
658 WRITE(msgBuf,'(A,A)')
659 & 'S/R INI_PARMS: Paramater "zonal_filt_lat" is',
660 & ' now read from file "data.zonfilt".'
661 CALL PRINT_ERROR( msgBuf , myThid)
662 ENDIF
663 IF ( gravitySign .NE. UNSET_RL ) THEN
664 nRetired = nRetired+1
665 WRITE(msgBuf,'(A,A)')
666 & 'S/R INI_PARMS: "gravitySign" is set according to vertical ',
667 & ' coordinate and is no longer allowed in file "data".'
668 CALL PRINT_ERROR( msgBuf , myThid)
669 ENDIF
670 IF ( tracerAdvScheme .NE. UNSET_I ) THEN
671 nRetired = nRetired+1
672 WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: "tracerAdvScheme" ',
673 & '(old passive tracer code) is no longer allowed in file "data"'
674 CALL PRINT_ERROR( msgBuf , myThid)
675 ENDIF
676 IF ( trac_EvPrRn .NE. UNSET_RL ) THEN
677 nRetired = nRetired+1
678 WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: "trac_EvPrRn" ',
679 & '(old passive tracer code) is no longer allowed in file "data"'
680 CALL PRINT_ERROR( msgBuf , myThid)
681 ENDIF
682 IF ( .NOT. tempDiffusion ) THEN
683 nRetired = nRetired+1
684 WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: "tempDiffusion" ',
685 & 'is no longer allowed in file "data"'
686 CALL PRINT_ERROR( msgBuf , myThid)
687 WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: to turn off diffusion',
688 & ' => set diffusivity to zero'
689 CALL PRINT_ERROR( msgBuf , myThid)
690 ENDIF
691 IF ( .NOT. saltDiffusion ) THEN
692 nRetired = nRetired+1
693 WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: "saltDiffusion" ',
694 & 'is no longer allowed in file "data"'
695 CALL PRINT_ERROR( msgBuf , myThid)
696 WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: to turn off diffusion',
697 & ' => set diffusivity to zero'
698 CALL PRINT_ERROR( msgBuf , myThid)
699 ENDIF
700 IF ( viscAstrain .NE. UNSET_RL ) THEN
701 nRetired = nRetired+1
702 WRITE(msgBuf,'(A,A)')
703 & 'S/R INI_PARMS: "viscAstrain" ',
704 & 'is no longer allowed in file "data"'
705 CALL PRINT_ERROR( msgBuf , myThid)
706 WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: to use Strain & Tension',
707 & ' formulation => set useStrainTensionVisc to TRUE'
708 CALL PRINT_ERROR( msgBuf , myThid)
709 ENDIF
710 IF ( viscAtension .NE. UNSET_RL ) THEN
711 nRetired = nRetired+1
712 WRITE(msgBuf,'(A,A)')
713 & 'S/R INI_PARMS: "viscAtension" ',
714 & 'is no longer allowed in file "data"'
715 CALL PRINT_ERROR( msgBuf , myThid)
716 WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: to use Strain & Tension',
717 & ' formulation => set useStrainTensionVisc to TRUE'
718 CALL PRINT_ERROR( msgBuf , myThid)
719 ENDIF
720 IF (useAnisotropicViscAgridMax .NEQV. .TRUE. ) THEN
721 nRetired = nRetired+1
722 WRITE(msgBuf,'(A,A)')
723 & 'S/R INI_PARMS: "useAnisotropicViscAgridMax" ',
724 & 'is not allowed in "data" substitute useAreaViscLength=true'
725 CALL PRINT_ERROR( msgBuf , myThid)
726 ENDIF
727
728 C-- Elliptic solver parameters
729 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; starts to read PARM02'
730 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
731 & SQUEEZE_RIGHT , 1)
732 READ(UNIT=iUnit,NML=PARM02) !,IOSTAT=errIO)
733 IF ( errIO .LT. 0 ) THEN
734 WRITE(msgBuf,'(A)')
735 & 'S/R INI_PARMS'
736 CALL PRINT_ERROR( msgBuf , 1)
737 WRITE(msgBuf,'(A)')
738 & 'Error reading numerical model '
739 CALL PRINT_ERROR( msgBuf , 1)
740 WRITE(msgBuf,'(A)')
741 & 'parameter file "data".'
742 CALL PRINT_ERROR( msgBuf , 1)
743 WRITE(msgBuf,'(A)')
744 & 'Problem in namelist PARM02'
745 CALL PRINT_ERROR( msgBuf , 1)
746 CALL MODELDATA_EXAMPLE( myThid )
747 STOP 'ABNORMAL END: S/R INI_PARMS'
748 ELSE
749 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM02 : OK'
750 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
751 & SQUEEZE_RIGHT , 1)
752 ENDIF
753
754 C-- Time stepping parameters
755 rCD = -1.D0
756 latBandClimRelax = UNSET_RL
757 deltaTtracer = 0. _d 0
758 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; starts to read PARM03'
759 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
760 & SQUEEZE_RIGHT , 1)
761 READ(UNIT=iUnit,NML=PARM03) !,IOSTAT=errIO)
762 IF ( errIO .LT. 0 ) THEN
763 WRITE(msgBuf,'(A)')
764 & 'S/R INI_PARMS'
765 CALL PRINT_ERROR( msgBuf , 1)
766 WRITE(msgBuf,'(A)')
767 & 'Error reading numerical model '
768 CALL PRINT_ERROR( msgBuf , 1)
769 WRITE(msgBuf,'(A)')
770 & 'parameter file "data"'
771 CALL PRINT_ERROR( msgBuf , 1)
772 WRITE(msgBuf,'(A)')
773 & 'Problem in namelist PARM03'
774 CALL PRINT_ERROR( msgBuf , 1)
775 CALL MODELDATA_EXAMPLE( myThid )
776 STOP 'ABNORMAL END: S/R INI_PARMS'
777 ELSE
778 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM03 : OK'
779 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
780 & SQUEEZE_RIGHT , 1)
781 ENDIF
782 C Process "timestepping" params
783 C o Time step size
784 IF ( deltaTtracer .NE. dTtracerLev(1) .AND.
785 & deltaTtracer .NE. 0. .AND. dTtracerLev(1) .NE. 0. ) THEN
786 WRITE(msgBuf,'(A)')
787 & 'S/R INI_PARMS: deltaTtracer & dTtracerLev(1) not equal'
788 CALL PRINT_ERROR( msgBuf , myThid)
789 STOP 'ABNORMAL END: S/R INI_PARMS'
790 ELSEIF ( dTtracerLev(1) .NE. 0. ) THEN
791 deltaTtracer = dTtracerLev(1)
792 ENDIF
793 IF ( deltaT .EQ. 0. ) deltaT = deltaTmom
794 IF ( deltaT .EQ. 0. ) deltaT = deltaTtracer
795 IF ( deltaTmom .EQ. 0. ) deltaTmom = deltaT
796 IF ( deltaTtracer .EQ. 0. ) deltaTtracer = deltaT
797 IF ( deltaTClock .EQ. 0. ) deltaTClock = deltaT
798 DO k=1,Nr
799 IF (dTtracerLev(k).EQ.0.) dTtracerLev(k)= deltaTtracer
800 ENDDO
801 C Note that this line should set deltaFreesurf=deltaTmom
802 C but this would change a lot of existing set-ups so we are
803 C obliged to set the default inappropriately.
804 C Be advised that when using asynchronous time stepping
805 C it is better to set deltaTreesurf=deltaTtracer
806 IF ( deltaTfreesurf .EQ. 0. ) deltaTfreesurf = deltaTmom
807 IF ( periodicExternalForcing ) THEN
808 IF ( externForcingCycle*externForcingPeriod .EQ. 0. ) THEN
809 WRITE(msgBuf,'(A)')
810 & 'S/R INI_PARMS: externForcingCycle,externForcingPeriod =0'
811 CALL PRINT_ERROR( msgBuf , 1)
812 STOP 'ABNORMAL END: S/R INI_PARMS'
813 ENDIF
814 IF ( INT(externForcingCycle/externForcingPeriod) .NE.
815 & externForcingCycle/externForcingPeriod ) THEN
816 WRITE(msgBuf,'(A)')
817 & 'S/R INI_PARMS: externForcingCycle <> N*externForcingPeriod'
818 CALL PRINT_ERROR( msgBuf , 1)
819 STOP 'ABNORMAL END: S/R INI_PARMS'
820 ENDIF
821 IF ( externForcingCycle.lt.externForcingPeriod ) THEN
822 WRITE(msgBuf,'(A)')
823 & 'S/R INI_PARMS: externForcingCycle < externForcingPeriod'
824 CALL PRINT_ERROR( msgBuf , 1)
825 STOP 'ABNORMAL END: S/R INI_PARMS'
826 ENDIF
827 IF ( externForcingPeriod.lt.deltaTclock ) THEN
828 WRITE(msgBuf,'(A)')
829 & 'S/R INI_PARMS: externForcingPeriod < deltaTclock'
830 CALL PRINT_ERROR( msgBuf , 1)
831 STOP 'ABNORMAL END: S/R INI_PARMS'
832 ENDIF
833 ENDIF
834 C o Convection frequency
835 IF ( cAdjFreq .LT. 0. ) THEN
836 cAdjFreq = deltaTClock
837 ENDIF
838 IF ( ivdc_kappa .NE. 0. .AND. cAdjFreq .NE. 0. ) THEN
839 WRITE(msgBuf,'(A,A)')
840 & 'S/R INI_PARMS: You have enabled both ivdc_kappa and',
841 & ' convective adjustment.'
842 CALL PRINT_ERROR( msgBuf , myThid)
843 STOP 'ABNORMAL END: S/R INI_PARMS'
844 ENDIF
845 IF (useCDscheme) THEN
846 C o CD coupling (CD scheme):
847 IF ( tauCD .EQ. 0.D0 ) tauCD = deltaTmom
848 IF ( rCD .LT. 0. ) rCD = 1. _d 0 - deltaTMom/tauCD
849 ENDIF
850 C o Temperature climatology relaxation time scale
851 IF ( tauThetaClimRelax .EQ. 0.D0 ) THEN
852 doThetaClimRelax = .FALSE.
853 ELSE
854 doThetaClimRelax = .TRUE.
855 ENDIF
856 C o Salinity climatology relaxation time scale
857 IF ( tauSaltClimRelax .EQ. 0.D0 ) THEN
858 doSaltClimRelax = .FALSE.
859 ELSE
860 doSaltClimRelax = .TRUE.
861 ENDIF
862 C o Tracer 1 climatology relaxation time scale
863 IF ( tauTr1ClimRelax .EQ. 0.D0 ) THEN
864 doTr1ClimRelax = .FALSE.
865 lambdaTr1ClimRelax = 0.D0
866 ELSE
867 doTr1ClimRelax = .TRUE.
868 lambdaTr1ClimRelax = 1./tauTr1ClimRelax
869 ENDIF
870
871 C o Base time
872 IF ( nIter0.NE.0 .AND. startTime.NE.0. .AND. baseTime.EQ.0. )
873 & baseTime = startTime - deltaTClock*float(nIter0)
874 C o Start time
875 IF ( nIter0 .NE. 0 .AND. startTime .EQ. 0. )
876 & startTime = baseTime + deltaTClock*float(nIter0)
877 C o nIter0
878 IF ( nIter0 .EQ. 0 .AND. startTime .NE. baseTime )
879 & nIter0 = NINT( (startTime-baseTime)/deltaTClock )
880
881 C o nTimeSteps 1
882 IF ( nTimeSteps .EQ. 0 .AND. nEndIter .NE. 0 )
883 & nTimeSteps = nEndIter-nIter0
884 C o nTimeSteps 2
885 IF ( nTimeSteps .EQ. 0 .AND. endTime .NE. 0. )
886 & nTimeSteps = NINT((endTime-startTime)/deltaTclock)
887 C o nEndIter 1
888 IF ( nEndIter .EQ. 0 .AND. nTimeSteps .NE. 0 )
889 & nEndIter = nIter0+nTimeSteps
890 C o nEndIter 2
891 IF ( nEndIter .EQ. 0 .AND. endTime .NE. 0. )
892 & nEndIter = NINT((endTime-baseTime)/deltaTclock)
893 C o End Time 1
894 IF ( endTime .EQ. 0. .AND. nTimeSteps .NE. 0 )
895 & endTime = startTime + deltaTClock*float(nTimeSteps)
896 C o End Time 2
897 IF ( endTime .EQ. 0. .AND. nEndIter .NE. 0 )
898 & endTime = baseTime + deltaTClock*float(nEndIter)
899
900 C o Consistent?
901 IF ( startTime .NE. baseTime+deltaTClock*float(nIter0) ) THEN
902 WRITE(msgBuf,'(A)')
903 & 'S/R INI_PARMS: startTime, baseTime and nIter0 are inconsistent'
904 CALL PRINT_ERROR( msgBuf , 1)
905 WRITE(msgBuf,'(A)')
906 & 'S/R INI_PARMS: Perhaps more than two were set at once'
907 CALL PRINT_ERROR( msgBuf , 1)
908 STOP 'ABNORMAL END: S/R INI_PARMS'
909 ENDIF
910 IF ( nEndIter .NE. nIter0+nTimeSteps ) THEN
911 WRITE(msgBuf,'(A)')
912 & 'S/R INI_PARMS: nIter0, nTimeSteps and nEndIter are inconsistent'
913 CALL PRINT_ERROR( msgBuf , 1)
914 WRITE(msgBuf,'(A)')
915 & 'S/R INI_PARMS: Perhaps more than two were set at once'
916 CALL PRINT_ERROR( msgBuf , 1)
917 STOP 'ABNORMAL END: S/R INI_PARMS'
918 ENDIF
919 IF ( nTimeSteps .NE. NINT((endTime-startTime)/deltaTClock) )
920 & THEN
921 WRITE(msgBuf,'(A)')
922 & 'S/R INI_PARMS: both endTime and nTimeSteps have been set'
923 CALL PRINT_ERROR( msgBuf , 1)
924 WRITE(msgBuf,'(A)')
925 & 'S/R INI_PARMS: but are inconsistent'
926 CALL PRINT_ERROR( msgBuf , 1)
927 STOP 'ABNORMAL END: S/R INI_PARMS'
928 ENDIF
929
930 C o Monitor (should also add CPP flag for monitor?)
931 IF (monitorFreq.LT.0.) THEN
932 monitorFreq=0.
933 IF (dumpFreq.NE.0.) monitorFreq=dumpFreq
934 IF (diagFreq.NE.0..AND.diagFreq.LT.monitorFreq)
935 & monitorFreq=diagFreq
936 IF (taveFreq.NE.0..AND.taveFreq.LT.monitorFreq)
937 & monitorFreq=taveFreq
938 IF (chkPtFreq.NE.0..AND.chkPtFreq.LT.monitorFreq)
939 & monitorFreq=chkPtFreq
940 IF (pChkPtFreq.NE.0..AND.pChkPtFreq.LT.monitorFreq)
941 & monitorFreq=pChkPtFreq
942 IF (monitorFreq.EQ.0.) monitorFreq=deltaTclock
943 ENDIF
944
945 C-- Grid parameters
946 C In cartesian coords distances are in metres
947 DO K =1,Nr
948 delZ(K) = UNSET_RL
949 delP(K) = UNSET_RL
950 delR(K) = UNSET_RL
951 ENDDO
952 C In spherical polar distances are in degrees
953 recip_rSphere = 1.D0/rSphere
954 dxSpacing = UNSET_RL
955 dySpacing = UNSET_RL
956 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; starts to read PARM04'
957 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
958 & SQUEEZE_RIGHT , 1)
959 READ(UNIT=iUnit,NML=PARM04,IOSTAT=errIO)
960 IF ( errIO .LT. 0 ) THEN
961 WRITE(msgBuf,'(A)')
962 & 'S/R INI_PARMS'
963 CALL PRINT_ERROR( msgBuf , 1)
964 WRITE(msgBuf,'(A)')
965 & 'Error reading numerical model '
966 CALL PRINT_ERROR( msgBuf , 1)
967 WRITE(msgBuf,'(A)')
968 & 'parameter file "data"'
969 CALL PRINT_ERROR( msgBuf , 1)
970 WRITE(msgBuf,'(A)')
971 & 'Problem in namelist PARM04'
972 CALL PRINT_ERROR( msgBuf , 1)
973 CALL MODELDATA_EXAMPLE( myThid )
974 STOP 'ABNORMAL END: S/R INI_PARMS'
975 ELSE
976 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM04 : OK'
977 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
978 & SQUEEZE_RIGHT , 1)
979 ENDIF
980
981 C Check for retired parameters still being used
982 IF ( rkFac .NE. UNSET_RL ) THEN
983 nRetired = nRetired+1
984 WRITE(msgBuf,'(A,A)')
985 & 'S/R INI_PARMS: "rkFac" has been replaced by -rkSign ',
986 & ' and is no longer allowed in file "data".'
987 CALL PRINT_ERROR( msgBuf , myThid)
988 ENDIF
989 IF ( groundAtK1 ) THEN
990 c nRetired = nRetired+1
991 WRITE(msgBuf,'(A,A)')
992 & 'S/R INI_PARMS: "groundAtK1" is set according to vertical ',
993 & ' coordinate and is no longer allowed in file "data".'
994 CALL PRINT_ERROR( msgBuf , myThid)
995 ENDIF
996
997 C X coordinate : Check for multiple definitions
998 goptCount = 0
999 IF ( delX(1) .NE. UNSET_RL ) goptCount = goptCount + 1
1000 IF ( dxSpacing .NE. UNSET_RL ) goptCount = goptCount + 1
1001 IF ( delXFile .NE. ' ' ) goptCount = goptCount + 1
1002 IF ( goptCount.GT.1 ) THEN
1003 WRITE(msgBuf,'(A,A)') 'Too many specifications for delX:',
1004 & 'Specify only one of delX, dxSpacing or delXfile'
1005 CALL PRINT_ERROR( msgBuf , 1)
1006 STOP 'ABNORMAL END: S/R INI_PARMS'
1007 ENDIF
1008 IF ( dxSpacing .NE. UNSET_RL ) THEN
1009 DO i=1,Nx
1010 delX(i) = dxSpacing
1011 ENDDO
1012 ENDIF
1013 C Y coordinate : Check for multiple definitions
1014 goptCount = 0
1015 IF ( delY(1) .NE. UNSET_RL ) goptCount = goptCount + 1
1016 IF ( dySpacing .NE. UNSET_RL ) goptCount = goptCount + 1
1017 IF ( delYFile .NE. ' ' ) goptCount = goptCount + 1
1018 IF ( goptCount.GT.1 ) THEN
1019 WRITE(msgBuf,'(A,A)') 'Too many specifications for delY:',
1020 & 'Specify only one of delY, dySpacing or delYfile'
1021 CALL PRINT_ERROR( msgBuf , 1)
1022 STOP 'ABNORMAL END: S/R INI_PARMS'
1023 ENDIF
1024 IF ( dySpacing .NE. UNSET_RL ) THEN
1025 DO j=1,Ny
1026 delY(j) = dySpacing
1027 ENDDO
1028 ENDIF
1029 C
1030 IF ( rSphere .NE. 0 ) THEN
1031 recip_rSphere = 1.D0/rSphere
1032 ELSE
1033 recip_rSphere = 0.
1034 ENDIF
1035 C-- Check for conflicting grid definitions.
1036 goptCount = 0
1037 IF ( usingCartesianGrid ) goptCount = goptCount+1
1038 IF ( usingSphericalPolarGrid ) goptCount = goptCount+1
1039 IF ( usingCurvilinearGrid ) goptCount = goptCount+1
1040 IF ( usingCylindricalGrid ) goptCount = goptCount+1
1041 IF ( goptCount .GT. 1 ) THEN
1042 WRITE(msgBuf,'(A)')
1043 & 'S/R INI_PARMS: More than one coordinate system requested'
1044 CALL PRINT_ERROR( msgBuf , myThid)
1045 STOP 'ABNORMAL END: S/R INI_PARMS'
1046 ENDIF
1047 IF ( goptCount .LT. 1 ) THEN
1048 C- No horizontal grid is specified => use Cartesian grid as default:
1049 WRITE(msgBuf,'(A)')
1050 & 'S/R INI_PARMS: No horizontal grid requested'
1051 CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
1052 & SQUEEZE_RIGHT , myThid)
1053 WRITE(msgBuf,'(A)')
1054 & 'S/R INI_PARMS: => Use Cartesian Grid as default'
1055 CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
1056 & SQUEEZE_RIGHT , myThid)
1057 usingCartesianGrid = .TRUE.
1058 ENDIF
1059 C-- Make metric term settings consistent with underlying grid.
1060 IF ( usingCartesianGrid ) THEN
1061 usingSphericalPolarMterms = .FALSE.
1062 metricTerms = .FALSE.
1063 useNHMTerms = .FALSE.
1064 mTFacMom = 0.
1065 useBetaPlaneF = .TRUE.
1066 ENDIF
1067 C-- Make metric term settings consistent with underlying grid.
1068 IF ( usingCylindricalGrid) THEN
1069 usingSphericalPolarMterms = .FALSE.
1070 metricTerms = .FALSE.
1071 useNHMTerms = .FALSE.
1072 mTFacMom = 1.
1073 useBetaPlaneF = .TRUE.
1074 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; Cylinder OK'
1075 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1076 & SQUEEZE_RIGHT , 1)
1077
1078 ENDIF
1079
1080 IF ( usingSphericalPolarGrid ) THEN
1081 useSphereF = .TRUE.
1082 usingSphericalPolarMterms = metricTerms
1083 ENDIF
1084 IF ( usingCurvilinearGrid ) THEN
1085 useSphereF = .TRUE.
1086 metricTerms = .FALSE.
1087 useNHMTerms = .FALSE.
1088 ENDIF
1089 C-- Set default for latitude-band where relaxation to climatology applies
1090 C note: done later (once domain size is known) if using CartesianGrid
1091 IF ( latBandClimRelax .EQ. UNSET_RL) THEN
1092 IF ( usingSphericalPolarGrid ) latBandClimRelax= 180. _d 0
1093 IF ( usingCurvilinearGrid ) latBandClimRelax= 180. _d 0
1094 ENDIF
1095 C-- set cell Center depth and put Interface at the middle between 2 C
1096 setCenterDr = .FALSE.
1097 IF (delRc(1).NE.UNSET_RL) setCenterDr=.TRUE.
1098 DO K=1,Nr+1
1099 IF (delRc(K).EQ.UNSET_RL) setCenterDr = .FALSE.
1100 ENDDO
1101 C-- p, z, r coord parameters
1102 DO K = 1, Nr
1103 IF ( delZ(K) .NE. UNSET_RL ) zCoordInputData = .TRUE.
1104 IF ( delP(K) .NE. UNSET_RL ) pCoordInputData = .TRUE.
1105 IF ( delR(K) .NE. UNSET_RL ) rCoordInputData = .TRUE.
1106 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delZ(K)
1107 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delP(K)
1108 IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delRDefault(K)
1109 IF (.NOT.setCenterDr .AND. delR(K).EQ.delRDefault(K) ) THEN
1110 WRITE(msgBuf,'(A,I4)')
1111 & 'S/R INI_PARMS: No value for delZ/delP/delR at K = ',K
1112 CALL PRINT_ERROR( msgBuf , 1)
1113 STOP 'ABNORMAL END: S/R INI_PARMS'
1114 ELSEIF ( setCenterDr .AND. delR(K).NE.delRDefault(K) ) THEN
1115 WRITE(msgBuf,'(2A,I4)') 'S/R INI_PARMS:',
1116 & ' Cannot specify both delRc and delZ/delP/delR at K=', K
1117 CALL PRINT_ERROR( msgBuf , 1)
1118 STOP 'ABNORMAL END: S/R INI_PARMS'
1119 ENDIF
1120 ENDDO
1121 C Check for multiple coordinate systems
1122 CoordsSet = 0
1123 IF ( zCoordInputData ) coordsSet = coordsSet + 1
1124 IF ( pCoordInputData ) coordsSet = coordsSet + 1
1125 IF ( rCoordInputData ) coordsSet = coordsSet + 1
1126 IF ( coordsSet .GT. 1 ) THEN
1127 WRITE(msgBuf,'(A)')
1128 & 'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
1129 CALL PRINT_ERROR( msgBuf , myThid)
1130 STOP 'ABNORMAL END: S/R INI_PARMS'
1131 ENDIF
1132
1133 C-- When using the dynamical pressure in EOS (with Z-coord.),
1134 C needs to activate specific part of the code (restart & exchange)
1135 c useDynP_inEos_Zc = .FALSE.
1136 useDynP_inEos_Zc = ( fluidIsWater .AND. usingZCoords
1137 & .AND. ( eosType .EQ. 'JMD95P' .OR.
1138 & eosType .EQ. 'UNESCO' .OR.
1139 & eosType .EQ. 'MDJWF' ) )
1140
1141 C-- Input files
1142 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; starts to read PARM05'
1143 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1144 & SQUEEZE_RIGHT , 1)
1145 READ(UNIT=iUnit,NML=PARM05) !,IOSTAT=errIO)
1146 IF ( errIO .LT. 0 ) THEN
1147 WRITE(msgBuf,'(A)')
1148 & 'Error reading numerical model '
1149 CALL PRINT_ERROR( msgBuf , 1)
1150 WRITE(msgBuf,'(A)')
1151 & 'parameter file "data"'
1152 CALL PRINT_ERROR( msgBuf , 1)
1153 WRITE(msgBuf,'(A)')
1154 & 'Problem in namelist PARM05'
1155 CALL PRINT_ERROR( msgBuf , 1)
1156 CALL MODELDATA_EXAMPLE( myThid )
1157 STOP 'ABNORMAL END: S/R INI_PARMS'
1158 ELSE
1159 WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM05 : OK'
1160 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1161 & SQUEEZE_RIGHT , 1)
1162 ENDIF
1163
1164 C-- Set factors required for mixing pressure and meters as vertical coordinate.
1165 C rkSign is a "sign" parameter which is used where the orientation of the vertical
1166 C coordinate (pressure or meters) relative to the vertical index (K) is important.
1167 C rkSign = -1 applies when K and the coordinate are in the opposite sense.
1168 C rkSign = 1 applies when K and the coordinate are in the same sense.
1169 C horiVertRatio is a parameter that maps horizontal units to vertical units.
1170 C It is used in certain special cases where lateral and vertical terms are
1171 C being combined and a single frame of reference is needed.
1172 rkSign = -1. _d 0
1173 horiVertRatio = 1. _d 0
1174 gravitySign = -1. _d 0
1175 IF ( usingPCoords ) THEN
1176 gravitySign = 1. _d 0
1177 horiVertRatio = Gravity * rhoConst
1178 ENDIF
1179 convertEmP2rUnit = rhoConstFresh*recip_rhoConst*horiVertRatio
1180 recip_horiVertRatio = 1./horiVertRatio
1181
1182 c-- gradually replacing debugMode by debugLevel
1183 IF ( debugMode ) debugLevel = debLevB
1184
1185 c-- flag for approximate adjoint
1186 IF ( inAdExact ) THEN
1187 inAdTrue = .FALSE.
1188 inAdFALSE = .FALSE.
1189 ELSE
1190 inAdTrue = .TRUE.
1191 inAdFALSE = .FALSE.
1192 ENDIF
1193 C
1194 CLOSE(iUnit)
1195
1196 C-- Check whether any retired parameters were found.
1197 C-- Stop if they were
1198 IF ( nRetired .GT. 0 ) THEN
1199 WRITE(msgBuf,'(A)')
1200 & 'Error reading parameter file "data"'
1201 CALL PRINT_ERROR( msgBuf , 1)
1202 WRITE(msgBuf,'(A)')
1203 & 'some out of date parameters were found in the namelist'
1204 CALL PRINT_ERROR( msgBuf , 1)
1205 STOP 'ABNORMAL END: S/R INI_PARMS'
1206 ENDIF
1207
1208 _END_MASTER(myThid)
1209
1210 C-- Everyone else must wait for the parameters to be loaded
1211 _BARRIER
1212 C
1213 RETURN
1214 END

  ViewVC Help
Powered by ViewVC 1.1.22