/[MITgcm]/MITgcm/model/src/ini_parms.F
ViewVC logotype

Diff of /MITgcm/model/src/ini_parms.F

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

revision 1.27 by cnh, Mon Aug 24 02:25:01 1998 UTC revision 1.73 by adcroft, Fri Nov 30 21:19:22 2001 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2    C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5    
6    CBOP
7    C     !ROUTINE: INI_PARMS
8    C     !INTERFACE:
9        SUBROUTINE INI_PARMS( myThid )        SUBROUTINE INI_PARMS( myThid )
10  C     /==========================================================\  C     !DESCRIPTION: \bv
11  C     | SUBROUTINE INI_PARMS                                     |  C     *==========================================================*
12  C     | o Routine to set model "parameters"                      |  C     | SUBROUTINE INI_PARMS                                      
13  C     |==========================================================|  C     | o Routine to set model "parameters"                      
14  C     | Notes:                                                   |  C     *==========================================================*
15  C     | ======                                                   |  C     | Notes:                                                    
16  C     | The present version of this routine is a place-holder.   |  C     | ======                                                    
17  C     | A production version needs to handle parameters from an  |  C     | The present version of this routine is a place-holder.    
18  C     | external file and possibly reading in some initial field |  C     | A production version needs to handle parameters from an  
19  C     | values.                                                  |  C     | external file and possibly reading in some initial field  
20  C     \==========================================================/  C     | values.                                                  
21    C     *==========================================================*
22    C     \ev
23    
24    C     !USES:
25          IMPLICIT NONE
26  C     === Global variables ===  C     === Global variables ===
27  #include "SIZE.h"  #include "SIZE.h"
28  #include "EEPARAMS.h"  #include "EEPARAMS.h"
29  #include "PARAMS.h"  #include "PARAMS.h"
30  #include "CG2D.h"  #include "GRID.h"
31    
32    C     !INPUT/OUTPUT PARAMETERS:
33  C     === Routine arguments ===  C     === Routine arguments ===
34  C     myThid - Number of this instance of INI_PARMS  C     myThid - Number of this instance of INI_PARMS
35        INTEGER myThid        INTEGER myThid
36    
37    C     !LOCAL VARIABLES:
38  C     === Local variables ===  C     === Local variables ===
39  C     dxSpacing, dySpacing - Default spacing in X and Y.  C     dxSpacing, dySpacing - Default spacing in X and Y.
40  C                            Units are that of coordinate system  C                            Units are that of coordinate system
# Line 37  C     errIO     - IO error flag Line 47  C     errIO     - IO error flag
47  C     iUnit - Work variable for IO unit number  C     iUnit - Work variable for IO unit number
48  C     record - Work variable for IO buffer  C     record - Work variable for IO buffer
49  C     K, I, J - Loop counters  C     K, I, J - Loop counters
50        REAL dxSpacing  C     xxxDefault - Default value for variable xxx
51        REAL dySpacing        _RL  dxSpacing
52          _RL  dySpacing
53          CHARACTER*(MAX_LEN_FNAM) delXfile
54          CHARACTER*(MAX_LEN_FNAM) delYfile
55        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
56        CHARACTER*(MAX_LEN_PREC) record        CHARACTER*(MAX_LEN_PREC) record
57        INTEGER goptCount        INTEGER goptCount
# Line 48  C     K, I, J - Loop counters Line 61  C     K, I, J - Loop counters
61        EXTERNAL IFNBLNK        EXTERNAL IFNBLNK
62        INTEGER  ILNBLNK        INTEGER  ILNBLNK
63        EXTERNAL ILNBLNK        EXTERNAL ILNBLNK
64    C     Default values for variables which have vertical coordinate system
65    C     dependency.
66          _RL viscArDefault
67          _RL diffKrTDefault
68          _RL diffKrSDefault
69          _RL hFacMinDrDefault
70          _RL delRDefault(Nr)
71          _RS rkFacDefault
72    C     zCoordInputData - These are used to select between different coordinate systems.
73    C     pCoordInputData   The vertical coordinate system in the rest of the model is
74    C     rCoordInputData   written in terms of r. In the model "data" file input data can
75    C     coordsSet         be interms of z, p or r.
76    C                       e.g. delZ or delP or delR for the vertical grid spacing.
77    C                       The following rules apply:
78    C                       All parameters must use the same vertical coordinate system.
79    C                       e.g. delZ and viscAz is legal but
80    C                            delZ and viscAr is an error.
81    C                       Similarly specifyinh delZ and delP is an error.
82    C                       zCoord..., pCoord..., rCoord... are used to flag when z, p or r are
83    C                       used. coordsSet counts how many vertical coordinate systems have been
84    C                       used to specify variables. coordsSet > 1 is an error.
85    C
86          LOGICAL zCoordInputData
87          LOGICAL pCoordInputData
88          LOGICAL rCoordInputData
89          INTEGER coordsSet
90    CEOP
91    
92  C--   Continuous equation parameters  C--   Continuous equation parameters
93        NAMELIST /PARM01/        NAMELIST /PARM01/
94       & gravity, gBaro, rhonil, tAlpha, sBeta, f0, beta,       & gravity, gBaro, rhonil, tAlpha, sBeta, f0, beta, omega,
95       &   viscAh,  viscAz,  viscA4,       & viscAh,  viscAz,  viscA4, cosPower, viscAstrain, viscAtension,
96       &  diffKhT, diffKzT, diffK4T,       & diffKhT, diffKzT, diffK4T,
97       &  diffKhS, diffKzS, diffK4S,       & diffKhS, diffKzS, diffK4S,
98       &  GMmaxslope,GMlength,GMalpha,GMdepth,GMkbackground,GMmaxval,       & tRef, sRef, eosType, Integr_GeoPot,
99       &  tRef, sRef, eosType,       & no_slip_sides,no_slip_bottom,
100       & momViscosity,  momAdvection, momForcing, useCoriolis,       & momViscosity,  momAdvection, momForcing, useCoriolis,
101       & momPressureForcing, metricTerms,       & momPressureForcing, metricTerms, vectorInvariantMomentum,
102       & tempDiffusion, tempAdvection, tempForcing,       & tempDiffusion, tempAdvection, tempForcing,
103       & saltDiffusion, saltAdvection, saltForcing,       & saltDiffusion, saltAdvection, saltForcing,
104         & implicSurfPress, implicDiv2DFlow,
105       & implicitFreeSurface, rigidLid, freeSurfFac, hFacMin, hFacMinDz,       & implicitFreeSurface, rigidLid, freeSurfFac, hFacMin, hFacMinDz,
106       & tempStepping, saltStepping, momStepping, implicitDiffusion,       & exactConserv,uniformLin_PhiSurf,nonlinFreeSurf,hFacInf,hFacSup,
107         & staggerTimeStep,
108         & tempStepping, saltStepping, momStepping, tr1Stepping,
109         & implicitDiffusion, implicitViscosity,
110       & viscAr, diffKrT, diffKrS, hFacMinDr,       & viscAr, diffKrT, diffKrS, hFacMinDr,
111       & rhoConst       & viscAp, diffKpT, diffKpS, hFacMinDp,
112         & rhoConst, buoyancyRelation, HeatCapacity_Cp,
113         & writeBinaryPrec, readBinaryPrec, writeStatePrec,
114         & nonHydrostatic, globalFiles,
115         & allowFreezing, ivdc_kappa,
116         & zonal_filt_lat, zonal_filt_sinpow, zonal_filt_cospow,
117         & bottomDragLinear,bottomDragQuadratic,
118         & usePickupBeforeC35, debugMode,
119         & readPickupWithTracer, writePickupWithTracer,
120         & tempAdvScheme, saltAdvScheme, tracerAdvScheme,
121         & multiDimAdvection, useEnergyConservingCoriolis,
122         & useJamartWetPoints
123    
124  C--   Elliptic solver parameters  C--   Elliptic solver parameters
125        NAMELIST /PARM02/        NAMELIST /PARM02/
126       & cg2dMaxIters, cg2dChkResFreq, cg2dTargetResidual, cg2dpcOffDFac       & cg2dMaxIters, cg2dChkResFreq, cg2dTargetResidual,
127         & cg2dTargetResWunit, cg2dpcOffDFac,
128         & cg3dMaxIters, cg3dChkResFreq, cg3dTargetResidual
129    
130  C--   Time stepping parammeters  C--   Time stepping parammeters
131        NAMELIST /PARM03/        NAMELIST /PARM03/
132       & nIter0, nTimeSteps, deltaT, deltaTmom, deltaTtracer, abEps, tauCD, rCD,       & nIter0, nTimeSteps, nEndIter, deltaT, deltaTmom, deltaTtracer,
133       & startTime, endTime, chkPtFreq, dumpFreq, taveFreq, deltaTClock,       & abEps, tauCD, rCD,
134       & pChkPtFreq, cAdjFreq, tauThetaClimRelax, tauSaltClimRelax,       & startTime, endTime, chkPtFreq,
135         & dumpFreq, taveFreq, deltaTClock, diagFreq,
136         & monitorFreq, pChkPtFreq, cAdjFreq,
137         & tauThetaClimRelax, tauSaltClimRelax, tauTr1ClimRelax,
138       & periodicExternalForcing, externForcingPeriod, externForcingCycle       & periodicExternalForcing, externForcingPeriod, externForcingCycle
139    
140  C--   Gridding parameters  C--   Gridding parameters
141        NAMELIST /PARM04/        NAMELIST /PARM04/
142       & usingCartesianGrid, delZ, dxSpacing, dySpacing, delX, delY,       & usingCartesianGrid, dxSpacing, dySpacing, delX, delY, delZ,
143       & usingSphericalPolarGrid, phiMin, thetaMin, rSphere,       & usingSphericalPolarGrid, phiMin, thetaMin, rSphere,
144       & l, m, n       & usingCurvilinearGrid,
145         & delP, delR, rkFac, Ro_SeaLevel, groundAtK1,
146         & delXfile, delYfile
147    
148  C--   Input files  C--   Input files
149        NAMELIST /PARM05/        NAMELIST /PARM05/
150       & bathyFile, hydrogThetaFile, hydrogSaltFile,       & bathyFile, topoFile, hydrogThetaFile, hydrogSaltFile,
151       & zonalWindFile, meridWindFile, thetaClimFile,       & zonalWindFile, meridWindFile,
152       & saltClimFile       & thetaClimFile, saltClimFile,
153         & surfQfile, EmPmRfile, surfQswfile,
154         & uVelInitFile, vVelInitFile, pSurfInitFile,
155         & dQdTFile
156    
157  C  C
158        _BEGIN_MASTER(myThid)        _BEGIN_MASTER(myThid)
159    
160    C     Defaults values for input parameters
161          CALL SET_DEFAULTS(
162         O   viscArDefault, diffKrTDefault, diffKrSDefault,
163         O   hFacMinDrDefault, delRdefault, rkFacDefault,
164         I   myThid )
165    
166    C--   Initialise "which vertical coordinate system used" flags.
167          zCoordInputData = .FALSE.
168          pCoordInputData = .FALSE.
169          rCoordInputData = .FALSE.
170          usingPCoords    = .FALSE.
171          usingZCoords    = .FALSE.
172          coordsSet       = 0
173    
174  C--   Open the parameter file  C--   Open the parameter file
175        OPEN(UNIT=scrUnit1,STATUS='SCRATCH')        OPEN(UNIT=scrUnit1,STATUS='SCRATCH')
176        OPEN(UNIT=scrUnit2,STATUS='SCRATCH')        OPEN(UNIT=scrUnit2,STATUS='SCRATCH')
177        OPEN(UNIT=modelDataUnit,FILE='data',STATUS='OLD',err=1,IOSTAT=errIO)        OPEN(UNIT=modelDataUnit,FILE='data',STATUS='OLD',
178        IF ( errIO .GE. 0 ) GOTO 2       &     IOSTAT=errIO)
179      1 CONTINUE        IF ( errIO .LT. 0 ) THEN
180         WRITE(msgBuf,'(A)')         WRITE(msgBuf,'(A)')
181       &  'S/R INI_PARMS'       &  'S/R INI_PARMS'
182         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
# Line 109  C--   Open the parameter file Line 188  C--   Open the parameter file
188         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
189         CALL MODELDATA_EXAMPLE( myThid )         CALL MODELDATA_EXAMPLE( myThid )
190         STOP 'ABNORMAL END: S/R INI_PARMS'         STOP 'ABNORMAL END: S/R INI_PARMS'
191      2 CONTINUE        ENDIF    
192    
193   1000 CONTINUE        DO WHILE ( .TRUE. )
194        READ(modelDataUnit,FMT='(A)',END=1001) RECORD         READ(modelDataUnit,FMT='(A)',END=1001) RECORD
195        IL = MAX(ILNBLNK(RECORD),1)         IL = MAX(ILNBLNK(RECORD),1)
196        IF ( RECORD(1:1) .NE. commentCharacter )         IF ( RECORD(1:1) .NE. commentCharacter )
197       &    WRITE(UNIT=scrUnit1,FMT='(A)') RECORD(:IL)       &     WRITE(UNIT=scrUnit1,FMT='(A)') RECORD(:IL)
198         WRITE(UNIT=scrUnit2,FMT='(A)') RECORD(:IL)          WRITE(UNIT=scrUnit2,FMT='(A)') RECORD(:IL)
199        GOTO 1000        ENDDO
200   1001 CONTINUE   1001 CONTINUE
201        CLOSE(modelDataUnit)        CLOSE(modelDataUnit)
202    
203  C--   Report contents of model parameter file  C--   Report contents of model parameter file
204        WRITE(msgBuf,'(A)')        WRITE(msgBuf,'(A)')
205       &'// ======================================================='       &'// ======================================================='
206        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit, SQUEEZE_RIGHT , 1)        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
207         &                    SQUEEZE_RIGHT , 1)
208        WRITE(msgBuf,'(A)') '// Model parameter file "data"'        WRITE(msgBuf,'(A)') '// Model parameter file "data"'
209        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit, SQUEEZE_RIGHT , 1)        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
210         &                    SQUEEZE_RIGHT , 1)
211        WRITE(msgBuf,'(A)')        WRITE(msgBuf,'(A)')
212       &'// ======================================================='       &'// ======================================================='
213        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
214       &  SQUEEZE_RIGHT , 1)       &  SQUEEZE_RIGHT , 1)
215        iUnit = scrUnit2        iUnit = scrUnit2
216        REWIND(iUnit)        REWIND(iUnit)
217   2000 CONTINUE        DO WHILE ( .TRUE. )
218         READ(UNIT=iUnit,FMT='(A)',END=2001) RECORD         READ(UNIT=iUnit,FMT='(A)',END=2001) RECORD
219         IL = MAX(ILNBLNK(RECORD),1)         IL = MAX(ILNBLNK(RECORD),1)
220         WRITE(msgBuf,'(A,A)') '>',RECORD(:IL)         WRITE(msgBuf,'(A,A)') '>',RECORD(:IL)
221         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit, SQUEEZE_RIGHT , 1)         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
222        GOTO 2000       &                    SQUEEZE_RIGHT , 1)
223          ENDDO
224   2001 CONTINUE   2001 CONTINUE
225        CLOSE(iUnit)        CLOSE(iUnit)
226        WRITE(msgBuf,'(A)') ' '        WRITE(msgBuf,'(A)') ' '
# Line 151  C--   Read settings from model parameter Line 233  C--   Read settings from model parameter
233        REWIND(iUnit)        REWIND(iUnit)
234    
235  C--   Set default "physical" parameters  C--   Set default "physical" parameters
236        DO K =1,Nr        viscAz   = UNSET_RL    
237         tRef(K) = 30.D0 - FLOAT(K)        viscAr   = UNSET_RL
238        ENDDO        viscAp   = UNSET_RL
239        gravity =   9.81 D0        diffKzT  = UNSET_RL
240        gBaro   = gravity        diffKpT  = UNSET_RL
241        rhoNil   = 999.8 D0        diffKrT  = UNSET_RL
242        rhoConst = 999.8 D0        diffKzS  = UNSET_RL
243        f0     = 1.D-4        diffKpS  = UNSET_RL
244        beta   = 1. _d -11        diffKrS  = UNSET_RL
245        viscAh = 1.d3        gBaro    = UNSET_RL
246        diffKhT= 1.0d3        rhoConst = UNSET_RL
247        diffKhS= 1.0d3        hFacMinDr           = UNSET_RL
248        viscAr = 1.d-3        hFacMinDz           = UNSET_RL
249        diffKrT= 1.d-5        hFacMinDp           = UNSET_RL
250        diffKrS= 1.d-5        READ(UNIT=iUnit,NML=PARM01) !,IOSTAT=errIO)
251        viscA4 = 0.        IF ( errIO .LT. 0 ) THEN
       diffK4T= 0.  
       diffK4S= 0.  
       GMmaxslope   = 1.d-2  
       GMlength     = 200.d3  
       GMalpha      = 0.  
       GMdepth      = 1000.  
       GMkbackground= 0.  
       GMmaxval     = 2500.d0  
       tAlpha=2.d-4  
       sBeta=7.4d-4  
       eosType='LINEAR'  
       buoyancyRelation='OCEANIC'  
       implicitFreeSurface = .TRUE.  
       rigidLid            = .FALSE.  
       freeSurfFac         = 1. _d 0  
       hFacMin             = 0. _d 0  
       hFacMinDr           = 0. _d 0  
       momViscosity        = .TRUE.  
       momAdvection        = .TRUE.  
       momForcing          = .TRUE.  
       useCoriolis         = .TRUE.  
       momPressureForcing  = .TRUE.  
       momStepping         = .TRUE.  
       tempStepping        = .TRUE.  
       saltStepping        = .TRUE.  
       metricTerms         = .TRUE.  
       implicitDiffusion   = .FALSE.  
       READ(UNIT=iUnit,NML=PARM01,IOSTAT=errIO,err=3)  
       IF ( errIO .GE. 0 ) GOTO 4  
     3 CONTINUE  
