WESTERN REGION TECHNICAL ATTACHMENT NO. 98-21

WESTERN REGION TECHNICAL ATTACHMENT
NO. 98-21
June 2, 1998

A DESCRIPTION OF THE FIFTH GENERATION
PENN STATE/NCAR MESOSCALE MODEL (MM5)

James A. Nelson, Jr. - WRH-SSD/SWSFO-SLC

With computer power across the NWS increasing, real-time mesoscale modeling will be possible at each forecast office in the near future. There is the potential for increasing the skill of the meteorologists as model resolution and better terrain representation are achieved. Mesoscale models can be used to study or forecast fronts, land-sea breezes, mountain-valley circulations, marine pushes, and other mesoscale phenomena.

Model Description

The Fifth-Generation Penn State/NCAR Mesoscale Model Version 5 or MM5 (Grell et. al. 1995) originated from a mesoscale model developed by Richard Anthes at Penn State University, documented in Anthes and Warner (1978). The MM5 has grown into a community model with users continually improving and adding features. These include:

(a) a multiple-nesting capability,

(b) non-hydrostatic dynamics, and

(c) a four-dimensional data-assimilation capability as well as additional physics options.

The MM5 is best described as a hydrostatic or non-hydrostatic, sigma-coordinate model designed to simulate and predict mesoscale and regional-scale atmospheric circulations. The model itself is supported by several complimentary programs. These are known collectively as the MM5 modeling system (Fig.1). These are described in detail later in this Technical Attachment (TA).

The sigma-coordinate was first introduced by Phillips (1957). Sigma surfaces near the ground follow the terrain, and the higher-level sigma surfaces tend to approximate isobaric surfaces.

Sigma is defined as:



where P is the pressure, Pt is a specified constant top pressure, and Ps is the surface pressure. Sigma is equal to one at the surface (Ps) and zero at the top of the model (Ps). Normally, the vertical resolution of the model is much finer in the boundary layer with the number of sigma levels and their location determined by the discretion of the user and the phenomenon being modeled. The NCAR program INTERP (detailed later in this TA) performs the vertical interpolation from pressure levels to the sigma coordinate system (Fig. 2) of MM5.

The sigma coordinate system permits more efficient use of computer resources (Pielke and Martin, 1981). The advantages are two-fold: 1) The vertical advection terms are exactly zero at the top and bottom of the atmosphere 2) The effect of orography can be introduced without leading to either a) uncentered horizontal differences in the vicinity of mountains, or b) the assumption that the hypothetical flow patterns obtained by reduction to sea level actually exists. This method also simplifies the system of equations for the lower boundary conditions. But this simplification comes at an expense. The expense being that over steep orography, care needs to be taken to ensure the cancellation between the two components which comprise the pressure gradient force in the momentum equation. This results in more complicated vorticity and divergence equations, which is one reason why the eta coordinate system was developed.

The model's horizontal grid structure is based on the staggered Arakawa-Lamb B grid. In this scheme, the velocity variables are staggered with respect to the mass variables. It is designed to minimize errors associated with geostrophic adjustment and topographic forcing that occur in other horizontal grid domains. A sample subset of the B grid can be seen in Fig. 3. Each X (or "cross" point) represents a mass variable point (such as temperature or moisture) and each DOT represents both components of the wind. Each grid box consists of mass points surrounded by four velocity points. The distance between X or DOT points represents the model grid spacing. When data are input to the model, the preprocessor does the necessary interpolations to assure consistency with the grid. Horizontal velocities are represented in the middle of each model sigma layer, while vertical velocities are defined on the sigma levels. Finite differencing schemes within the model are crucially dependent on the grid staggering wherever gradients or averaging is required.

The model's topography is represented in a terrain file generated at NCAR using a USGS (United States Geological Survey) database for terrain and land use. The resolution of the source terrain and land use data are: 1) 111km, 2) 56km, 3) 19km, 4) 9km terrain and 19km land use, 5) .9km terrain over US and 19km land use. These terrain files are interpolated to the desired model grid spacing through the use of NCAR's TERRAIN program. This program will be detailed later in this TA.

Nesting

