!################################################################## !################################################################## !###### ###### !###### MODULE GLOBCST ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## MODULE globcst !----------------------------------------------------------------------- ! ! PURPOSE: ! ! This module contains most of the control parameters used by ARPS and ! related programs. ! ! It is based on include globcst.inc used by earlier version of ARPS. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! 10/01/1991 ! ! MODIFICATION HISTORY: ! ! 5/03/1992 (M. Xue and K. Droegemeier) ! Major restructuring. Some parameter names changed. ! ! 7/21/1992 (M. Xue) ! Added parameters for energy statistics calculations, initial ! random perturbation option, model run controls and additional ! history dump format. ! ! 3/10/1993 (M. Xue) ! Added parameters for advection options, and for terrain ! specifications. ! ! 2/12/94 (Y. Liu) ! Added parameters for the surface energy budget model. ! ! 9/10/94 (Weygandt & Y. Lu) ! Cleaned up documentation. ! ! 1/28/95 (G. Bassett) ! Added new parameter, buoyopt ! ! 02/07/1995 (Yuhe Liu) ! Deleted the table data array, veg(14) and added an input ! parameter, veg0. ! ! 2/2/96 (Donghai Wang & yuhe Liu) ! Added parameters for the map projection factor. ! ! 03/26/96 (Yuhe Liu) ! Added parameters for radiation computing. ! ! 4/2/96 (Donghai Wang, X. Song and M. Xue) ! Added parameters for the vertically implicit mixing. ! ! 5/7/96 (Donghai Wang and M. Xue) ! Added a new option for Rayleigh damping. ! ! 10/21/97 (Donghai Wang) ! Added two new options, buoy2nd and rhofctopt. ! ! 7/30/1998 (M. Xue) ! Converted to modeul from include file globcst.inc. ! !----------------------------------------------------------------------- IMPLICIT NONE SAVE INTEGER, PARAMETER :: mgrdmax = 100 ! Maximum of mgrid INTEGER :: mgrid ! Grid identifier INTEGER :: nestgrd ! Nested grid option ! = 0, no grid nesting ! = 1, grid nesting !----------------------------------------------------------------------- ! runname: A character string containing the pertinent information ! for this run. The initial characters before a blank space, a comma ! or otherwise the first 6 characters will be used to construct ! history data and restart data dump file names. This string is ! written into the output data file. !----------------------------------------------------------------------- CHARACTER :: runname*80 ! Name of this run INTEGER :: lfnkey ! Length of string (no. of characters) ! to be used for naming files (file name keys). INTEGER :: nocmnt ! Number of comment lines CHARACTER :: cmnt(50)*80 ! String of comments on this model run !----------------------------------------------------------------------- ! Model initialization control parameters. !----------------------------------------------------------------------- INTEGER :: runmod ! Model configuration option. ! = 1, model run in 3-D mode ! = 2, model run in 2-D mode (x-z plane) ! = 3, model run in 2-D mode (y-z plane) ! = 4, model run in 1-D mode (z direction) !----------------------------------------------------------------------- ! Model initialization control parameters. !----------------------------------------------------------------------- INTEGER :: initopt ! Model initialization option. ! = 1, self initialization (e.g. analytical field). ! = 2, restart run, initialize from a restart file. ! = 3, initialize from an external data set. INTEGER :: inibasopt ! Initialization option for base state fields. ! = 1, Initialize the base state using a single sounding. ! = 2, isentropic atmosphere ! = 3, isothermal atmosphere ! = 4, atmosphere with a constant static stability. ! = 5, an analytical thermodynamic sounding profile ! after Weisman and Klemp (1982, WMR) ! For options rather than 1, additional ! parameters are hardwired inside subroutine INIBASE. INTEGER :: viniopt ! Initialization option for base state ! velocity when inibasopt.ne.1. ! = 1, ubar= ubar0, vbar= vbar0 ! = 2, ubar0 and vbar0 are specified by user ! inside code. By default, ubar0=vbar0=0.0. INTEGER :: pt0opt ! Initial potential temperature perturbation ! option for initopt=1. ! = 1, initial perturbation is bubble-shaped, ! = 2, initial perturbation has a random distribution. ! = 3, random initial perturbation symmetric about x and y axis. ! = 4, initial perturbation as given in Skamarock and Klemp (1994) INTEGER :: inifmt ! Format option for the history dump type restart file ! See hdmpfmt. CHARACTER :: inifile*80 ! Filename of t-dependent variable initialization for initopt=3. CHARACTER :: inigbf*80 ! Filename of base&grid variable initialization for initopt=3. INTEGER :: inifunt ! Fortran I/O unit number for initial data input. INTEGER :: buoyopt ! Flag for turning off buoyancy terms ! (buoyopt = 0 dissables buoyancy). INTEGER :: buoy2nd ! Option for the second order terms in the linerized ! buoyancy terms ! = 1, including the 2nd-order terms ! = 0, only the 1st-order terms INTEGER :: rhofctopt ! Option for removing the density approximation ! in the pressure gradient force(PGF) terms ! = 1, removing the approximation ! = 0, using the approximation INTEGER :: bsnesq ! Bousinessq approximation? ! = 1, yes ! = 0, no REAL :: ptpert0 ! The magnitude of initial potential ! temperature perturbation (K) REAL :: pt0radx ! The radius of bubble in x-direction. REAL :: pt0rady ! The radius of bubble in y-direction. REAL :: pt0radz ! The radius of bubble in z-direction. REAL :: pt0ctrx ! The x coordinate of the bubble center. REAL :: pt0ctry ! The y coordinate of the bubble center. REAL :: pt0ctrz ! The z coordinate of the bubble center. REAL :: ubar0 ! User specified base state wind in x-dir. REAL :: vbar0 ! User specified base state wind in y-dir. !----------------------------------------------------------------------- ! Coordinate transforamtion and terrain option parameters: !----------------------------------------------------------------------- INTEGER :: crdtrns ! Coordinate transformation option ! = 0, no coordinate transforamtion at all ! = 1, with coordinate transforamtion. INTEGER :: ternopt ! Model terrain option. ! = 0, no terrain, the ground is flat; ! = 1, bell-shaped mountain; ! = 2, terrain data read from terrain data ! base (not implemented yet) INTEGER :: mntopt ! Option for choosing idealized mountain type. REAL :: hmount ! The mountain height (m) REAL :: mntwidx ! The half-width of bell-shaped mountain in x-dir. REAL :: mntwidy ! The half-width of bell-shaped mountain in y-dir. REAL :: mntctrx ! The x coordinate of the bell-shaped mountain center. REAL :: mntctry ! The y coordinate of the bell-shaped mountain center. CHARACTER :: terndta*80 ! Name of the terrain data file !----------------------------------------------------------------------- ! Environmental sounding file name and Fortran I/O input unit. !----------------------------------------------------------------------- CHARACTER :: sndfile*80 ! Name of the input sounding file. INTEGER :: sfunit ! Fortran I/O unit used by sounding input. !----------------------------------------------------------------------- ! Model grid parameters. !----------------------------------------------------------------------- REAL :: dx ! Grid spacing in x-direction in ! computational as well as in physical space (m) REAL :: dy ! Grid spacing in y-direction in ! computational as well as in physical space (m) INTEGER :: strhopt ! Vertical grid stretching option. ! = 0, no stretching in vectical. ! = 1, vertical stretching with function ! dzp=a+b*k**3 being used to specify ! dz as a function of k. ! = 2, vertical stretching with function ! dzp=c+a*tanh(b*k) being used to ! specify dz as a function of k. REAL :: dz ! Averaged grid spacing in vertical ! direction in transformed computataional space (m). REAL :: dzmin ! Minimun grid spacing in vertical ! direction in physcal space (m). REAL :: zrefsf! ! The reference height of the surface (ground level) (m) REAL :: dlayer1 ! The depth of the lower layer with uniform ! (dz=dzmin) vertical spacing (m) REAL :: dlayer2 ! The depth of the mid layer with stetched vertical spacing (m) REAL :: strhtune ! Tuning parameter for stretching option 2 ! A Value between 0.2 and 5.0 is appropriate. ! A larger value gives a more linear stretching. REAL :: zflat ! The height at which the grid levels ! becomes flat in the terrain-following ! coordinate transformation (m). REAL :: dxinv ! Reciprocal of dx (i.e. 1/dx) REAL :: dyinv ! Reciprocal of dy (i.e. 1/dy) REAL :: dzinv ! Reciprocal of dz (i.e. 1/dz) REAL :: xl ! Length of the physical model domain in ! x-direction (m). In uniform spacing case ,xl = (nx-3)*dx. REAL :: yl ! Length of the physical model domain in y-direction (m). ! In uniform spacing case ,yl = (ny-3)*dy. REAL :: zh ! Length of the physical model domain in ! z-direction (m). In uniform spacing case, zh = (nz-3)*dz. REAL :: ctrlat ! Latitude of the southwest corner of the ! physical model domain (deg. N) REAL :: ctrlon ! Longitude of the southwest corner of the ! physical model domain (deg. E) REAL :: xorig ! x-coordinate origin of the model physical domain REAL :: yorig ! y-coordinate origin of the model physical domain REAL :: zorig ! z-coordinate origin of the model physical domain REAL :: xgrdorg ! x-coordinate of model grid in absolute space REAL :: ygrdorg ! y-coordinate of model grid in absolute space !----------------------------------------------------------------------- ! Options for solution used in w equation ! Acoustic/sound wave speed definition used in the pressure equation. !----------------------------------------------------------------------- INTEGER :: vimplct ! Vertical implicit option for the w and p equations. ! = 1, implicit solution ! = 0, explicit solution INTEGER :: ptsmlstp ! Option to solve potential temperature eqaution. ! = 0, solve potential temperature eq. outside small time steps ! = 1, solve potential temperature eq. within small time steps INTEGER :: peqopt ! Option parameter of p-equation formulations ! = 1 or 2 (1 is as in 4.0 User's guide). INTEGER :: csopt ! Sound wave speed option used in the pressure equation. ! = 1, csound = cp/cv * rd * T ! = 2, csound = cp/cv * rd * T * csfactr ! = 3, csound = specfied constant. REAL :: csfactr ! A multiplication factor for the sound speed. REAL :: csound ! User specified constant sound speed. !----------------------------------------------------------------------- ! Model integration time step and integration control parameters. !----------------------------------------------------------------------- REAL :: dtbig ! Model integration time step. REAL :: dtsml ! Small time step for the acoustic wave integrations. REAL :: tstart ! Time when the model is initialized. REAL :: tstop ! Time when the integration is to be stoped. REAL :: curtim ! Current model time. REAL :: tacoef ! Time average weighting coefficient used ! in the vertically implicit solver. INTEGER :: nstep ! Current number of time steps taken. INTEGER :: nsteps ! The total number of time steps in this run. INTEGER :: nsmstp ! The number of small acoustic wave time ! integration steps per large time step. !----------------------------------------------------------------------- ! Parameters of the reference calendar day and real time !----------------------------------------------------------------------- INTEGER :: year ! Year of the reference calendar day INTEGER :: month ! Month of the reference calendar day INTEGER :: day ! Day of the reference calendar day INTEGER :: hour ! Hours of the reference real time INTEGER :: minute ! Minutes of the reference real time INTEGER :: second ! Seconds of the reference real time INTEGER :: jday ! Julian day starting from Jan. 1 of the year !----------------------------------------------------------------------- ! Robert-Asselin time filter coefficient for leapfrog time ! integration: !----------------------------------------------------------------------- REAL :: flteps ! Robert-Asselin time filter coefficient (non-dimensional). !----------------------------------------------------------------------- ! Option parameters for momentum and scalar advections. !----------------------------------------------------------------------- INTEGER :: madvopt ! Option parameter for momentum advection ! = 1, 2nd order centerd; ! = 2, 4th order centerd; ! = 3, ... INTEGER :: sadvopt ! Option parameter for scalar advection ! = 1, 2nd order centerd; ! = 2, 4th order centerd; ! = 3, Flux-corrected transport (FCT) INTEGER :: fctorderopt ! Option of the spactial order of accuracy of ! the FCT advection (sadvopt=4) and MPDCD ! advection schemes (sadvopt=5) ! = 1 2nd order ! = 2 4th order INTEGER :: fctadvptprt ! Option for FCT advection for potential temperature. ! Used only when sadvopt=4. ! = 0, FCT scheme is applied to ptbar+ptprt. ! = 1, FCT scheme is applied to ptprt only. ! = 2, FCT scheme is applied to ! ptbar+ptprt-min(ptbar+ptprt). !----------------------------------------------------------------------- ! Subgrid-scale turbulent mixing parameters. !----------------------------------------------------------------------- INTEGER :: tmixopt ! Control parameter for turbulent mixing options. ! = 0, zero turbulent mixing ! = 1, constant mixing coefficient, km = tmixcst ! = 2, Smagorinsky mixing coefficient. ! = 3, Smagorinsky mixing coefficient ! plus a constant coeffcient, tmixcst INTEGER :: trbisotp ! Option for isotropic subgrid scale turbulence. ! = 0, the turbulence is assumed to be anisotropic; ! = 1, the turbulence is assumed to be isotropic. INTEGER :: tkeopt ! Option for 1.5 order TKE formulation, used by tmixopt=4 ! = 1, after Wyngaard; ! = 2, after Deardroff; ! = 3, after Sun and Chang. INTEGER :: trbvimp ! Option parameter for the vertically implicit mixing ! = 0, vertical explicit ! = 1, vertical implicit REAL :: alfcoef ! Time average weighting coefficient used ! in the vertically implicit mixing REAL :: prantl ! Turbulent prantl number for Smagorinsky turb. REAL :: tmixcst ! Value of the constant mixing coefficient REAL :: kmlimit ! Constant used to set a limit of km !----------------------------------------------------------------------- ! Computational mixing parameters. !----------------------------------------------------------------------- INTEGER :: cmix2nd ! 2nd order computational mixing option. ! = 0, mixing off ! = 1, mixing on REAL :: cfcm2h ! 2nd order horizontal computational mixing ! coefficient caled by horizontal grid spacing (1/s). REAL :: cfcm2v ! 2nd order vertical computational mixing ! coefficient scaled by vertical grid spacing (1/s). REAL :: cfcmh2 ! Dimensional 2nd order horizontal ! computational mixing coefficient (m**2/s) REAL :: cfcmv2 ! Dimensional 2nd order vertical ! computational mixing coefficient (m**2/s) INTEGER :: cmix4th ! 4th order computational mixing option ! = 0, mixing off ! = 1, mixing on REAL :: cfcm4h ! 4th order horizontal computational mixing ! coefficient scaled by horizontal grid spacing (1/s). REAL :: cfcm4v ! 4th order vertical computational mixing ! coefficient scaled by vertical grid spacing (1/s). REAL :: cfcmh4 ! Dimensional 4th order horizontal ! computational mixing coefficient (m**2/s) REAL :: cfcmv4 ! Dimensional 4th order vertical ! computational mixing coefficient (m**2/s) REAL :: scmixfctr ! Reduction factor of the computational mixing ! coefficients for scalars relative to those of velocities. !----------------------------------------------------------------------- ! Upper level Rayleigh damping parameters. !----------------------------------------------------------------------- INTEGER :: raydmp ! Rayleigh damping option. ! = 0, damping off. ! = 1, damping the difference between the ! total and base state fields. ! = 2, damping the difference between the total and ! external fields defined in the EXBC file. REAL :: zbrdmp ! Height of the bottom of Rayleigh damping (sponge) layer (m). REAL :: cfrdmp ! Rayleigh damping coefficient (1/second). !----------------------------------------------------------------------- ! Acoustic wave divergence damping parameters. !----------------------------------------------------------------------- INTEGER :: divdmp ! Acoustic wave divergence damping option. ! = 0, divergence damping off ! = 1, divergence damping on REAL :: cdvdmph ! Divergence damping coefficient (m**2/s) in horizontal REAL :: cdvdmpv ! Divergence damping coefficient (m**2/s) in vertical REAL :: divdmpndh ! Non-dimensional divergence damping coefficient in horizontal REAL :: divdmpndv ! Non-dimensional divergence damping coefficient in vertical !----------------------------------------------------------------------- ! Moisture/ice microphysics parameters. !----------------------------------------------------------------------- INTEGER :: mphyopt ! Microphysics option. ! = 0, no microphysics. ! = 1, Kessler warm rain microphysics. ! = 2, Ice microphysics. INTEGER :: moist ! Moist processes option. ! = 0, dry run, all processes ! related to moisture are turned off ! = 1, moist processes are activated INTEGER :: cnvctopt ! Option for convective cumulus parameterization ! = 0, no convective parameterization ! = 1, Kuo scheme ! = 2, Kuo scheme and warm rain microphysics ! = 3, Kain and Fristch INTEGER :: ice ! Ice microphysics option. ! = 0, ice processes are turned off ! = 1, ice processes are activated. !----------------------------------------------------------------------- ! Cumulus physics parameters. !----------------------------------------------------------------------- REAL :: kffbfct ! A feedback factor for K-F scheme. REAL :: wcldbs ! Vertical motion needed at cloud base for convection. REAL :: confrq ! Frequency of conv param. updates in sec. REAL :: qpfgfrq ! Frequency of grid param. updates in sec. INTEGER :: idownd ! Downdraft flag ! = 0, no downdrafts. ! = 1, simple downdraft model. !----------------------------------------------------------------------- ! Microphysics constants. !----------------------------------------------------------------------- REAL, PARAMETER :: autotr = 1.0e-3 ! Cloud water mixing ratio threshold for ! autoconversion from cloud to rain water (kg/kg). REAL, PARAMETER :: autort = 1.0E-3 ! Autoconversion rate (1/seconds). REAL, PARAMETER :: accrrt = 2.2 ! Accretion rate (1/seconds). REAL, PARAMETER :: rho0 = 1.225 ! Reference density at surface (kg/meter**3). ) !----------------------------------------------------------------------- ! Saturation specific humidity parameters used in enhanced Teten's ! formula. (See A. Buck, JAM 1981) !----------------------------------------------------------------------- REAL, PARAMETER :: satfwa = 1.0007 REAL, PARAMETER :: satfwb = 3.46e-8 ! for p in Pa REAL, PARAMETER :: satewa = 611.21 ! es in Pa REAL, PARAMETER :: satewb = 17.502 REAL, PARAMETER :: satewc = 32.18 REAL, PARAMETER :: satfia = 1.0003 REAL, PARAMETER :: satfib = 4.18e-8 ! for p in Pa REAL, PARAMETER :: sateia = 611.15 ! es in Pa REAL, PARAMETER :: sateib = 22.452 REAL, PARAMETER :: sateic = 0.6 !----------------------------------------------------------------------- ! The x and y component of the domain translation speed ! (ground-relative). !----------------------------------------------------------------------- REAL :: umove ! Model domain translation speed in x direction. REAL :: vmove ! Model domain translation speed in y direction. !----------------------------------------------------------------------- ! Surface physics parameters. !----------------------------------------------------------------------- INTEGER :: sfcphy ! Surface physics option. ! = 0, no surface physics parameterization ! = 1, simple bulk aerodynamic formulation INTEGER :: landwtr ! Flag to indicate if land and water will be ! distinquished when the surface roughness ! length is calculated. If yes,the ! vagetation data (in sfcdtfl) has to be ! read in for all active options of sfcopt ! = 0, no distinquishing. ! = 1, surface roughness length is ! calculated differently for land and water. INTEGER :: cdhwtropt ! Option to set constant value of cdh (cdq) over water. INTEGER :: wtrexist ! Flag to check if any grid point is covered by water. REAL :: cdmlnd ! Bulk aerodynamic momentum drag coefficient. REAL :: cdmwtr ! Bulk aerodynamic momentum drag coefficient. REAL :: cdhlnd ! Bulk aerodynamic heat flux drag coefficient. REAL :: cdhwtr ! Bulk aerodynamic heat flux drag coefficient. REAL :: cdqlnd ! Bulk aerodynamic moisture flux drag coefficient. REAL :: cdqwtr ! Bulk aerodynamic moisture flux drag coefficient. ! cdm, cdh and cdq are all non-dimensional INTEGER :: pbldopt ! Option for determining PBL depth. REAL :: pbldpth0 ! User-specified PBL depth. REAL :: lsclpbl0 ! Length scale INTEGER :: sflxdis ! Option for distributing fluxes in PBL ! = 0, without flux distribution; ! = 1, with flux distribution for heat, ! moist, and momentum; ! = 2, with flux distribution for heat, moist. INTEGER :: tqflxdis ! Option for distributing heat and moisture ! fluxes qudratically in a specified depth dtqflxdis ! = 0, no distribution; ! = 1, with distribution REAL :: dtqflxdis ! Depth of flux distribution for tqflxdis=1, 200 m recommended. INTEGER :: sfcdiag ! Flag of surface diagnotic calculation INTEGER :: sfcdat ! Surface data input flag. INTEGER :: soilinit ! Soil model varialble initialization flag. INTEGER :: sfcunit ! Fortran I/O unit number for surface data input. REAL :: dtsfc ! Time step for surface energy budget model INTEGER :: nsfcst ! Number of steps for surface energy budget model INTEGER :: nstyp ! Number of soil types INTEGER :: styp ! Soil type INTEGER :: vtyp ! Vegetation type REAL :: lai0 ! Leaf Area Index REAL :: roufns0 ! Surface roughness REAL :: veg0 ! Vegetation fraction REAL :: ptslnd0 ! Initial ground-level soil potential temperature over land REAL :: ptswtr0 ! Initial ground-level soil potential temperature over water REAL :: tsoil0 ! Initial deep soil temperature REAL :: wetsfc0 ! Initial surface soil moisture REAL :: wetdp0 ! Initial deep soil moisture REAL :: wetcanp0 ! Initial canopy moisture REAL :: tsprt ! Offset of tsfc from sfc air temperature REAL :: t2prt ! Offset of tsoil from sfc air temperature REAL :: wgrat ! Saturation rate of sfc soil moisture REAL :: w2rat ! Saturation rate of deep soil moisture CHARACTER :: sfcdtfl*80 ! File name of surface data. CHARACTER :: soilinfl*80 ! File name of surface initialization data. !----------------------------------------------------------------------- ! Coriolis parameters !----------------------------------------------------------------------- INTEGER :: coriopt ! Coriolis term option. ! = 0, no coriolis terms ! = 1, coriolis terms involving horizontal wind ! = 2, coriolis terms involving both ! horizontal and vertical wind INTEGER :: coriotrm ! A flag indicating if total or perturbation ! wind is used in the Coriolis terms ! in the horizontal momemtum equations. ! = 1, total wind, ! = 2, perturbation wind. REAL :: latitud ! Model center latitude (degrees). REAL :: longitud ! Model center longitude (degrees). !----------------------------------------------------------------------- ! The name of the directory into which output files will be written: !----------------------------------------------------------------------- CHARACTER :: dirname*80 ! The name of output directory INTEGER :: ldirnam ! The length of the directory name string !----------------------------------------------------------------------- ! Model output parameters. !----------------------------------------------------------------------- INTEGER :: nfmtprt ! time steps between formatted print out REAL :: tfmtprt ! time interval between formatted print out INTEGER :: nmaxmin ! time steps between max/min statistics calc. REAL :: tmaxmin ! time interval between max/min statistics calc. INTEGER :: nenergy ! time steps between energy statistics calc. REAL :: tenergy ! time interval between energy statistics calc. INTEGER :: nrstout ! time steps between restart data dumps REAL :: trstout ! time interval between restart data dumps INTEGER :: totout ! Flag to indicate that the history dump ! contains total values instead of perturbation ! = 0, Perturbations ! = 1, total values INTEGER :: grdout ! Grid output option. ! = 0, no grid output ! = 1, grid output INTEGER :: basout ! Base state field output option. ! = 0, no base state fields output ! = 1, base state fields output INTEGER :: varout ! Mass & perturbation wind output option. ! = 0, no mass or perturbation wind output ! = 1, mass & perturbation wind output INTEGER :: mstout ! Moist variable output option. ! = 0, no moisture variables output ! = 1, qv, qc, qr, qi, qs and qh output. INTEGER :: rainout ! Rain variable output option. ! = 0, no rain variables output ! = 1, raing and rainc output. INTEGER :: prcout ! Precipitation rates output option. ! = 0, no rain variables output ! = 1, prcrate output. INTEGER :: iceout ! Ice variable ouptu option. ! = 0, no ice variables output ! = 1, qi, qs and qh output INTEGER :: tkeout ! TKE output option. ! = 0, no tke output ! = 1, output INTEGER :: trbout ! Turbulence field (km) output option. ! = 0, no km output ! = 1, km output INTEGER :: sfcout ! Surface variable output option. ! = 0, no surface variables output ! = 1, tsfc,tsoil,wetsfc,wetdp,wetsfc output INTEGER :: landout ! Surface property array output option. ! = 0, no surface property array output ! = 1, soiltyp,vegtyp,lai,roufns,veg output INTEGER :: radout ! Radiation arrays output option. ! = 0, no radiation arrays output ! = 1, radfrc, radsw, rnflx output INTEGER :: flxout ! Surface fluxes output option. ! = 0, no surface fluxes output ! = 1, usflx, vsflx, ptsflx, qvsflx output INTEGER :: exbcdmp ! Flag to dump the ARPS fields into an ! additional file to make the ARPS ! external boundary data set. ! = 0, no EXBC dumping ! = 1, EXBC dumping if lbcopt=2 INTEGER :: extdadmp ! Flag to dump the fields that contains ! external data array interpolated to the ! current model time. When lbcopt.ne.2, ! reset extdadmp to 0. ! = 0, no dumping ! = 1, dump. INTEGER :: qcexout ! Flag to dump qc into EXBC file ! = 0, no EXBC dumping ! = 1, EXBC dumping if exbcdmp=1 INTEGER :: qrexout ! Flag to dump qr into EXBC file ! = 0, no EXBC dumping ! = 1, EXBC dumping if exbcdmp=1 INTEGER :: qiexout ! Flag to dump qi into EXBC file ! = 0, no EXBC dumping ! = 1, EXBC dumping if exbcdmp=1 INTEGER :: qsexout ! Flag to dump qs into EXBC file ! = 0, no EXBC dumping ! = 1, EXBC dumping if exbcdmp=1 INTEGER :: qhexout ! Flag to dump qh into EXBC file ! = 0, no EXBC dumping ! = 1, EXBC dumping if exbcdmp=1 !----------------------------------------------------------------------- ! History data dump parameters. !----------------------------------------------------------------------- INTEGER :: hdmpopt ! History data dump option. ! = 1, dump with constant time interval ! = 2, dump at the time specified by user INTEGER :: hdmpfmt ! History data dump format option. ! = 0, No data dump is produced ! = 1, Unformatted binary data dump ! = 2, Formatted ascii data dump ! = 3, NCSA HDF format data dump ! = 4, Compressed binary data dump ! = 5, Dump for Savi3D visualization package ! = 6, binary allowing data point skip ! = 7, Uncompressed NetCDF format dump. ! = 8, Compressed NetCDF format dump. ! = 9, GrADS data dump ! = 10, GRIB data dump INTEGER :: grbpkbit ! Bit length to pack floating data in GRIB format REAL :: thisdmp ! time interval between history data dumps INTEGER :: nhisdmp ! time steps between history data dumps REAL :: tstrtdmp ! time to start history data dump INTEGER :: nstrtdmp ! time steps from beginning to tstrtdmp INTEGER :: numhdmp ! number of history dumps specified by user INTEGER, PARAMETER :: hdmpmax=100 ! maximum number of history dumps allowed in user specifying option REAL :: hdmptim(hdmpmax) ! times of model specified by user to dump out history data INTEGER :: hdmpstp(hdmpmax) ! steps from beginning of model, specified ! by user to dump out history data ! hdmpstp = nint( hdmptim/dtbig ) INTEGER :: nchdmp ! Fortran I/O channel for history data dump. CHARACTER :: hdmpfn*80 ! History data dump file name. INTEGER :: ldmpf ! Length of the history data dump file name !----------------------------------------------------------------------- ! More output control parameters !----------------------------------------------------------------------- INTEGER :: imgopt ! HDF image dumping option, 0 or 1 INTEGER :: nimgdmp ! time steps between HDF image dumps REAL :: timgdmp ! time interval between HDF image dumps INTEGER :: pltopt ! Graphic plotting option, 0 or 1 INTEGER :: nplots ! time steps between graphic plotting calls REAL :: tplots ! time interval between graphic plotting calls INTEGER :: filcmprs ! Option for output file compression, 0 or 1. !----------------------------------------------------------------------- ! Parameters for grid translation. !----------------------------------------------------------------------- INTEGER :: cltkopt ! Option for cell tracking. ! = 0, no cell-tracking, ! = 1, cell-tracking is turned on. INTEGER :: grdtrns ! Domain translation option. ! = 0, no change in domain translation. ! = 1, domain translation based on cell tracking. ! = 2, automatic domain translation via optimal procedure. REAL :: chkdpth ! Depth (m) over which the optimal domain translation ! speed is calculated for the automatic domain translation. REAL :: twindow ! Time window (s) between update of grid translation ! components for automatic domain translation. INTEGER :: nceltrk ! time steps between cell tracking calls REAL :: tceltrk ! time interval between cell tracking calls REAL :: tcrestr ! time required to restore cell center to domain center !----------------------------------------------------------------------- ! Debug information printing level. !----------------------------------------------------------------------- INTEGER :: lvldbg ! Level of debug information printing. ! =0, no printing; ! =1, Print variables in big t-step; ! =2, Print forcings in big t-step; ! =3, Print variables in small t-step; ! =4, Print forcings in small t-step; ! =5, Print individual forcing terms and misc. info; !----------------------------------------------------------------------- ! Model restart data dump parameters. !----------------------------------------------------------------------- INTEGER :: restrt ! Restart run option. ! = 0, non-restart run; ! = 1, restart run; INTEGER :: rstiunt ! Fortran I/O unit for restart data read. CHARACTER :: rstinf*80 ! File name of restart input data. INTEGER :: rstount ! Fortran I/O unit for restart data dump. CHARACTER :: rstoutf*80 ! File name of restart data dump. !----------------------------------------------------------------------- ! The Fortran I/O channel (unit) used for max/min statistics output. !----------------------------------------------------------------------- INTEGER :: nchmax ! Fortran I/O channel for max/min statistics output. !----------------------------------------------------------------------- ! The Fortran I/O channel (unit) used by energy statistics output !----------------------------------------------------------------------- INTEGER :: ncheng ! Fortran I/O channel for energy statistics output. !----------------------------------------------------------------------- ! Parameters used by Savi3D. !----------------------------------------------------------------------- INTEGER :: grafh ! Grafhandle used by savi3D INTEGER :: gridh ! Gridhandle used by savi3D INTEGER :: dsindex INTEGER :: gridid !----------------------------------------------------------------------- ! Parameters for map projections !----------------------------------------------------------------------- REAL, PARAMETER :: eradius=6371000.0 ! Mean earth radius in meters INTEGER :: mapproj ! Type of map projection in the model grid: ! modproj = 0 No projection. ! = 1 Polar Stereographic ! projection. ! = 2 Mercator projection. ! = 3 Lambert projection. REAL :: trulat1 ! 1st real true latitude of map projection REAL :: trulat2 ! 2nd real true latitude of map projection REAL :: trulon ! Real true longitude of map projection REAL :: sclfct ! Map scale factor. ! Distance on map, between two latitudes ! trulat1 and trulat2, is = (Distance on earth)*sclfct ! For ARPS model runs, generally this is 1.0 INTEGER :: mpfctopt ! Option parameter for map factor ! = 0, map factor set to 1 ! = 1, map factor calculated according to mapproj INTEGER :: mptrmopt ! Option parameter for map factor terms in momentum advection ! = 0, ignore the terms ! = 1, include the terms REAL :: swlatu ! latitude at southwest corner for u-grid REAL :: swlonu ! longitude at southwest corner for u-grid REAL :: nelatu ! latitude at northeast corner for u-grid REAL :: nelonu ! longitude at northeast corner for u-grid REAL :: swlatv ! latitude at southwest corner for v-grid REAL :: swlonv ! longitude at southwest corner for v-grid REAL :: nelatv ! latitude at northeast corner for v-grid REAL :: nelonv ! longitude at northeast corner for v-grid REAL :: swlats ! latitude at southwest corner for s-grid REAL :: swlons ! longitude at southwest corner for s-grid REAL :: nelats ! latitude at northeast corner for s-grid REAL :: nelons ! longitude at northeast corner for s-grid !----------------------------------------------------------------------- ! Parameters for radiation !----------------------------------------------------------------------- INTEGER :: radopt ! Option to switch on/off radiation physics ! = 0, No radiation physics; ! = 1, Radiation physics; INTEGER :: radstgr ! Option to switch on/off radiation ! computing at staggering points ! = 0, No staggering ! = 1, staggering INTEGER :: rlwopt ! Option for longwave schemes ! = 0, for high = .false. ! = 1, for high = .true. INTEGER :: raddiag ! Option to dump variables calculated in ! radiation transfer equations ! = 0, no dump ! = 1, dump to GrADS format INTEGER :: nradstp ! Time steps to update radiation calculation REAL :: dtrad ! Time interval to update radiation calculation REAL, PARAMETER :: cldh2m=5500.0 ! Height in meter that separates high cloud and middle cloud REAL, PARAMETER :: cldm2l=5500.0 ! Height in meter that separates middle cloud and low cloud INTEGER :: ict ! Vertical index of height cldh2m INTEGER :: icb ! Vertical index of height cldm2l CHARACTER :: arpsversion*20 ! A string containing the version number of ARPS. !----------------------------------------------------------------------- ! Accumulated CPU time usage by various processes or terms in the ! equations. !----------------------------------------------------------------------- REAL :: total_cpu_time, init_cpu_time,output_cpu_time, & rad_cpu_time,sfcphy_cpu_time,soil_cpu_time,tinteg_cpu_time, & tmix_cpu_time,cmix_cpu_time,raydmp_cpu_time,advuvw_cpu_time, & advs_cpu_time,coriol_cpu_time,buoy_cpu_time,tkesrc_cpu_time, & bc_cpu_time,qpfgrd_cpu_time,kuocum_cpu_time,kfcum_cpu_time, & smlstp_cpu_time,warmrain_cpu_time,linice_cpu_time, & nemice_cpu_time,qfall_cpu_time,misc_cpu_time END MODULE globcst