252         WRITE(msgBuf,'(A)')         WRITE(msgBuf,'(A)')
253       &  'S/R INI_PARMS'       &  'S/R INI_PARMS'
254         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
# Line 211  C--   Set default "physical" parameters Line 263  C--   Set default "physical" parameters
263         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
264         CALL MODELDATA_EXAMPLE( myThid )         CALL MODELDATA_EXAMPLE( myThid )
265         STOP 'ABNORMAL END: S/R INI_PARMS'         STOP 'ABNORMAL END: S/R INI_PARMS'
266     4  CONTINUE        ELSE
267        IF ( implicitFreeSurface ) freeSurfFac = 1. _d 0         WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM01 : OK'
268        IF ( rigidLid            ) freeSurfFac = 0. _d 0         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
269         &                    SQUEEZE_RIGHT , 1)
270          ENDIF
271          IF ( implicitFreeSurface ) freeSurfFac = 1.D0
272          IF ( rigidLid            ) freeSurfFac = 0.D0
273          IF ( gBaro .EQ. UNSET_RL ) gBaro=gravity
274          IF ( rhoConst .EQ. UNSET_RL ) rhoConst=rhoNil
275    C--   Momentum viscosity on/off flag.
276        IF ( momViscosity        ) THEN        IF ( momViscosity        ) THEN
277         vfFacMom = 1. _d 0         vfFacMom = 1.D0
278        ELSE        ELSE
279         vfFacMom = 0. _d 0         vfFacMom = 0.D0
280        ENDIF        ENDIF
281    C--   Momentum advection on/off flag.
282        IF ( momAdvection        ) THEN        IF ( momAdvection        ) THEN
283         afFacMom = 1. _d 0         afFacMom = 1.D0
284        ELSE        ELSE
285         afFacMom = 0. _d 0         afFacMom = 0.D0
286        ENDIF        ENDIF
287    C--   Momentum forcing on/off flag.
288        IF ( momForcing ) THEN        IF ( momForcing ) THEN
289         foFacMom = 1. _d 0         foFacMom = 1.D0
290        ELSE        ELSE
291         foFacMom = 0. _d 0         foFacMom = 0.D0
292        ENDIF        ENDIF
293    C--   Coriolis term on/off flag.
294        IF ( useCoriolis ) THEN        IF ( useCoriolis ) THEN
295         cfFacMom = 1. _d 0         cfFacMom = 1.D0
296        ELSE        ELSE
297         cfFacMom = 0. _d 0         cfFacMom = 0.D0
298        ENDIF        ENDIF
299    C--   Pressure term on/off flag.
300        IF ( momPressureForcing ) THEN        IF ( momPressureForcing ) THEN
301         pfFacMom = 1. _d 0         pfFacMom = 1.D0
302        ELSE        ELSE
303         pfFacMom = 0. _d 0         pfFacMom = 0.D0
304        ENDIF        ENDIF
305    C--   Metric terms on/off flag.
306        IF ( metricTerms ) THEN        IF ( metricTerms ) THEN
307         mTFacMom = 1. _d 0         mTFacMom = 1.D0
308        ELSE        ELSE
309         mTFacMom = 1. _d 0         mTFacMom = 0.D0
310          ENDIF
311    C--   z,p,r coord input switching.
312          IF ( viscAz .NE. UNSET_RL ) zCoordInputData = .TRUE.
313          IF ( viscAp .NE. UNSET_RL ) pCoordInputData = .TRUE.
314          IF ( viscAr .NE. UNSET_RL ) rCoordInputData = .TRUE.
315          IF ( viscAr .EQ. UNSET_RL )          viscAr = viscAz
316          IF ( viscAr .EQ. UNSET_RL )          viscAr = viscAp
317          IF ( viscAr .EQ. UNSET_RL )          viscAr = viscArDefault
318    
319          IF ( diffKzT .NE. UNSET_RL ) zCoordInputData  = .TRUE.
320          IF ( diffKpT .NE. UNSET_RL ) pCoordInputData  = .TRUE.
321          IF ( diffKrT .NE. UNSET_RL ) rCoordInputData  = .TRUE.
322          IF ( diffKrT .EQ. UNSET_RL )          diffKrT = diffKzT
323          IF ( diffKrT .EQ. UNSET_RL )          diffKrT = diffKpT
324          IF ( diffKrT .EQ. UNSET_RL )          diffKrT = diffKrTDefault
325    
326          IF ( diffKzS .NE. UNSET_RL ) zCoordInputData  = .TRUE.
327          IF ( diffKpS .NE. UNSET_RL ) pCoordInputData  = .TRUE.
328          IF ( diffKrS .NE. UNSET_RL ) rCoordInputData  = .TRUE.
329          IF ( diffKrS .EQ. UNSET_RL )          diffKrS = diffKzS
330          IF ( diffKrS .EQ. UNSET_RL )          diffKrS = diffKpS
331          IF ( diffKrS .EQ. UNSET_RL )          diffKrS = diffKrSDefault
332    
333          IF ( hFacMinDz .NE. UNSET_RL ) zCoordInputData = .TRUE.
334          IF ( hFacMinDp .NE. UNSET_RL ) pCoordInputData = .TRUE.
335          IF ( hFacMinDr .NE. UNSET_RL ) rCoordInputData = .TRUE.
336          IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr       = hFacMinDz
337          IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr       = hFacMinDp
338          IF ( hFacMinDr .EQ. UNSET_RL ) hFacMinDr       = hFacMinDrDefault
339    
340          IF ( ivdc_kappa .NE. 0. .AND. .NOT. implicitDiffusion ) THEN
341           WRITE(msgBuf,'(A,A)')
342         &  'S/R INI_PARMS: To use ivdc_kappa you must enable implicit',
343         &  ' vertical diffusion.'
344           CALL PRINT_ERROR( msgBuf , myThid)
345           STOP 'ABNORMAL END: S/R INI_PARMS'
346        ENDIF        ENDIF
347    
348        IF ( implicitFreeSurface .AND.  rigidLid ) THEN        IF ( implicitFreeSurface .AND.  rigidLid ) THEN
349           WRITE(msgBuf,'(A,A)')
350         & 'S/R INI_PARMS: Cannot select both implicitFreeSurface',
351         & ' and rigidLid.'
352           CALL PRINT_ERROR( msgBuf , myThid)
353           STOP 'ABNORMAL END: S/R INI_PARMS'
354          ENDIF
355          IF ( (implicSurfPress.NE.1. .OR. implicDiv2DFlow.NE.1.)
356         &    .AND. nonHydrostatic ) THEN
357           WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: nonHydrostatic',
358         & ' NOT SAFE with non-fully implicit Barotropic solver'
359           CALL PRINT_ERROR( msgBuf , myThid)
360           WRITE(msgBuf,'(A,A)') 'S/R INI_PARMS: To by-pass this',
361         &    'STOP, comment this test and re-compile ini_params'
362           CALL PRINT_ERROR( msgBuf , myThid)
363           STOP 'ABNORMAL END: S/R INI_PARMS'
364          ENDIF
365    
366          coordsSet = 0
367          IF ( zCoordInputData ) coordsSet = coordsSet + 1
368          IF ( pCoordInputData ) coordsSet = coordsSet + 1
369          IF ( rCoordInputData ) coordsSet = coordsSet + 1
370          IF ( coordsSet .GT. 1 ) THEN
371           WRITE(msgBuf,'(A)')
372         &  'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
373           CALL PRINT_ERROR( msgBuf , myThid)
374           STOP 'ABNORMAL END: S/R INI_PARMS'
375          ENDIF
376          IF ( rhoConst .LE. 0. ) THEN
377           WRITE(msgBuf,'(A)')
378         &  'S/R INI_PARMS: rhoConst must be greater than 0.'
379           CALL PRINT_ERROR( msgBuf , myThid)
380           STOP 'ABNORMAL END: S/R INI_PARMS'
381          ELSE
382           recip_rhoConst = 1.D0 / rhoConst
383          ENDIF
384          IF ( rhoNil .LE. 0. ) THEN
385         WRITE(msgBuf,'(A)')         WRITE(msgBuf,'(A)')
386       &  'S/R INI_PARMS: Cannot select implicitFreeSurface and rigidLid.'       &  'S/R INI_PARMS: rhoNil must be greater than 0.'
387         CALL PRINT_ERROR( msgBuf , myThid)         CALL PRINT_ERROR( msgBuf , myThid)
388         STOP 'ABNORMAL END: S/R INI_PARMS'         STOP 'ABNORMAL END: S/R INI_PARMS'
389          ELSE
390           recip_rhoNil = 1.D0 / rhoNil
391          ENDIF
392          IF ( HeatCapacity_Cp .LE. 0. ) THEN
393           WRITE(msgBuf,'(A)')
394         &  'S/R INI_PARMS: HeatCapacity_Cp must be greater than 0.'
395           CALL PRINT_ERROR( msgBuf , myThid)
396           STOP 'ABNORMAL END: S/R INI_PARMS'
397          ELSE
398           recip_Cp = 1.D0 / HeatCapacity_Cp
399          ENDIF
400          IF ( gravity .LE. 0. ) THEN
401           WRITE(msgBuf,'(A)')
402         &  'S/R INI_PARMS: gravity must be greater than 0.'
403           CALL PRINT_ERROR( msgBuf , myThid)
404           STOP 'ABNORMAL END: S/R INI_PARMS'
405          ELSE
406           recip_gravity = 1.D0 / gravity
407        ENDIF        ENDIF
408    C     Set globalFiles flag for READ_WRITE_FLD package
409          CALL SET_WRITE_GLOBAL_FLD( globalFiles )
410    C     Set globalFiles flag for READ_WRITE_REC package
411          CALL SET_WRITE_GLOBAL_REC( globalFiles )
412    C     Set globalFiles flag for READ_WRITE_REC package
413          CALL SET_WRITE_GLOBAL_PICKUP( globalFiles )
414    
415  C--   Elliptic solver parameters  C--   Elliptic solver parameters
416        cg2dMaxIters   = 150        READ(UNIT=iUnit,NML=PARM02) !,IOSTAT=errIO)
417        cg2dTargetResidual = 1. _d -7        IF ( errIO .LT. 0 ) THEN
       cg2dChkResFreq = 1  
       cg2dpcOffDFac  = 0.51 _d 0  
       READ(UNIT=iUnit,NML=PARM02,IOSTAT=errIO,err=5)  
       IF ( errIO .GE. 0 ) GOTO 6  
     5 CONTINUE  