The MM5 allows for multiple nesting (Fig. 4) of up to nine domains at any given time. The resolution ratio is always 3:1 for two-way nesting. Two-way nesting means that the nest's input from the coarse mesh comes via its boundaries, while the feedback to the coarser mesh occurs over the inner domain. Domains can be moved during simulation as long as they are not the mother domain. The mother domain is defined as the model's outer-most domain. A two-way nest can be initialized from interpolating its data from the mother domain or from a high-resolution data set (e.g., Eta-Model). One-way nesting is another option. This method is not restricted to the 3:1 ratio. In this method, the model is first run to create output, which is then interpolated using any ratio to a new inner domain. In this process, a boundary file is also created. The boundary file typically is represented by the output frequency of the outer domain, and these data are then time-interpolated to supply the nest. The model is then rerun with these boundary conditions. The major disadvantages for the one-way nest versus the two-way nest are: a) no feedback occurs from the inner nest to the coarser mesh, b) coarser temporal resolution at the inner nest's boundaries.

Dynamics

Typically, models are run with the hydrostatic approximation because horizontal grid scales are comparable with or greater than the vertical scales. The hydrostatic approximation assumes that pressure is determined solely by the air mass above. Non-hydrostatic dynamics are needed when the horizontal scale is of the same magnitude on the vertical scale. This leads to an additional term, vertical acceleration, which determines the vertical pressure gradient such that hydrostatic balance is no longer in tact. This creates pressure perturbations from a reference state.

In a non-hydrostatic model, a reference state is needed for the calculation of these perturbations. The reference state is an idealized temperature profile in hydrostatic equilibrium. Usually, this is set to a typical sounding in the domain. An accurate sounding will reduce the pressure gradient force error associated with sloped coordinate surfaces over terrain.

Prognostic Variables

The primary prognostic variables in the MM5 are temperature, specific humidity, horizontal wind components, surface pressure, and turbulent kinetic energy. A split-explicit integration approach, which is also used in the Eta model, produces forecasts based on these quantities. This means that during each time step, after each process is computed, each of the primary variables is updated and the integration continues to the next advective time step when both the primary variables and any advective processes are updated. The fundamental time step in the MM5 is about three times the model resolution (e.g., 20km resolution requires at least a 60 second time step), with twice that time step required for advective processes.

Four-Dimensional Data Assimilation (FDDA)

FDDA allows the model to run with data that is entered into the model over an extended period of time which "nudges" the output toward the observations. This has a benefit that after a period of nudging, the model has been tuned to the data, while maintaining its dynamical balance. This method has been used in two primary ways; dynamical initialization and four-dimensional data sets. Dynamical initialization is where FDDA is used over a pre-forecast period to optimize the initial conditions for a real-time forecast. The second process, using four-dimensional data sets, is a method of producing dynamically balanced analyses that have a variety of uses from budget to tracer studies.

Meteorological data are horizontally and vertically interpolated via DATAGRID and INTERP. This interpolation does not provide mesoscale detail. Therefore, it may be necessary to supplement the data with surface or rawinsonde data. GOES-9, snow cover, and SST data are also available for integration into the initialization using the FDDA process.

The Modeling System

The MM5 modeling system is a group of programs run in order to prepare the data for model processing. The terrain program, DATAGRID, RAWINS (optional), and finally INTERP must be run before the MM5 is run.

This program horizontally interpolates the latitude/longitude terrain elevation and land use for the mesoscale domains chosen. These terrain files will be used by the DATAGRID program.

Land use is separated into 13 different categories (Table 1) for each grid point. These 13 land-use types determine surface properties such as; albedo, roughness length, long-wave emissivity, heat capacity, and moisture availability. Snow cover may also be added to this model from the NCEP snow cover data. Terrain is run only once each time a model domain is set up.

  • DATAGRID

This program horizontally and vertically interpolates gridded data to the mesoscale domains chosen by the user. These data are written onto pressure surfaces for interpolation by the programs INTERP or RAWINS. The interpolation is achieved by a bi-linear interpolation or Cressman analysis.

  • INTERP

The INTERP program transforms the data from the analysis program DATAGRID to the mesoscale model grids. This entails interpolating data onto the vertical and horizontal grid from pressure surfaces to sigma surfaces; diagnostic computation and data reformatting. This generation is needed for the model to receive its initial and boundary files. INTERP can also be used to place forecast data on pressure surfaces and interpolate these data to the inner domains for a one-way nest.

