/[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.180 - (show annotations) (download)
Tue Mar 7 15:28:02 2006 UTC (18 years, 2 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint58b_post
Changes since 1.179: +16 -2 lines
change forcing_In_AB to affects both T,S forcing and Momentum forcing
(allow to differentiate between forcing components using new integer flags:
 momForcingOutAB=1/0 & tracForcingOutAB=1/0)
and add new flag to put Dissipation tendency out of Adams-Bashforth.

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

  ViewVC Help
Powered by ViewVC 1.1.22