418         WRITE(msgBuf,'(A)')         WRITE(msgBuf,'(A)')
419       &  'S/R INI_PARMS'       &  'S/R INI_PARMS'
420         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
# Line 274  C--   Elliptic solver parameters Line 429  C--   Elliptic solver parameters
429         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
430         CALL MODELDATA_EXAMPLE( myThid )         CALL MODELDATA_EXAMPLE( myThid )
431         STOP 'ABNORMAL END: S/R INI_PARMS'         STOP 'ABNORMAL END: S/R INI_PARMS'
432     6  CONTINUE        ELSE
433           WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM02 : OK'
434           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
435         &                    SQUEEZE_RIGHT , 1)
436          ENDIF    
437    
438  C--   Time stepping parameters  C--   Time stepping parameters
439        startTime      = 0.        rCD               = -1.D0
440        nTimeSteps     = 0        READ(UNIT=iUnit,NML=PARM03) !,IOSTAT=errIO)
441        endTime        = 0.        IF ( errIO .LT. 0 ) THEN
       nIter0         = 0  
       deltaT         = 0.  
       deltaTClock    = 0.  
       deltaTtracer   = 0.  
       deltaTMom      = 0.  
       abEps          = 0.01  
       pchkPtFreq     = 0.  
       chkPtFreq      = 3600.*25  
       dumpFreq       = 3600.*100  
       taveFreq       = 0.  
       writeStatePrec = precFloat32  
       nCheckLev      = 1  
       checkPtSuff(1) = 'ckptA'  
       checkPtSuff(2) = 'ckptB'  
       cAdjFreq       = -1. _d 0  
       rCD            = -1. _d 0  
       tauCD          = 0. _d 0  
       tauThetaClimRelax = 0. _d 0  
       doThetaClimRelax  = .FALSE.  
       tauSaltClimRelax  = 0. _d 0  
       doSaltClimRelax   = .FALSE.  
       periodicExternalForcing = .FALSE.  
       externForcingPeriod = 0.  
       externForcingCycle = 0.  
       READ(UNIT=iUnit,NML=PARM03,IOSTAT=errIO,err=7)  
       IF ( errIO .GE. 0 ) GOTO 8  
     7 CONTINUE  