THE MM5 Model Physics Suite

The basic equations of the MM5 in their non-hydrostatic form momentum(x, y, and sigma), and thermodynamic equations are described by Anthes and Warner (1978). Prognostic equations also exist for water vapor and microphysical variables such as; cloud and precipitation.

Precipitation Physics

Precipitation schemes are usually divided into two groups; explicit and implicit. Explicit schemes handle resolved (grid scale) precipitation physics, while implicit schemes deal with non-resolved (sub-grid scale) precipitation physics. Each of these processes can occur in the same grid box throughout the model run and at the same time. Other options include a dry and "fake" dry run. The "fake" dry run has non-convective precipitation. Large-scale saturation is removed and rained out. There is no rain evaporation or explicit cloud prediction. The only effects removed are those of the latent heat process. The dry run has no moisture or water vapor.

Explicit Schemes

These schemes are usually activated whenever grid-scale saturation is reached. In other words, they treat resolved precipitation processes. In simplest terms, which is still used on large scales, remove super-saturation as precipitation and add latent heat to the thermodynamic processes. This is usually referred to as the warm rain scheme. More sophisticated schemes add additional variables such as cloud, rainwater, ice, and snow. In all these processes, the Marshall-Palmer size distributions are assumed for rain and snow.

Simple Ice (Dudhia)

This scheme adds ice phase processes to the warm rain scheme without adding computing time or memory. The explicit treatment of cloud water, rainwater, snow and ice allows for ice-phase processes below 0C, where the cloud water is treated as cloud ice and rain is treated as snow. As snow falls through the 0C, it is set to immediately melt into rain. Advection of ice or snow downward or of rain or cloud upward through this level also melts or freezes the particles. In both cases, the 0C isotherm is taken to be at a full sigma level. Melting occurs at the level immediately below this boundary and freezing above it.

Mixed-Phase (Reisner)

In the mixed-phase ice scheme, water and ice do not immediately freeze or melt above or below the 0C isotherm. Supercooled water can exist below 0C, as can unmelted snow exist above 0C. Separate arrays are used to store the vapor, cloud, rain, cloud ice, and snow variables. Immediate freezing of cloud water to cloud ice occurs at -40C and cloud ice immediately melts above 0C. This scheme adds supercooled water to the above scheme and allows for the slow melting of snow. Memory is added for cloud ice and snow. No graupel or riming processes are present in this scheme.

Goddard Microphysics

Included in the scheme are additional equations for the prediction of ice number concentrations and graupel.

Resiner Graupel

This scheme is based on the mixed-phase scheme but adds graupel and ice number concentration prediction equations.

Implicit Cumulus Parameterization Schemes

Anthes-Kuo Scheme

In this scheme, the amount of convection is determined by the vertically integrated moisture convergence. The feedback to the larger scale is determined with the help of vertical profiles of convective heating and moistening. A vertical eddy-flux divergence of water vapor is associated with cumulus convection. A portion of the vertically integrated moisture convergence is assumed to condense and precipitate, while a remaining portion is assumed to moisten the grid column. This portion is a function of the mean relative humidity of the column. If this mean relative humidity is equal to or less than 50%, this process always results in moistening the column. This scheme is mostly applicable to larger grid sizes (>30km). It tends to produce too much convective rainfall.

Arakawa-Schubert

This scheme allows for the interaction of a cumulus cloud ensemble with the large-scale environment. The large-scale environment is divided into the sub-cloud mixed layer and the Region above. The time changes of the environment are governed by the heat and moisture budget equations for the sub-cloud mixed layer and for the Region above, and by a prognostic equation for the depth of the mixed layer. In the environment above the mixed layer, the cumulus convection affects the temperature and moisture fields through downdrafts and detrainment of saturated air which evaporates in the environment. In the sub-cloud mixed layer, the cumulus convection does not act directly on the temperature and moisture fields, but it affects the depth of the mixed layer downdrafts. This scheme is mostly applicable to larger grid scales (>30km).

Grell (Modified Arakawa-Schubert) Scheme

