/[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.168 - (show annotations) (download)
Tue Sep 20 21:01:30 2005 UTC (18 years, 8 months ago) by baylor
Branch: MAIN
Changes since 1.167: +2 -1 lines
Adding clean version of mom_calc_visc and viscAhRemax and viscA4Remax params.

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

  ViewVC Help
Powered by ViewVC 1.1.22