442         WRITE(msgBuf,'(A)')         WRITE(msgBuf,'(A)')
443       &  'S/R INI_PARMS'       &  'S/R INI_PARMS'
444         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
# Line 321  C--   Time stepping parameters Line 453  C--   Time stepping parameters
453         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
454         CALL MODELDATA_EXAMPLE( myThid )         CALL MODELDATA_EXAMPLE( myThid )
455         STOP 'ABNORMAL END: S/R INI_PARMS'         STOP 'ABNORMAL END: S/R INI_PARMS'
456     8  CONTINUE        ELSE
457           WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM03 : OK'
458           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
459         &                    SQUEEZE_RIGHT , 1)
460          ENDIF  
461  C     Process "timestepping" params  C     Process "timestepping" params
462  C     o Time step size  C     o Time step size
463        IF ( deltaT       .EQ. 0. ) deltaT       = deltaTmom        IF ( deltaT       .EQ. 0. ) deltaT       = deltaTmom
# Line 360  C     o Convection frequency Line 496  C     o Convection frequency
496        IF ( cAdjFreq .LT. 0. ) THEN        IF ( cAdjFreq .LT. 0. ) THEN
497         cAdjFreq = deltaTClock         cAdjFreq = deltaTClock
498        ENDIF        ENDIF
499          IF ( ivdc_kappa .NE. 0. .AND. cAdjFreq .NE. 0. ) THEN
500           WRITE(msgBuf,'(A,A)')
501         &  'S/R INI_PARMS: You have enabled both ivdc_kappa and',
502         &  ' convective adjustment.'
503           CALL PRINT_ERROR( msgBuf , myThid)
504           STOP 'ABNORMAL END: S/R INI_PARMS'
505          ENDIF
506  C     o CD coupling  C     o CD coupling
507        IF ( tauCD .EQ. 0. _d 0 ) THEN        IF ( tauCD .EQ. 0.D0 ) THEN
508          tauCD = deltaTmom          tauCD = deltaTmom
509        ENDIF        ENDIF
510        IF ( rCD .LT. 0. ) THEN        IF ( rCD .LT. 0. ) THEN
511         rCD = 1. - deltaTMom/tauCD         rCD = 1. - deltaTMom/tauCD
512        ENDIF        ENDIF
513  C     o Temperature climatology relaxation time scale  C     o Temperature climatology relaxation time scale
514        IF ( tauThetaClimRelax .EQ. 0. _d 0 ) THEN        IF ( tauThetaClimRelax .EQ. 0.D0 ) THEN
515         doThetaClimRelax     = .FALSE.         doThetaClimRelax     = .FALSE.
516         lambdaThetaClimRelax = 0. _d 0         lambdaThetaClimRelax = 0.D0
517        ELSE        ELSE
518         doThetaClimRelax     = .TRUE.         doThetaClimRelax     = .TRUE.
519         lambdaThetaClimRelax = 1./tauThetaClimRelax         lambdaThetaClimRelax = 1./tauThetaClimRelax
520        ENDIF        ENDIF
521  C     o Salinity climatology relaxation time scale  C     o Salinity climatology relaxation time scale
522        IF ( tauSaltClimRelax .EQ. 0. _d 0 ) THEN        IF ( tauSaltClimRelax .EQ. 0.D0 ) THEN
523         doSaltClimRelax     = .FALSE.         doSaltClimRelax     = .FALSE.
524         lambdaSaltClimRelax = 0. _d 0         lambdaSaltClimRelax = 0.D0
525        ELSE        ELSE
526         doSaltClimRelax     = .TRUE.         doSaltClimRelax     = .TRUE.
527         lambdaSaltClimRelax = 1./tauSaltClimRelax         lambdaSaltClimRelax = 1./tauSaltClimRelax
528        ENDIF        ENDIF
529  C     o Time step count  C     o Tracer 1 climatology relaxation time scale
530        IF ( endTime .NE. 0 ) THEN        IF ( tauTr1ClimRelax .EQ. 0.D0 ) THEN
531         IF ( deltaTClock .NE. 0 ) nTimeSteps =         doTr1ClimRelax     = .FALSE.
532       &  INT((endTime-startTime)/deltaTClock)         lambdaTr1ClimRelax = 0.D0
533          ELSE
534           doTr1ClimRelax     = .TRUE.
535           lambdaTr1ClimRelax = 1./tauTr1ClimRelax
536          ENDIF
537    
538    C     o Start time
539          IF ( nIter0 .NE. 0 .AND. startTime .EQ. 0. )
540         &   startTime = deltaTClock*float(nIter0)
541    C     o nIter0
542          IF ( nIter0 .EQ. 0 .AND. startTime .NE. 0. )
543         &   nIter0 = INT( startTime/deltaTClock )
544    
545    C     o nTimeSteps 1
546          IF ( nTimeSteps .EQ. 0 .AND. nEndIter .NE. 0 )
547         &     nTimeSteps = nEndIter-nIter0
548    C     o nTimeSteps 2
549          IF ( nTimeSteps .EQ. 0 .AND. endTime .NE. 0. )
550         &     nTimeSteps = int(0.5+(endTime-startTime)/deltaTclock)
551    C     o nEndIter 1
552          IF ( nEndIter .EQ. 0 .AND. nTimeSteps .NE. 0 )
553         &     nEndIter = nIter0+nTimeSteps
554    C     o nEndIter 2
555          IF ( nEndIter .EQ. 0 .AND. endTime .NE. 0. )
556         &     nEndIter = int(0.5+endTime/deltaTclock)
557    C     o End Time 1
558          IF ( endTime .EQ. 0. .AND. nTimeSteps .NE. 0 )
559         &     endTime = startTime + deltaTClock*float(nTimeSteps)
560    C     o End Time 2
561          IF ( endTime .EQ. 0. .AND. nEndIter .NE. 0 )
562         &     endTime = deltaTClock*float(nEndIter)
563    
564    C     o Consistent?
565          IF ( nEndIter .NE. nIter0+nTimeSteps ) THEN
566           WRITE(msgBuf,'(A)')
567         & 'S/R INI_PARMS: nIter0, nTimeSteps and nEndIter are inconsistent'
568           CALL PRINT_ERROR( msgBuf , 1)
569           WRITE(msgBuf,'(A)')
570         & 'S/R INI_PARMS: Perhaps more than two were set at once'
571           CALL PRINT_ERROR( msgBuf , 1)
572           STOP 'ABNORMAL END: S/R INI_PARMS'
573          ENDIF
574          IF ( nTimeSteps .NE. int(0.5+(endTime-startTime)/deltaTClock) )
575         & THEN
576            WRITE(msgBuf,'(A)')
577         &  'S/R INI_PARMS: both endTime and nTimeSteps have been set'
578            CALL PRINT_ERROR( msgBuf , 1)
579            WRITE(msgBuf,'(A)')
580         &  'S/R INI_PARMS: but are inconsistent'
581            CALL PRINT_ERROR( msgBuf , 1)
582            STOP 'ABNORMAL END: S/R INI_PARMS'
583          ENDIF
584    
585    C     o Monitor (should also add CPP flag for monitor?)
586          IF (monitorFreq.LT.0.) THEN
587           monitorFreq=0.
588           IF (dumpFreq.NE.0.) monitorFreq=dumpFreq
589           IF (diagFreq.NE.0..AND.diagFreq.LT.monitorFreq)
590         &         monitorFreq=diagFreq
591           IF (taveFreq.NE.0..AND.taveFreq.LT.monitorFreq)
592         &         monitorFreq=taveFreq
593           IF (chkPtFreq.NE.0..AND.chkPtFreq.LT.monitorFreq)
594         &         monitorFreq=chkPtFreq
595           IF (pChkPtFreq.NE.0..AND.pChkPtFreq.LT.monitorFreq)
596         &         monitorFreq=pChkPtFreq
597           IF (monitorFreq.EQ.0.) monitorFreq=deltaTclock
598        ENDIF        ENDIF
 C     o Finish time  
       IF ( endTime .EQ. 0. ) endTime = FLOAT(nTimeSteps)*deltaTClock  