This scheme is based on the rate of destabilization of a single cloud with updraft and downdraft fluxes and compensating motion determining the heating and moistening profile. In contrast to the original scheme, this scheme includes moist convective-scale downdrafts and allows for detrainment below the cloud. This cold downdraft outflow below the cloud changes the origins of the updraft air from a mixture of downdraft and the surrounding environment to one of air slightly above the outflow. Maximum buoyancy for both the updraft and downdraft is allowed in this scheme. This scheme is useful for smaller grid scales (10-30km) over mid-latitudes and tends to allow more resolved scale rainfall than convective rainfall.

Fritsch-Chappell

The formulation is based on the hypothesis that the buoyant energy available to a parcel, in combination with a prescribed period of time for the convection to remove that energy, can be used to regulate the amount of convection in a mesoscale numerical model grid element. This scheme is based on relaxation to a profile due to an updraft, downdraft, and subsidence region properties of a single cloud. Individual clouds are represented as the entraining moist updraft and downdraft plumes. The fraction of updraft condensate evaporated in moist downdrafts is determined from an empirical relationship between the vertical shear of the horizontal wind and precipitation efficiency. Vertical transports of horizontal momentum and warming are included in the parameterization. The influence of subsidence heating in a volume containing convection is computed explicitly by balancing the mass continuity equation. The convective mass flux removes 50% of the available buoyant energy in the relaxation time (approximately the Brunt-Vaisala Frequency). This scheme is useful for smaller grid scales (20-30km) due to the single-cloud assumptions and local subsidence.

Kain-Fritsch

This is similar to the Fritsch-Chappell scheme. However, a sophisticated cloud-mixing scheme modulates the two-way exchange of mass between cloud and environment (i.e. entrainment/detrainment) as a function of the buoyancy characteristics of various mixtures of clear and cloudy air. This scheme was formulated to ensure conservation of mass, thermal energy, total moisture, and momentum. In larger scale models over longer time frames, this assurance of conservation is needed. Grid scale considerations should be on the same order as the Fritsch-Chappell scheme.

Betts-Miller

In this scheme, non-precipitating shallow convection carries moisture upward and maintains low-level temperature inversions. Deep convection transports heat and moisture upward and produces precipitation. Model profiles of temperature and specific humidity are compared to profiles derived from numerous observations based on actual tropical convection. The model profiles are then relaxed (adjusted) toward the reference profiles. The precipitation is then deduced from the net negative change of specific humidity in the model's deep convective cloud. If the net change is positive (i.e., net evaporation occurred rather than condensation), no adjustment of the variables is made at that grid point. This scheme is useful at grid scales large than 30km, but not in severe weather simulations.

Parameterization of Shallow Convection

The shallow convection scheme serves to parameterize the planetary boundary layer (PBL) forced shallow non-precipitating convection as well as mid-tropospheric shallow convection caused by other sub-grid scale effects. The former is expected to transport boundary layer moisture to the level just above the boundary layer via convective bubbles. The bubbles are produced either by surface heating or moisture fluxes in a layer of strong lateral mixing. The bubbles are expected to rise without producing precipitation. Due to the strong lateral mixing, they usually only rise somewhere between 50-75mb. The latter process is much the same as the first, but the forcing mechanism is different. In this scheme, a bubble rises through several model layers forced by PBL fluxes or radiational cooling. The bubbles usually have large mixing, are non-precipitating, and have no convective downdrafts.

Planetary Boundary Layer Parameterization Schemes

Surface Flux Schemes

Ground temperatures are set from a choice of three different schemes. These schemes range from no scheme at all to a multi-soil layer-based scheme. The first option is one of no ground temperature prediction and the surface temperature is fixed.

Force/Restore (Blackadar) Scheme

The Force/Restore (Blackadar) scheme has a single slab (~10-20cm) and a fixed-temperature substrate. The slab temperature is based on the energy budget and depth assumed to represent the depth of the diurnal temperature variation. In this scheme, the slab temperature is set to be identical to the surface temperature of a real soil layer. In clear air, the amount of solar radiation absorbed by the slab, including multiple reflections of short waves, is proportional to the solar constant, albedo, zenith angle, and the short-wave transmissivity. All clear air transmissivities and backscattering coefficients are determined as a function of path length and precipitable water from a look-up table. Transmissivities are then adjusted for surface pressure. For cloudy skies, the effects of clouds on short-wave and downward long-wave radiation are incorporated into the parameterization scheme. Groups of sigma levels are chosen to correspond to low-, middle-, and upper-cloud layers. The attenuation of short-wave radiation by cloud is parameterized with absorption and scattering transmissivities. Cloud fraction is based upon relative humidity. The transfer of heat due to molecular conduction is defined by a heat transfer coefficient, angular momentum of the earth, and the temperature of the substrate.

