/[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.146 - (show annotations) (download)
Thu Feb 10 05:25:37 2005 UTC (19 years, 3 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint57d_post
Changes since 1.145: +13 -3 lines
Add flag(s) inAdExact to help D. grow seaice in the Sahara.

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

  ViewVC Help
Powered by ViewVC 1.1.22