599    
600  C     o If taveFreq is finite, then we must make sure the diagnostics  C     o If taveFreq is finite, then we must make sure the diagnostics
601  C       code is being compiled  C       code is being compiled
602  #ifndef ALLOW_DIAGNOSTICS  #ifndef ALLOW_TIMEAVE
603        IF (taveFreq.NE.0.) THEN        IF (taveFreq.NE.0.) THEN
604          WRITE(msgBuf,'(A)')          WRITE(msgBuf,'(A)')
605       &  'S/R INI_PARMS: taveFreq <> 0  but you have'       &  'S/R INI_PARMS: taveFreq <> 0  but you have'
# Line 401  C       code is being compiled Line 607  C       code is being compiled
607          WRITE(msgBuf,'(A)')          WRITE(msgBuf,'(A)')
608       &  'not compiled the model with the diagnostics routines.'       &  'not compiled the model with the diagnostics routines.'
609          CALL PRINT_ERROR( msgBuf , 1)          CALL PRINT_ERROR( msgBuf , 1)
610          WRITE(msgBuf,'(A)')          WRITE(msgBuf,'(A,A)')
611       &  'Re-compile with:  #define ALLOW_DIAGNOSTICS  or  -DALLOW_DIAGNOSTICS'       &  'Re-compile with:  #define ALLOW_TIMEAVE',
612         &  '              or  -DALLOW_TIMEAVE'
613          CALL PRINT_ERROR( msgBuf , 1)          CALL PRINT_ERROR( msgBuf , 1)
614          STOP 'ABNORMAL END: S/R INI_PARMS'          STOP 'ABNORMAL END: S/R INI_PARMS'
615        ENDIF        ENDIF
# Line 410  C       code is being compiled Line 617  C       code is being compiled
617    
618  C--   Grid parameters  C--   Grid parameters
619  C     In cartesian coords distances are in metres  C     In cartesian coords distances are in metres
620        usingCartesianGrid = .TRUE.        rkFac = UNSET_RS
621        DO K =1,Nr        DO K =1,Nr
622         delZ(K) = 100. _d 0         delZ(K) = UNSET_RL
623        ENDDO         delP(K) = UNSET_RL
624        dxSpacing = 20. _d 0 * 1000. _d 0         delR(K) = UNSET_RL
       dySpacing = 20. _d 0 * 1000. _d 0  
       DO i=1,Nx  
        delX(i) = dxSpacing  
       ENDDO  
       DO j=1,Ny  
        delY(j) = dySpacing  