The temperature is predicted in 1, 2, 4, 8, 16 cm layers with a fixed substrate below using a vertical diffusion equation. Thermal inertia is resolved in the same way as the force/restore scheme, but this scheme vertically resolves diurnal temperature variations. This allows for a more rapid response in surface temperature. Sensible-heat flux and surface-moisture flux are determined depending upon what PBL scheme is used. The available PBL schemes are listed below:

Bulk PBL

This scheme is suitable for coarse vertical resolution in the boundary layer (e.g. >250km vertical grid sizes). It is a very inexpensive choice. The surface-heat fluxes are defined by density and potential temperature at the lowest model layer and exchange coefficients. The surface moisture flux is defined by a moisture availability parameter which varies from 1 for a wet surface to 0 for a surface with no potential for evaporation. This parameter is determined via land-use. The surface momentum flux is determined as a function of the drag coefficient, density, and wind velocity. It has two stability regimes (stable and unstable).

High-resolution Blackadar PBL

This scheme is used to forecast the vertical mixing of horizontal winds, potential temperature, mixing ratio, cloud water, and ice. The surface heat and moisture fluxes are computed from similarity theory. The surface heat flux is a function of surface roughness, frictional velocity, potential temperature, density, and non-dimensional stability parameters, which are a function of the Bulk Richardson Number. This scheme is suitable for high-resolution PBL. It has five layers in the lowest km and the surface layer is <100m thick. There are four stability regimes (stable, mechanically driven turbulence, forced convection, and free-convection).

Burk-Thompson PBL

This scheme adds the ability to predict turbulent kinetic energy for use in vertical mixing to counteract the inadequate capability of previous schemes in a well-mixed atmosphere. The scheme vertically nests a high-resolution grid within the framework of a more coarsely spaced regional model. The physics computations are performed on the fine grid, while the atmospheric dynamic calculations are done at the coarser resolution. The two grids interact with each time step.

MRF

In order to explain the MRF scheme, one must know the history of vertical diffusion modeling at NCEP. The vertical diffusion scheme, based on local gradients of wind and potential temperature (local-K approach), was in operation at NCEP during the 1980s and early 1990s. In this scheme, the diffusivity coefficients are parameterized as functions of the local Richardson number. A shortcoming of this scheme is that the transport of mass and momentum in the PBL is mostly accomplished by the largest eddies, which should be modeled by the bulk properties and not the local ones. Thus, the scheme cannot handle conditions when the atmosphere is well mixed because of the "countergradient fluxes". In order to contain this deficiency, the nonlocal K approach was introduced to represent large eddy turbulence within a well-mixed boundary layer.

Radiation Schemes

Atmospheric radiation options in the model range from a scheme that considers no mean tendency being applied to the atmospheric temperature to a scheme which incorporates multiple spectral bands in shortwave and longwave. The atmospheric radiation scheme preferred by most modelers is represented via a longwave and shortwave scheme that interacts with the atmosphere, cloud, and precipitation fields and with the surface.

Longwave Radiation

Longwave absorption by water vapor is strongly spectral in character, and the method employed uses a precalculated emissivity function which represents the frequency-integrated absorption spectrum of water vapor, weighted by a suitable function (upward and emissivity as a function of the water vapor path). Cloud water is assumed to have a constant absorption coefficient which is slightly different for upward and downward radiation. Ice cloud is assumed to be a black body and has an absorption rate defined using the cross-sectional area. This is also done for rain and snow. These effects are combined with the water vapor absorption and the transmissivities are multiplied. Carbon dioxide is less easily treated. The spectrum is divided into a carbon-dioxide band and a non-carbon-dioxide Region. The former requires overlapping of the carbon dioxide transmissivity function while the latter does not. The relative weights of these two regions are slightly temperature dependent. The upward and downward longwave radiation are defined on model full sigma levels.