625        ENDDO        ENDDO
626  C     In spherical polar distances are in degrees  C     In spherical polar distances are in degrees
627        usingSphericalPolarGrid = .FALSE.        recip_rSphere  = 1.D0/rSphere
628        phiMin    = -5.0        dxSpacing = UNSET_RL
629        thetaMin  = 0.        dySpacing = UNSET_RL
630        rSphere   = 6370. * 1. _d 3        delXfile = ' '
631        recip_rSphere  = 1. _d 0/rSphere        delYfile = ' '
632        IF ( usingSphericalPolarGrid ) THEN        READ(UNIT=iUnit,NML=PARM04) !,IOSTAT=errIO)
633         dxSpacing = 1.        IF ( errIO .LT. 0 ) THEN
        dySpacing = 1.  
        DO I=1,Nx  
         delX(I) = dxSpacing  
        ENDDO  
        DO J=1,Ny  
         delY(J) = dySpacing  
        ENDDO  
       ENDIF  
   
       READ(UNIT=iUnit,NML=PARM04,IOSTAT=errIO,err=9)  
       IF ( errIO .GE. 0 ) GOTO 10  
     9 CONTINUE  
634         WRITE(msgBuf,'(A)')         WRITE(msgBuf,'(A)')
635       &  'S/R INI_PARMS'       &  'S/R INI_PARMS'
636         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
# Line 456  C     In spherical polar distances are i Line 645  C     In spherical polar distances are i
645         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
646         CALL MODELDATA_EXAMPLE( myThid )         CALL MODELDATA_EXAMPLE( myThid )
647         STOP 'ABNORMAL END: S/R INI_PARMS'         STOP 'ABNORMAL END: S/R INI_PARMS'
648    10  CONTINUE        ELSE
649           WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM04 : OK'
650           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
651         &                    SQUEEZE_RIGHT , 1)
652          ENDIF    
653    
654    C     X coordinate
655          IF ( delXfile .NE. ' ' ) THEN
656           IF ( delX(1) .NE. UNSET_RL .OR. dxSpacing .NE. UNSET_RL ) THEN
657             WRITE(msgBuf,'(A,A)') 'Too many specifications for delX:',
658         &   'Specify only one of delX, dxSpacing or delXfile'
659            CALL PRINT_ERROR( msgBuf , 1)
660            STOP 'ABNORMAL END: S/R INI_PARMS'
661           ELSE
662            _BEGIN_MASTER( myThid )
663            IF (readBinaryPrec.EQ.precFloat32) THEN
664             OPEN(37,FILE=delXfile,STATUS='OLD',FORM='UNFORMATTED',
665         &        ACCESS='DIRECT',RECL=WORDLENGTH*Nx)
666             READ(37,rec=1) delX
667    #ifdef _BYTESWAPIO
668                call MDS_BYTESWAPR4( Nx, delX )
669    #endif
670             CLOSE(37)
671            ELSEIF (readBinaryPrec.EQ.precFloat64) THEN
672             OPEN(37,FILE=delXfile,STATUS='OLD',FORM='UNFORMATTED',
673         &        ACCESS='DIRECT',RECL=WORDLENGTH*2*Nx)
674             READ(37,rec=1) delX
675    #ifdef _BYTESWAPIO
676                call MDS_BYTESWAPR8( Nx, delX )
677    #endif
678             CLOSE(37)
679            ENDIF
680            _END_MASTER(myThid)
681           ENDIF
682          ENDIF
683          IF ( dxSpacing .NE. UNSET_RL ) THEN
684           DO i=1,Nx
685            delX(i) = dxSpacing
686           ENDDO
687          ENDIF
688    C     Y coordinate
689          IF ( delYfile .NE. ' ' ) THEN
690           IF ( delY(1) .NE. UNSET_RL .OR. dySpacing .NE. UNSET_RL ) THEN
691             WRITE(msgBuf,'(A,A)') 'Too many specifications for delY:',
692         &   'Specify only one of delY, dySpacing or delYfile'
693            CALL PRINT_ERROR( msgBuf , 1)
694            STOP 'ABNORMAL END: S/R INI_PARMS'
695           ELSE
696            _BEGIN_MASTER( myThid )
697            IF (readBinaryPrec.EQ.precFloat32) THEN
698             OPEN(37,FILE=delYfile,STATUS='OLD',FORM='UNFORMATTED',
699         &        ACCESS='DIRECT',RECL=WORDLENGTH*Ny)
700             READ(37,rec=1) delY
701    #ifdef _BYTESWAPIO
702                call MDS_BYTESWAPR4( Ny, delY )
703    #endif
704             CLOSE(37)
705            ELSEIF (readBinaryPrec.EQ.precFloat64) THEN
706             OPEN(37,FILE=delYfile,STATUS='OLD',FORM='UNFORMATTED',
707         &        ACCESS='DIRECT',RECL=WORDLENGTH*2*Ny)
708             READ(37,rec=1) delY
709    #ifdef _BYTESWAPIO
710                call MDS_BYTESWAPR8( Ny, delY )
711    #endif
712             CLOSE(37)
713            ENDIF
714            _END_MASTER(myThid)
715           ENDIF
716          ENDIF
717          IF ( dySpacing .NE. UNSET_RL ) THEN
718           DO i=1,Ny
719            delY(i) = dySpacing
720           ENDDO
721          ENDIF
722    C
723        IF ( rSphere .NE. 0 ) THEN        IF ( rSphere .NE. 0 ) THEN
724         recip_rSphere = 1. _d 0/rSphere         recip_rSphere = 1.D0/rSphere
725        ELSE        ELSE
726         recip_rSphere = 0.         recip_rSphere = 0.
727        ENDIF        ENDIF
728    C--   Initialize EOS coefficients (3rd order polynomial)
 C     Initialize EOS coefficients (3rd order polynomial)  
729        IF (eostype.eq.'POLY3') THEN        IF (eostype.eq.'POLY3') THEN
730         OPEN(37,FILE='POLY3.COEFFS',STATUS='OLD',FORM='FORMATTED')         OPEN(37,FILE='POLY3.COEFFS',STATUS='OLD',FORM='FORMATTED')
731         READ(37,*) I         READ(37,*) I
732         IF (I.NE.Nr) THEN         IF (I.NE.Nr) THEN
733          WRITE(0,*) 'ini_parms: attempt to read POLY3.COEFFS failed'          WRITE(msgBuf,'(A)')
734          WRITE(0,*) '           because bad # of levels in data'       &  'ini_parms: attempt to read POLY3.COEFFS failed'
735            CALL PRINT_ERROR( msgBuf , 1)
736            WRITE(msgBuf,'(A)')
737         &  '           because bad # of levels in data'
738            CALL PRINT_ERROR( msgBuf , 1)
739          STOP 'Bad data in POLY3.COEFFS'          STOP 'Bad data in POLY3.COEFFS'
740         ENDIF         ENDIF
741         READ(37,*) (eosRefT(K),eosRefS(K),eosSig0(K),K=1,Nr)         READ(37,*) (eosRefT(K),eosRefS(K),eosSig0(K),K=1,Nr)
742         DO K=1,Nr         DO K=1,Nr
743          READ(37,*) (eosC(I,K),I=1,9)          READ(37,*) (eosC(I,K),I=1,9)
         write(0,'(i3,13f8.3)') K,eosRefT(K),eosRefS(K),eosSig0(K),  
      &                (eosC(I,K),I=1,9)  