Shortwave Radiation

The downward component of shortwave radiation is a function of the solar zenith angle, clouds, and clear air. The solar zenith angle affects the path length while the clouds affect the albedo and absorption, and the clear air flux is determined by scattering and water-vapor absorption.

The radiation schemes available in model simulations are listed below:

None

No mean tendency is applied to atmospheric temperature and is considered unrealistic in long-term simulations.

Simple Cooling

Atmospheric cooling rates depend just on temperature. There is no cloud interaction or diurnal cycle.

Surface Radiation

This is used with the above two options. It provides diurnally varying shortwave and longwave flux at the surface for use in the ground energy budget. These fluxes are calculated based on atmospheric column-intergrated water vapor and low/middle/high cloud fraction estimated from relative humidity.

Cloud-Radiation Scheme

This scheme is sophisticated enough to account for longwave and shortwave interactions with explicit cloud and clear-air. As well as atmospheric temperature tendencies, this scheme provides surface radiation fluxes.

CCM2 Radiation Scheme

This scheme has multiple spectral bands in shortwave and longwave, but closed treated simply based on relative humidity. This is suitable for larger grid scales, and probably more accurate for long time integrations. This scheme also provides radiative fluxes at the surface.

References

Anthes, R. A., 1977: A cumulus parameterization scheme utilizing a one-dimensional cloud model. Mon. Wea. Rev., 105, 270-286.

______, and T. T. Warner, 1978: Development of hydrodynamic models suitable for air pollution and other meso-meteorological studies. Mon. Wea. Rev., 106, 1045-1078.

Arakawa, A., and W. H. Schubert, 1974: Interaction of a cumulus cloud ensemble with the large-scale environment. Part I. J. Atmos. Sci., 31, 674-701.

Betts, A. K., and M.J. Miller, 1993: The Betts-Miller scheme. The representation of cumulus convection in numerical models, K. A. Emanuel and D. J. Raymond, Eds., Amer. Meteor. Soc., 246pp.

Brankovic, C., 1981: A transformed isentropic coordinate and its use in an atmospheric model. Mon Wea. Rev., 109, 2029-2039.

Burk, S. D., and W. T. Thompson, 1989: A vertically nested regional numerical prediction model with second-order-closure physics. Mon. Wea. Rev., 117, 2305-2324.

Dudhia, J., 1989: A nonhydrostatic version of the Penn State-NCAR mesoscale model. Mon. Wea. Rev., 121, 1493-1513.

Fritsch, J. M., and C. F. Chappell, 1980: Numerical prediction of convectively driven mesoscale pressure systems. Part I: Convective parameterization. J. Atmos. Sci., 37, 1722-1733.

Grell, G. A., 1993: Prognostic evaluation of assumptions used by cumulus parameterizations. Mon. Wea. Rev., 121, 764-787.

Hong, S. -Y., and H.-L. Pan, 1996: Nonlocal boundary layer vertical diffusion in a medium-range forecast model. Mon. Wea. Rev., 124, 2322-2339.

Kain, J. S., and J. M. Fritsch, 1993: Convective parameterization for mesoscale models: The Kain-Fritsch scheme. The representation of cumulus convection in numerical models, K. A. Emanuel and D. J. Raymond, Eds., Amer. Meteor. Soc., 246pp.

Kasahara, A., 1974: Various vertical coordinate systems used for numerical weather prediction, Mon. Wea. Rev., 102, 509-522.

Molinari, J., and M. Dudek, 1992: Parameterization of convective precipitation in mesoscale numerical models: A critical review. Mon. Wea. Rev., 120, 326-344.

Phillips, N. A., 1957: A coordinate system having some special advantages for numerical forecasting. J. Meteor., 14, 184-185.

Pielke, R. A., 1984: Mesoscale Meteorological Modeling, Academic Press, Inc., San Diego, 455pp.

______, and C. L. Martin, 1981: The derivation of a terrain-following coordinate system for use is a hydrostatic model. J. Atmos. Sci., 38, 1707-1713.

Zhang, D. L., J. S. Kain, J. M. Fritsch, and K. Gao, 1994: Comments on "Parameterization of convective precipitation in mesoscale numerical models: A critical review". Mon. Wea. Rev., 122, 2222-2231.