744         ENDDO         ENDDO
745         CLOSE(37)         CLOSE(37)
746        ENDIF        ENDIF
747    C--   Check for conflicting grid definitions.
748        goptCount = 0        goptCount = 0
749        IF ( usingCartesianGrid )      goptCount = goptCount+1        IF ( usingCartesianGrid )      goptCount = goptCount+1
750        IF ( usingSphericalPolarGrid ) goptCount = goptCount+1        IF ( usingSphericalPolarGrid ) goptCount = goptCount+1
751        IF ( goptCount .NE. 1 ) THEN        IF ( usingCurvilinearGrid )    goptCount = goptCount+1
752          IF ( goptCount .GT. 1 ) THEN
753         WRITE(msgBuf,'(A)')         WRITE(msgBuf,'(A)')
754       &  'S/R INI_PARMS: More than one coordinate system requested'       &  'S/R INI_PARMS: More than one coordinate system requested'
755         CALL PRINT_ERROR( msgBuf , myThid)         CALL PRINT_ERROR( msgBuf , myThid)
756         STOP 'ABNORMAL END: S/R INI_PARMS'         STOP 'ABNORMAL END: S/R INI_PARMS'
757        ENDIF        ENDIF
758          IF ( goptCount .LT. 1 ) THEN
759           WRITE(msgBuf,'(A)')
760         &  'S/R INI_PARMS: No coordinate system requested'
761           CALL PRINT_ERROR( msgBuf , myThid)
762           STOP 'ABNORMAL END: S/R INI_PARMS'
763          ENDIF
764    C--   Make metric term settings consistent with underlying grid.
765        IF ( usingCartesianGrid ) THEN        IF ( usingCartesianGrid ) THEN
766         usingSphericalPolarMterms = .FALSE.         usingSphericalPolarMterms = .FALSE.
767         metricTerms = .FALSE.         metricTerms = .FALSE.
768         mTFacMom = 0         mTFacMom = 0.
769         useBetaPlaneF = .TRUE.         useBetaPlaneF = .TRUE.
770        ENDIF        ENDIF
771        IF ( usingSphericalPolarGrid ) THEN        IF ( usingSphericalPolarGrid ) THEN
772         useConstantF  = .FALSE.         useConstantF  = .FALSE.
773         useBetaPlaneF = .FALSE.         useBetaPlaneF = .FALSE.
774         useSphereF    = .TRUE.         useSphereF    = .TRUE.
775         omega         = 2. _d 0 * PI / ( 3600. _d 0 * 24. _d 0 )         usingSphericalPolarMterms = metricTerms
776         usingSphericalPolarMterms = .TRUE.        ENDIF
777         metricTerms = .TRUE.        IF ( usingCurvilinearGrid ) THEN
778         mTFacMom = 1         useSphereF    = .TRUE.
779          ENDIF
780    C--   p, z, r coord parameters
781          DO K = 1, Nr
782           IF ( delZ(K) .NE. UNSET_RL ) zCoordInputData = .TRUE.
783           IF ( delP(K) .NE. UNSET_RL ) pCoordInputData = .TRUE.
784           IF ( delR(K) .NE. UNSET_RL ) rCoordInputData = .TRUE.
785           IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delZ(K)
786           IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delP(K)
787           IF ( delR(K) .EQ. UNSET_RL ) delR(K) = delRDefault(K)
788           IF ( delR(K) .EQ. 0. ) THEN
789             WRITE(msgBuf,'(A,I4)')
790         &  'S/R INI_PARMS: No value for delZ/delP/delR at K = ',K
791             CALL PRINT_ERROR( msgBuf , 1)
792             STOP 'ABNORMAL END: S/R INI_PARMS'
793           ENDIF
794          ENDDO
795    C     Check for multiple coordinate systems
796          CoordsSet = 0
797          IF ( zCoordInputData ) coordsSet = coordsSet + 1
798          IF ( pCoordInputData ) coordsSet = coordsSet + 1
799          IF ( rCoordInputData ) coordsSet = coordsSet + 1
800          IF ( coordsSet .GT. 1 ) THEN
801           WRITE(msgBuf,'(A)')
802         &  'S/R INI_PARMS: Cannot mix z, p and r in the input data.'
803           CALL PRINT_ERROR( msgBuf , myThid)
804           STOP 'ABNORMAL END: S/R INI_PARMS'
805        ENDIF        ENDIF
806    
807  C--   Input files  C--   Input files
808        bathyFile       = ' '        READ(UNIT=iUnit,NML=PARM05) !,IOSTAT=errIO)
809        hydrogSaltFile  = ' '        IF ( errIO .LT. 0 ) THEN    
       hydrogThetaFile = ' '  
       zonalWindFile   = ' '  
       meridWindFile   = ' '  
       thetaClimFile   = ' '  
       saltClimFile    = ' '  
       READ(UNIT=iUnit,NML=PARM05,IOSTAT=errIO,err=11)  
       IF ( errIO .GE. 0 ) GOTO 12  
    11 CONTINUE  
        WRITE(msgBuf,'(A)')  
      &  'S/R INI_PARMS'  
        CALL PRINT_ERROR( msgBuf , 1)  
810         WRITE(msgBuf,'(A)')         WRITE(msgBuf,'(A)')
811       &  'Error reading numerical model '       &  'Error reading numerical model '
812         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
# Line 533  C--   Input files Line 818  C--   Input files
818         CALL PRINT_ERROR( msgBuf , 1)         CALL PRINT_ERROR( msgBuf , 1)
819         CALL MODELDATA_EXAMPLE( myThid )         CALL MODELDATA_EXAMPLE( myThid )
820         STOP 'ABNORMAL END: S/R INI_PARMS'         STOP 'ABNORMAL END: S/R INI_PARMS'
821    12  CONTINUE        ELSE
822           WRITE(msgBuf,'(A)') 'S/R INI_PARMS ; read PARM05 : OK'
823           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
824         &                    SQUEEZE_RIGHT , 1)
825          ENDIF      
826    
827    C
828    C--   Set factors required for mixing pressure and meters as vertical coordinate.
829    C     rkFac is a "sign" parameter which is used where the orientation of the vertical
830    C     coordinate (pressure or meters) relative to the vertical index (K) is important.
831    C     rkFac =  1 applies when K and the coordinate are in the opposite sense.
832    C     rkFac = -1 applies when K and the coordinate are in the same sense.
833    C     horiVertRatio is a parameter that maps horizontal units to vertical units.
834    C     It is used in certain special cases where lateral and vertical terms are
835    C     being combined and a single frame of reference is needed.
836          IF ( zCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
837           rkFac       = 1.D0
838           horiVertRatio = 1.D0
839          ENDIF
840          IF ( pCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
841           rkFac = -1.D0
842           horiVertRatio = Gravity * rhoConst
843          ENDIF
844          IF ( rCoordInputData .AND. rkFac .EQ. UNSET_RS ) THEN
845           rkFac =  1.D0
846           horiVertRatio = 1.D0
847          ENDIF
848          IF (buoyancyRelation.EQ.'ATMOSPHERIC')
849         &   horiVertRatio = Gravity * rhoConst                                    
850          IF ( rkFac .EQ. UNSET_RS ) rkFac=rkFacDefault
851          recip_rkFac = 1.D0 / rkFac
852          recip_horiVertRatio = 1./horiVertRatio
853          IF ( zCoordInputData ) usingZCoords = .TRUE.
854          IF ( pCoordInputData ) usingPCoords = .TRUE.
855    
856  C  C
857        CLOSE(iUnit)        CLOSE(iUnit)
Legend:
Removed from v.1.27  
changed lines
  Added in v.1.73
  ViewVC Help
Powered by ViewVC 1.1.22