Partners header

ECOSIM Telematics Applications Project:
ECOSIM Deliverable D05.02

ECOSIM User Documentation

Programme name Telematics Applications Programme
Sector Environment
Project Number EN1006
Project Acronym ECOSIM
Project title Ecological and environmental monitoring and simulation system for management decision support in urban areas
Deliverable number D05.02
Deliverable title ECOSIM USer Documentation
Workpackage contributing to Deliverable WP 5
Authors P. Mieth, S. Unger, GMD, N. Moussiopoulos, K. Karatzas, AUT, K. Fedra, ESS, D. Fatta, M. Loizidou, TUA, L. Perivoliotis, A. Lascaratos, UA
Project Manager  
Name Dr. Kurt FEDRA
Organisation Environmental Software and Services GmbH
Address Kalkgewerk 1
Country-Code City A-2352 Gumpoldskirchen, AUSTRIA
Telephone +43 2253 633 050
Fax +43 2253 655 059


Table of Contents

    Technical abstract

    Executive summary

    1 Introduction

      1.1 General
      1.2 Objective
      1.3 Contents

    2 General features of the ECOSIM system and conceptual overview

    3 The ECOSIM main server

      3.1 Geographical information system (GIS)
      3.2 AirWare
      3.3 Systems architecture
      3.4 Hardware requirements and operating system
      3.5 Networking
      3.6 Auxiliary tools
      3.7 Data requirements
      3.7.1 GIS data
      3.7.2 Object data

    4 The Models - technical description

      4.1 The Mesoscale meteorological model MEMO
      4.1.1 Model equations
      4.1.2 Transformation to terrain-influenced coordinates
      4.1.3 Numerical solution of the equation system
      4.1.4 Parameterisations
      4.1.5 Initial and boundary conditions
      4.1.6 Grid definition
      4.1.7 Nesting facility
      4.2 The Dispersion and chemical model DYMOS
      4.2.1 The Chemical mechanism
      4.2.2 Treatment of surface removal processes
      4.3 The Forecast system
      4.3.1 The Ozone model REGOZON
      4.3.2 The Multilayer dispersion model MUSE
      4.4 The Groundwater model
      4.4.1 Simulation of the hydraulic conditions
      4.4.2 Simulation of the pollutants transport
      4.5 The Ocean water model
      4.5.1 Main features
      4.5.2 Boundary conditions

    5 Model input requirements (air domain)

      5.1 Digital elevation map and land use data
      5.2 Emission inventory
      5.3 Initial and background data
      5.3.1 Meteorological data
      5.3.2 Background concentrations

    6 Model output (air domain)

    7 Computer user notes and installation guide

    8 Demonstrator integration

      8.1 REGOZON forecasting model and the BLUME network
      8.1.1 REGOZON forecasting model: scenario editing
      8.1.2 REGOZON forecasting model: scenario display
      8.1.3 BLUME monitoring data acquisition
      8.2 The MEMO/MUSE meteorological and air quality models
      8.2.1 MEMO meteorological model: scenario editing
      8.2.2 MEMO/MUSE: scenario display
      8.3 The POM (Princeton Ocean Model) coastal water quality model
      8.4 The ECOSIM groundwater quality model

    9 Bibliography

Technical Abstract

ECOSIM is a support system to investigate and forecast pollution levels in urban areas. By allowing the effects of industrial developments or new roads to be quickly and easily examined, public authorities can use ECOSIM to ensure that their urban plans fully consider environmental impacts. ECOSIM is also able to forecast pollution levels - typically over the next 24 hours.

ECOSIM is based on sophisticated computer models such as those which allow ozone levels to be calculated from road traffic emissions. It combines these models with up-to-date measurements of current meteorological conditions and pollution by linking directly into local databases and pollution measurement stations. ECOSIM also links the various domains such as surface water, coastal water and air to ensure that as complete a picture as possible can be predicted of environmental conditions.

ECOSIM makes use of high performance computers whenever it needs to and uses the latest methods in handling maps and similar data to ensure that its results can be easily translated into practical measures for pollution control.

Executive Summary

This document, the user documentation, will enable an end-user of the ECOSIM system to decide if the application of ECOSIM is useful for solving his problems and will enable him to use the ECOSIM system through detailed guidance.

This deliverable describes the main model parts, its theoretical background, its capabilities and requirements as well as input and output data and their interfaces. In addition, it introduces the main server and explains his basic features.

The document gives an introduction and shows the principle simulation system design of ECOSIM. The ECOSIM main server is described which includes a Geographical Information System (GIS) and a data base, enables networking and remote control of simulation runs. It connects monitoring information with model data for initialisation, background information or validation purposes. The design of the models and the theoretical description is given. It includes opportunities of ozone forecasting with different complexity of the models. The document presents the set of input data which is required for the use of ECOSIM. The main data set for the air domain consists of a digital elevation map, a land use map, an emission inventory and initial and boundary values. The document describes the main features of the ocean and groundwater models (POM, MODFLOW) and its input and output facilities.

1. Introduction

1.1 General

This report has been prepared by GMD-Forschungszentrum Informationstechnik GmbH (GMD) with main input from Aristotle University of Thessaloniki (AUT) and Environmental Software & Services GmbH (ESS). It presents the user documentation, deliverable number D05.02, for the ECOSIM project (Project No EN 1006).

1.2 Objective

This document, the user documentation, will enable an end-user of the ECOSIM system to decide if the application of ECOSIM is useful for solving his problems and will enable him to use the ECOSIM system through detailed guidance. This first issue gives a principle overview of the system. It describes the main model parts, its theoretical background, its capabilities and requirements. In addition, it introduces the main server and explains its basic features.

1.3 Contents

Section 2 gives an introduction and shows the principle simulation system design of ECOSIM. In Section 3 the ECOSIM main server is described which includes a Geographical Information System (GIS) and a data base, enables networking and remote control of simulation runs. It connects monitoring information with model data for initialisation, background information or validation purposes. In Section 4 the design of the models and the theoretical description is given. It includes opportunities of ozone forecasting with different complexity of the models. Water domain models, such as the ground water model MODFLOW and the ocean water model POM, are described. Section 5 presents the set of input data which is required for the use of ECOSIM. The main data set in the air domain consists of a digital elevation model, a land use map, an emission inventory and initial and boundary values. Section 6 describes the quantities computed by the different models. Some computer related features and requirements, general hardware and software requirements as well as run time considerations are given in Section 7. Section 8 describes the main features of the demonstrator for the different application domains and it's usage.

2 .General Features of the ECOSIM System and Conceptual Overview

The ECOSIM project builds up a prototype of an environmental simulation and monitoring system including advanced visualisation and analysing capabilities for urban areas. The main features are designed for planning, analysing and forecast purposes in the field of air quality management and decision support. It is based on the idea of connecting state of the art air pollution models devised for complex terrain with graphical information systems (GIS) and monitoring networks. In order to take advantage of the computational power of different high performance computation centres and to simplify the implementation for the users, the different submoduls of the software system are using a client-server architecture connected via a http server.

Air pollutant dispersion models are extremely valuable for environmental planning and management, as they may be used to study complex source/receptor relationships. In this sense, these models may be used:

  • for the optimization of air pollution abatement strategies either locally or in a larger region;
  • for short-term predictions of smog episodes in highly polluted areas.

Coupling of model computations with on-line measuring networks and a data base allows analysing and evaluation of the model results as well as initialisation and determination of the actual background values. The use of high performance hardware platforms (e.g. parallel computers) allowed the formulation of more accurate air pollutant dispersion models which has been used up to now for research purposes only.

The ECOSIM simulation system is designed in a client-server mode. In general, the different model components run at different locations on different computers. They are connected with each other via the Internet and controlled by a main server which must be located at the user`s office. The main server provides the access to all parts of the software system and allows

  • model configuration and set up of the model domain
  • generating and modifying the model input data
  • operating and controlling simulation processes
  • advanced visualization with a geographic information system (GIS)
  • data base access
  • data analysis
  • on-line coupling with a monitoring network
  • data transfer between different simulation models
  • administration of the emission inventory.

The model system itself consists of two main parts:

  • The grid-based (Eulerian) mesoscale wind field model MEMO [Moussiopoulos 1989] calculates the wind field and all important meteorological parameters (e.g. temperatures, stability, pressure, solar radiation). Its output form one necessary input data set for the computation of the dispersion of air pollutants. Due to the use of advanced parameterisations and approximations of the governing equations MEMO is able to describe wind fields in rather complex terrain, e.g. sea breeze effects.

  • The dispersion model DYMOS. Using the same grid DYMOS calculates advective transport, turbulent diffusion and surface uptake of air pollutants. The incorporation of a gas phase chemical module enables the model to simulate reactions leading to photochemical smog situations with high near surface ozone levels.

    Figure 1 shows the main co mponents of the ECOSIM simulation system and the data flow between the different components.

Schematic overview of the ECOSIM system

Despite strong efforts in the design, the choice of the numerical methods and the implementation and optimization of the software system the use of the MEMO/DYMOS is still not suitable for forecast purposes because of its tremendous need of computation time. Even with the today`s fastest available computers it requires at least several hours computation time in order to simulate a 24-hours-period.

Although this shortcoming seems to be temporary because of the rapid growing computer power there is a urgent need for forecast purposes today. In order to give environmental agencies the possibility to run a photochemical dispersion model in a forecast mode two additional models have been included in the ECOSIM modelling system. These models are slim versions of the advanced MEMO/DYMOS system and are able to forecast ozone under restricted circumstances in an appropriate time, usually within two hours for a 24-hours-period. Nevertheless, the MEMO/DYMOS system is in principle able to operate in a forecast mode as soon as the necessary computer power is available. The ECOSIM system deals with urban scale air pollution problems and phenomena.

3. The ECOSIM Main Server

The ECOSIM server will be based on a UNIX workstation with a Xlib (X Windows) graphical user interface. Its basic user interface, display and analysis functionality will be based on the ACA ToolKit. The ACA ToolKit is a set of functions, utility programs, and objects, designed to support the efficient building of interactive information and decision support systems in the domain of environmental planning and management. Please note that the ToolKit is ESS proprietary software and is introduced into ECOSIM as Background. A public-domain alternative, based on HTML, and the implications of its use are described in a separate document. The ECOSIM main server provides the main user interface, coordination of functional components, and DSS part of the system. An example of the user interface and the model/GIS integration is provided by the AirWare air quality assessment and modelling system. Server functionality includes:

  • the top level start-up screen, showing systems identification and release information including the main icon menu, which offers an entry to the hypertext helpand explain system, which provides a limited functionality HTML browser, based on the CERN libraries used for the Mosaic browser. The browser is fully integrated into the ECOSIM server, and uses standard URL referencing of hypertext content files or cgi scripts.

  • the geographical information system (GIS) (please note that there are several other links and entry points to the GIS, in particular for/from the display of spatial model results or the database interface and (spatial) data analysis functions). The GIS is primarily a repository of all spatially referenced information, i.e., maps, together with interactive display and basic analysis functionality.

  • the basic objects data base that represents the main elements of the urban environmental systems and their attributes and functions, and its embedded rule-based expert system.

  • any case specific functionality, including a hierarchical decomposition of the spatial domain, and the various models or other specific information processing tools.

The main server is implemented on a (SUN) UNIX machine, and provides an X Windows (Xlib) based graphical user interface. This can also support PC clients in a local network running an X Windows emulation (under DOS, MS Windows, Windows NT, or an Apple/Macintosh), for example, using Hummingbird`s eXceed software. Based on its client-server architecture the server communicates with the various distributed information resources either through its internal functions and local file system, or through the http (hypertext transfer) protocol to access distributed resources like data bases, monitoring systems, or simulation models.

3.1 Geographical Information System (GIS)

Based on the ACA ToolKit, GIS functionality is integrated with data bases, models, expert systems, and DSS tools, as well as http servers supporting remote access through the Internet/WWW or appropriate client/server software. The GIS tools integrate vector and raster data as well as digital elevation data (DEMs), arbitrary zooming, multiple resolutions, tiling, hierarchical map layers, arbitrary overlay stacking and zooming, multiple windows, interactive colour editing and 3D display with arbitrary rotation and user-defined Gouraud shading.

Displaying output from spatial simulation and optimisation models, including the spatial interpretation of data as topical maps, is an important function of the GIS tools, furthermore the animation of interactive, dynamic models and statistical analysis in combination with any other map overlay. Spatial object data bases are accessible from the GIS through the display of corresponding symbols. These spatially referenced objects can represent, for example observation time series, compound objects (e.g., emission sources, industries, settlements, water bodies, parks and nature reserves, etc.). The main visualization features are:

  • display of overlays (points, lines, polygons, cell grids, rasters, DEM, spatial model results) from a (hierarchical) list of available overlays;
  • toggling on/off of features within an overlay;
  • restacking of overlays displayed;
  • display of satellite imagery or scanned maps as a background to the map overlays;
  • arbitrary zooming (arbitrary selection of location and size of the zooming window);
  • multiple resolutions with transparent switching between resolution within a map set;
  • support for tiling of large maps for increased performance;
  • highlighting and colour editing of map (legend) features;
  • switching from one to a four window mode, with four synchronized map windows;
  • query and read-back function for map attributes;
  • cross hair cursor synchronized in the four parallel windows;
  • 3D display of digital terrain models (or other numerical data such as model results or spatially interpolated measurement data) in arbitrary zooming states:
    • 3D display of elevation data over a 2D map background with Gouraud shading;
    • arbitrary 3D rotation of the display window (roll, pitch, yaw);
    • vertical stretching and positioning of the 3D display;
    • interactive selection of 3D overlays;
    • satellite imagery display in 3D (draping);
    • control of lighting source location and intensity for the shading;
  • animation of time series of maps (cell grids, lines/polygons, attribute time series);
  • extensive data import and export facilities;
  • dynamic linkage to the (optionally) embedded public-domain GRASS GIS for special data capture and analysis functions such as digitizing, map conversion and geo-correction, etc. - hardcopy output (HP-GL, PostScript);
  • support of WWW clients and Internet access (Netscape, Java);
  • hypertext (multimedia) on-line user manual, help functions, and metadata display;
  • access to spatially referenced object data such as observation time series, compound, objects (e.g., emission sources, industries, water bodies, settlement, parks and nature reserves).

Objects in turn provide linkages to (optional) environmental simulation models and a rule-based expert system for environmental assessment and data interpretation. Map sets in the GIS can be hierarchically structured, i.e., a regional GIS data set can contain several local, more detailed case study sub-regions, each with their individual map sets.

3.2 AirWare

AirWare is an interactive, integrated, model-based information and decision support system for air quality assessment and management. Designed for major cities and industrial areas, the flexibility of the modular software components makes it easy to adapt and apply the system to a broad range of physiographic and meteorological, socio-economic, and environmental conditions. The system is designed to provide easy access to advanced tools of data analysis and the design and evaluation of air quality control strategies, using a suite of simulation and optimization models together with integrated data base management, a geographical information system, and expert systems functionality. These components are integrated with an interactive and graphical user interface designed for users with little or no computer experience. The functionality of the AirWare system includes:

  • a GIS and object DB framework for emission sources and meteorology/monitoring data (with time series data display and analysis, including real-time coupling to monitoring equipment;

  • a hypertext system, with on-line help and explain functions, but also references to regulations, standard definitions and guidelines;

  • an embedded expert system (coupled to source inventory, used for the estimation of source terms);

  • integrated models:
    • US EPA ISC long- and short-term model;
    • US EPA PBM ozone box model;
    • an optimization routine for emission reduction strategies.

3.3 System Architecture

The ECOSIM demonstrator is based on a client-server architecture, taking advantage of the http (hypertext transfer) protocol. The main server provides the basic user interface and controls the user dialogue, displays information, and connects to external information resources (monitoring data, data bases, simulation models) as required. This communication is based on the public http protocol, and can be based on the Internet or dedicated connections (such as ISDN phone lines) for the physical communication layer. This protocol also forms the basis of World Wide Web browsers like Mosaic or Netscape.

3.4 Hardware Requirements and Operating System

The SERVER requires a UNIX workstation with a high-resolution graphics system (1280 by 1024), 256 simultaneous colours, i.e., an 8 bit graphics frame buffer) with 32 MB RAM and, depending on the application, 1 GB of disk space or more. The system will be supported on SUN Sparc architectures, under SUN Solaris 2.4 or higher. Support for alternative architectures and operating systems (HP UX, IBM AIX, Intel Linux) can be arranged, but would require additional resources not foreseen in the project.

In addition to the operating system, X Windows run-time libraries (usually part of the basic systems software) are required. Under SUN Solaris, the necessary libraries and tools are available under /usr/openwin/lib and /usr/openwin/bin. Window managers currently supported are the publish domain window managers uwm and ctwm as well as SUN OpenWindows. The system will also require a http daemon (httpd) for the client-server integration and Internet access. For more information on http and setting up a local Web server, see, for example [Liu 1994].

3.5 Networking

The host of the DEMONSTRATOR will require Internet access for some of its intended functionality. To access remote computing resources, a minimum of 64Kb will be required (for example through an ISDN phone line).

3.6 Auxiliary Tools

As external tools for data preparation, editing etc, the following tools are recommended:

  • basic screen editor like vi;
  • graphical editors and image conversion programs, including xv, xpaint, giftrans;
  • a standard GIS for geographical data capture and editing, for example:
    • Arc/Info (powerful but expensive);
    • IDRISI (PC based, reasonably cheap),
    • GRASS (public domain).

3.7 Data Requirements

The following data should be provided:

  • GIS data including spatially distributed model input;
  • object data including observation time series;
  • descriptor definitions and rules for the expert system;
  • hypertext content including metadata and user manuals;
  • auxiliary data including model interface descriptions.

3.7.1 GIS Data

The case study is defined for a geographical domain that includes the entire area of interest, i.e., the city boundaries and the surrounding area that is considered explicitly. For reasons of screen layout, this domain should either be square, or have a landscape (aspect ratio: 4x3) or portrait (aspect ratio: 3x4) layout. Other aspect ratios are in principle possible, but require extensive editing of the layout configuration files and may not lead to aesthetically pleasing" results. Please note that the screen layout is based on a 1280*1024 pixel resolution (also required for any PC client). The domain should be defined in terms of longitude/latitude (upper left and lower right corner) or the coordinate system used for the background maps. The GIS can handle any number of sub-areas within this domain (so-called case-study areas and sites), but cannot go beyond the domain boundaries (without re-implementing the GIS data base). For the domain and any sub-areas within the domain, for example, the region around a city and the city proper, the following map data sets (overlays) are recommended:

  • Basic topographical and land-use map (usually 1:50000 or 1:25000.

  • For the city proper, again depending on size, maps of up to 1:5000 could also be used to support zooming into small neighbourhoods, e.g., around a major source of pollution, individual street canyon, waste dump, etc.) in vector format (e.g., Arc/Info export format or any GIS format that can be converted with Arc/Info or GRASS) depending on the size of the area. Data should include at least features like administrative boundaries, transportation lines (rail, roads), water bodies, built-up areas.

  • Scanned maps (with names) for better orientation can be used and

  • satellite imagery, where available, LANDSAT TM and/or SPOT, for a recent land-use picture, and digital elevation model (gridded, typically 50x50 m spatial resolution and a minimum of 1 m vertical resolution). Where coastal waters are included this should also include the corresponding bathymetry.

  • Additional topical maps, e.g., geology, soils, vegetation, population density (in particular for the city area proper) can enhance the visualization facilities.

  • For the coupling of traffic emissions, a transportation graph (see also road segments in the objects data) is also required.

Please note that at least ONE map coverage is required to run the GIS system at all, but this could, in principle, be as simple as a single line (vector) overlay with the city boundaries. All other overlays (including scanned maps, satellite imagery, and DEM) are optional, but of course will limit the functionality if absent. Furthermore, several of the models require (gridded) input data sets that can also be included in form of a map.

3.7.2 Object Data

The main data base (top layer) of the ECOSIM server is using an object structure. Objects describe important functional elements (see the class examples below), and the different (model) scenarios that represent alternative management options. The following object classes are available or are under construction in parallel projects:

  • Climate_stations CS (including temperature, wind etc.)
  • Air_Quality AQ measurement station
  • Water_quality WQ measurement station
  • Flow_measurements (rivers) FS measurement station
  • Settlements SE hierarchy of administrative units
  • City district CD
  • City blocks CB
  • Point_source PTS (atmospheric, generic)
  • Area_source ARS (atmospheric, generic)
  • Line_source LNS (atmospheric, generic)
  • Thermal_Power_Station TPS
  • Road segments RDS
  • Road intersections RDI
  • Industries (generic) IND
  • Waste_water_treatment_plants WWTP
  • Process_Plants PPLANT
  • Chemicals_Storage CST
  • Waste_Dumps WDUMP
  • Waste_Processing WPROC -
  • Waste_Incinerators WINC
  • Port_Facilities PORTF
  • Airports AIRP
  • Marshalling_Yards MYARD
  • Railway_Loading_Docks RLDOCK
  • Animal_farms_and_feedlots AFF
  • Irrigation_districts IRR
  • Sub-catchments SC
  • Dams_and_Reservoirs DAM
  • Natural_Lakes LAKE
  • Abstractions_(river) ABS
  • Aquifers AQUF
  • Wells WELL
  • Scenic_sites SS

Each object class is defined by a header, specifying

  • NA : a name of the object;

  • ID: the unique ID which by convention consists of the class label (see above) and an integer number;

  • location, defined by LO, LA, EL longituted, latitude, and elevation (meters above sea level) or x,y, and z coordinates in the respective domain coordinate system;

  • HY: hypertext link, pointing to an explanation file with meta data and background information;

  • PI: hypertext link, pointing to the pictorial and textual (multi-media) description of the object;

  • geographical reference, a table of symbolic labels describing (hierarchical) location references like city, district, block, with the associated overlay element.

This provides the linkage to the GIS layer. Depending on the type of objects, this is followed by a list (set) of parameters (DE Descriptors), known to the expert system and defined in its knowledge base in terms of name, unit, legal range, and method of retrieval/estimation/deduction). Special object-specific data (like tables, graphs, time series etc. and the methods for their retrieval), are defined in a set of table display functions, and links to other objects. As a special feature, in particular for all observation station objects and objects that include time series data, specific time series display and analysis tools, including a number of statistical tests and spatial interpolation, are available as methods to these objects through their display functions.

Please note that the retrieval methods for any of the object`s data can include access to the local file system, (networked) data bases (embedded SQL), or remote (Internet) data sources (URLs). Once the (sub)set of objects required for each case study is selected, the detailed parameter lists for each object class can be defined.

4 The Models - Technical Description

4 The Models - technical description
4.1 The Mesoscale meteorological model MEMO

The original version of the nonhydrostatic mesoscale model MEMO was developed at the Universität Karlsruhe. In the last years MEMO has been increasingly installed and utilized at several research institutions throughout Europe. Recently, MEMO was selected as one of the core models of the EURAD Zooming Model (EZM) to be used for the refined modelling of transport and chemical transformation of pollutants in selected European regions in the frame of the EUROTRAC project [EUROTRAC 1992]. Further development of MEMO is currently undertaken at both the Universität Karlsruhe and the Aristotle University of Thessaloniki. In the following, a complete description of MEMO is given in brief outlines. More details can be found in [Moussiopoulos 1989], [Flassak 1990], [Moussiopoulos 1993].

4.1.1 Model equations

The prognostic mesoscale model MEMO describes the dynamics of the atmospheric boundary layer. In the present model version air is assumed to be unsaturated. The model solves the continuity equation, the momentum equations and several transport equations for scalars (including the thermal energy equation and, as options, transport equations for water vapour, the turbulent kinetic energy and pollutant concentrations).

The conservation equations solved in MEMO take the form:

x', y' and z' represent Cartesian coordinates and u, v and w are wind velocity components in the x', y' and z' direction, respectively. y is any scalar (e.g. the potential temperature q or the turbulent kinetic energy E).

In the conservation equations Ru , Rv , Rw and Ry represent turbulent diffusion which is discussed later, while Cu , Cv and Cw represent volumetric Coriolis force components ( C = 2rW¥v ). The source/sink terms Qy depend on the transported scalar quantity: For the potential temperature this term includes anthropogenic heat emission and the divergence of radiative fluxes, for the turbulent kinetic energy it contains the shear and buoyancy production rates as well as the dissipation rate.

Following the common practice in mesoscale models, variables are split into base-state parts (denoted by overbars) and mesoscale perturbations (marked by a hat). By definition, the base-state parts of the wind velocity components are taken as zero. For the thermodynamic variables the separation yields

In order to accelerate the convergence of the solution of the elliptic equation for pressure (see below), the mesoscale pressure perturbation is split into three components:

The first term on the right-hand-side of this equation represents the large-scale horizontal pressure gradient in the considered mesoscale domain. The hydrostatic part, i.e. the second term on the right-hand-side of Eq. (9), is obtained by integrating the hydrostatic equation:

where r follows from the ideal gas law. It should be noted that the vertical derivative of ph and the buoyancy term in Eq. (3) cancel each other out.

The last term on the right-hand-side of Eq. (9), the nonhydrostatic part, is computed implicitly by solving the elliptic pressure equation.

4.1.2 Transformation to terrain-influenced coordinates

The lower boundary of the model domain coincides with the ground. Because of the inhomogeneity of the terrain, it is not possible to impose boundary conditions at that boundary with respect to Cartesian coordinates. Therefore, a transformation of the vertical coordinate to a terrain-following one is performed:

H and h(x',y') are the (temporally and spatially constant) height of the model top boundary and the altitude at the location (x',y'), respectively. To allow for nonequidistant meshsizes, for example to achieve a better resolution near the ground, the additional transformation x = x(x') , y = y(y') , z = z(h) is employed, where x(x'), y(y') and z(h) represent arbitrary monotonic functions [Schumann 1984]. Hence, the originally irregularly bounded physical domain is mapped onto one consisting of unit cubes.

4.1.3 Numerical solution of the equation system

The discretized equations are solved numerically on a staggered grid, i.e. the scalar quantities r, p and q are defined at the cell centre while the velocity components u, v and w are defined at the centre of the appropriate interfaces.

Temporal discretization of the prognostic equations is based on the explicit second order Adams-Bashforth scheme. There are two deviations from the Adams-Bashforth scheme: The first refers to the implicit treatment of the nonhydrostatic part of the mesoscale pressure perturbation pnh. To ensure non-divergence of the flow field an elliptic equation is solved. The elliptic equation is derived from the continuity equation wherein velocity components are expressed in terms of pnh. It should be noted that since the elliptic equation is derived from the discrete form of the continuity equation and the discrete form of the pressure gradient conservativity is guaranteed [Flassak 1988]. The discrete pressure equation is solved numerically with a fast elliptic solver in conjunction with a generalized conjugate gradient method. The fast elliptic solver is based on fast Fourier analysis in both horizontal directions and Gaussian elimination in the vertical direction [Moussiopoulos 1989].

The second deviation from the explicit treatment is related to the turbulent diffusion in vertical direction. In case of an explicit treatment of this term, the stability requirement may necessitate an unacceptable abridgement of the time increment. To avoid this, vertical turbulent diffusion is treated using the second order Crank-Nicolson method.

On principle, advective terms can be computed using any suitable advection scheme. In the present version of MEMO a 3-D second-order total-variation-diminishing (TVD) scheme is implemented which is based on the 1-D scheme proposed by Harten [Harten 1986]. It achieves a fair (but not any) reduction of numerical diffusion, the solution being independent of the magnitude of the scalar (i.e. preserving transportivity).

4.1.4 Parameterizations

Turbulence and radiative transfer are the most important physical processes that have to be parameterized in a prognostic mesoscale model. In the MEMO model radiative transfer is calculated with an efficient scheme based on the emissivity method for longwave radiation and an implicit multilayer method for shortwave radiation [Moussiopoulos 1987].

The diffusion terms in Eqs (1-3) and (5) may be represented as the divergence of the corresponding fluxes. For turbulence parameterizations K-theory is applied. In case of MEMO turbulence can be treated either with a zero-, a one- or a two-equation turbulence model.

For most applications a one-equation model is used, where a conservation equation for the turbulent kinetic energy E is solved. The eddy viscosity Km follows from

where a length scale (the so-called mixing length) is needed. An appropriate length scale is given by Blackadar [Blackadar 1962]:

(k : von Kármán constant). In MEMO the limiting value l· is given by .

The eddy diffusivity KY can be calculated from the turbulent Prandtl number Prt Km / KY. The turbulent Prandtl number is given by

(typically zD = 1000 m).

4.1.5 Initial and boundary conditions

In the prognostic model MEMO initialization is performed with suitable diagnostic methods: A mass-consistent initial wind field is formulated using an objective analysis model; scalar fields are initialized using appropriate interpolating techniques [Kunz 1991]. Data needed to apply the diagnostic methods may be derived either from observations or from larger scale simulations.

Suitable boundary conditions have to be imposed for the wind velocity components u , v, w, the potential temperature q and pressure p at all boundaries. At open boundaries, wave reflection and deformation may be minimized by the use of so-called 'radiation conditions' [Orlanski 1976]. As the original formulation of the radiation conditions merely allows disturbances to propagate out through the boundary, but does not allow information from outside to be imposed at the boundary expanded radiation conditions at lateral boundaries

are used [Carpenter 1982], where n is the direction perpendicular to the boundary and C the phase velocity that includes wave propagation and advection. C is calculated for every variable f from known values from the interior of the computational domain. Eq. (15) differs from the standard formulation of the radiation condition [Orlanski 1976] by the last two terms describing the temporal and spatial change of the undisturbed large scale environment. This information can be specified to the model using either measurements or results of larger scale simulations.

According to the experience gained so far with the model MEMO, neglecting large scale environmental information might result in instabilities in case of simulations over longer time periods.

For the nonhydrostatic part of the mesoscale pressure perturbation, homogeneous Neumann boundary conditions are used at lateral boundaries. With these conditions the wind velocity component perpendicular to the boundary remains unaffected by the pressure change.

At the upper boundary Neumann boundary conditions are imposed for the horizontal velocity components and the potential temperature. To ensure non-reflectivity, a radiative condition is used for the hydrostatic part of the mesoscale pressure perturbation ph at that boundary. Hence, vertically propagating internal gravity waves are allowed to leave the computational domain [Klemp 1983]. For the nonhydrostatic part of the mesoscale pressure perturbation, homogeneous staggered Dirichlet conditions are imposed. Being justified by the fact that nonhydrostatic effects are negligible at large heights, this condition is necessary, if singularity of the elliptic pressure equation is to be avoided in view of the Neumann boundary conditions at all other boundaries.

The lower boundary coincides with the ground (or, more precisely, a height above ground corresponding to its aerodynamic roughness). For the nonhydrostatic part of the mesoscale pressure perturbation, inhomogeneous Neumann conditions are imposed at that boundary. All other conditions at the lower boundary follow from the assumption that the Monin-Obukhov similarity theory is valid. With the exception of water surfaces (where the temperature is specified) the surface temperature is calculated from the nonlinear heat balance equation

( R : longwave , S : shortwave radiative fluxes (Ro ~To4); Qs : heat flux to the soil ; Qo , Lo : sensible and latent heat to the atmosphere ; Qa : anthropogenic heat flux), which is solved using a Newton iteration technique. Radiative fluxes follow from the above mentioned radiation scheme [Moussiopoulos 1987]. For the calculation of the soil temperature, an one-dimensional heat conduction equation for the soil is solved. The specific humidity at the surface (needed even in the case that the transport equation for water vapor is not included in the model) is computed from

where Y is the evaporation parameter varying between 0 and 1 (water: Y = 1 , dry soil: Y = 0), qs is the saturation humidity and q2 is the humidity computed at the lowermost grid level above ground.

4.1.6 Grid definition

The governing equations are solved numerically on a staggered grid. Scalar quantities as the temperature, pressure, density and also the cell volume V are defined at the center of a grid cell and the velocity components u,v,w at the center of the appropriate interface. Turbulent fluxes are defined at different locations: Shear fluxes are defined at the center of the appropriate edges of a grid cell and normal stress fluxes at scalar points. With this definitions the outgoing fluxes of momentum, mass, heat and also turbulent fluxes of a grid cell are identical to incoming flux of the adjacent grid cell. So the numerical method is conservative.

4.1.7 Nesting facility

In MEMO Vs. 5.0 a one-way interactive nesting scheme is implemented. With this nesting scheme a coarse grid (CG) and a fine grid (FG) simulation can be nested. During the CG simulation data is interpolated and written to a file. A consecutive FG simulation uses this data (FT-UNIT 42) as lateral boundary values.

Compiling: When building up the FORTRAN-library for the CG simulation an INCLUDE-file has to be provided containing the maximum number of grid cells of the FG domain (IMN, JMN and KMN).

Controlling of CG run: Control of the nesting facility is performed from the FT-UNIT 5 file. If the current simulation is used as a CG simulation INESO must be set to 1. Furthermore the origin of the FG domain within the CG domain (XZERO and YZERO in [m]) must be given. The nesting factor, i. e. the relation (XNEST=DXCG/DXFG) of the grid spacing of both grids must be specified. The nesting factor XNEST should not exceed the value of 3.

The update time interval is given in the CG simulation by DTNESS. It should be chosen as small as possible (best results with DTNESS=DTCG) but should not exceed 3*DTCG.

Controlling of FG run: Within the control file of the FG simulation INESI can be either set to 1 or to 2. INESI=1 uses expanded radiation conditions for the nesting facility while INESI=2 applies a relaxation scheme. The latter being more efficient at inflow boundaries whereas expanded radiation conditions seem to be superior at outflow boundaries.

The FG simulation itself can also be used as CG simulation for a finer grid and so forth.

IMPORTANT: The altitude of the model top (ZT) and the minimal vertical grid spacing (HHMIN) must be identical in both simulations.

4.2 The Dispersion and chemical model DYMOS

The DYMOS submodel is a three-dimensional grid model designed to compute the concentration of air pollutants in the troposphere. It includes transport processes of advection and turbulent diffusion as well as chemical changes and surface uptake processes. The mathematical formulation is given by the advection-diffusion equation, which represents the mass balance for a species:

ci concentration of species i

u, v, w wind speed components

Kh horizontal diffusion coefficient

Kv vertical diffusion coefficient

Qi emission rate of the species i

Ri source or sink term due to chemical reactions

Di surface dry deposition of species i

The solution of this equation has to be computed for every substance. The model employs numerical differencing techniques on the discretized grid given by the MEMO output. An efficient solution technique is the method of fractional steps so that the terms representing different atmospheric processes are solved separately using the most suitable scheme.

The variables u, v, w, Kh and Kv are computed by the mesoscale meteorological model MEMO and read in given time intervals. The anthropogenic partition of emission rate Qi of the substance i could be obtained in different ways:

    • determined by the emission data base for the area, modified by a time dependent emission factor describing seasonal, weekly or diurnal variation of the emission strength,
    • model computation by an emission model (e.g. traffic flow model),
    • on-line emission measurements.

Only a fraction of all substances are primary species. They are emitted by anthropogenic activities or by the ecosystem (biogenic emissions). Others may appear as a result of the chemical reaction scheme. The deposition rate Di of a substance is approached by a resistance model (see below).

4.2.1 The Chemical mechanism

High surface near ozone concentrations are among the harzardous environmental problems. Ozone is formed in the troposphere through chemical reactions between nitrogen oxides (NOx) and volatile organic compounds (VOC) under presence of solar radiation. Hundreds of substances and thousands of reactions participate in the ozone formation process. That is why the explicit treatment of the chemistry is impossible in a Eulerian dispersion model for the following two reasons:

    • limitation on computer resources,
    • lack of input data.

Thus the photochemical mechanism must be compressed significantly. The chosen chemical reaction scheme is the Carbon Bond Mechanism (CBM-IV). Reactions with minor importance in the given time scale and typical concentration range were removed from the explicit scheme. The organic compounds are disaggregated based on the carbon bonds of the organic compounds. For example, butene would be split into one olefinic bond (OLE) and two paraffinic bonds (PAR). Some organic species remain explicit because of their special role in the ozone chemistry. CBM IV contains the following primary classes of nonmethane VOC`s:

    • PAR 1 C-atom carbon single bond s
    • ETH 2 C-atoms ethene
    • OLE 2 C-atoms carbon double bonds, olefine
    • TOL 7 C-atoms toluole, monoalkylbenzene
    • XYL 8 C-atoms xylole, C8 aromatic species
    • FORM 1 C-atom formaldehyde
    • ALD2 2 C-atoms acetaldehyde, higher aldehydes
    • ISOP 5 C-atoms isoprene
    • NR 1 C-atom nonreactive VOC`s

VOC`s are any compounds of carbon excluding carbon monoxide, carbon dioxide, carbon acid, metallic carbides or carbonates and ammonium carbonates. Methane is also excluded from this VOC definition because of its low importance on the ozone production in the urban scale.

The complete reaction scheme of the lumped chemical mechanism CBM VI is given in Table 1. The reaction rates could be constant, temperature dependent or radiation dependent (photolysis rates). These variable photolysis rates depend on light intensity and spectral distribution. The parameters are computed internally as a function of longitude and latitude, day of the year, time and fraction of cloud cover.

The mathematical formulation for the chemical reaction scheme leads to a stiff system of ordinary differential equations (ODES). The equations contain a wide variation in the reaction rate constants. The ODES must be solved for every grid box. No matter what numerical algorithm is selected, the solution of such systems remains computationally expensive. On the other hand, a parallel implementation of this code guarantees a relatively high acceleration (see installation) and allows a much faster execution in comparison with most serial high performance computer systems.

In addition, it is not always necessary to compute chemical changes at every dispersion time step. This speeds up the whole procedure tremendously because of the self-adjusting character of the numerical algorithm.

The chemical kinetic mechanism CBM IV



Rate Constant at 298 K

Activation Energy



NO2 = NO + O



O = O3




O3 + NO = NO2




O + NO2 = NO




O + NO2 = NO3




O + NO = NO2




NO2 + O3 = NO3




O3 = O



O3 = O1D



O1D = O




O1D + H2O = 2 OH




O3 + OH = HO2




O3 + HO2 = OH




NO3 = 0.89 NO2 + 0.89 O

+ 0.11 NO



NO3 + NO = 2 NO2




NO3 + NO2 = NO + NO2




NO3 + NO2 = N2O5




N2O5 + H2O = 2 HNO3




N2O5 = NO3 + NO2




NO + NO = 2 NO2




NO + NO2 + H2O = 2 HNO2




NO + OH = HNO2




HNO2 = NO + OH



OH + HNO2 = NO2




HNO2 + HNO2 = NO + NO2




NO2 + OH = HNO3




OH + HNO3 = NO3




HO2 + NO = OH + NO2




HO2 + NO2 = PNA




PNA = HO2 + NO2




OH + PNA = NO2




HO2 + HO2 = H2O2




HO2 + HO2 + H2O = H2O2




H2O2 = 2 OH



OH + H2O2 = HO2




OH + CO = HO2




FORM + OH = HO2 + CO




FORM = 2 HO2 + CO






FORM + O = OH + HO2

+ CO




FORM + NO3 = HNO3 + HO2

+ CO




ALD2 + O = C2O3 + OH




ALD2 + OH = C2O3




ALD2 + NO3 = C2O3 + HNO3




ALD2 = FORM + 2 HO2

+ CO + XO2



C2O3 + NO = FORM + NO2

+ HO2 + XO2




C2O3 + NO2 = PAN




PAN = C2O3 + NO2




C2O3 + C2O3 = 2 FORM + 2 XO2

+ 2 HO2




C2O3 + HO2 = 0.79 FORM + 0.79 XO2

+ 0.79 HO2 + 0.79 OH





+ HO2




PAR + OH = 0.87 XO2 + 0.13 XO 2 N

+ 0.11 HO2 + 0.11 ALD2

- 0.11 PAR + 0.76 ROR




ROR = 0.96 XO2 + 1.1 ALD2

+ 0.94 HO2 - 2.1 PAR

+ 0.04 XO2N + 0.02 ROR








ROR + NO2 = Produkte




O + OLE = 0.63 ALD2 + 0.38 HO2

+ 0.28 XO2 + 0.3 CO

+ 0.2 FORM + 0.02 XO2N

+ 0.22 PAR + 0.2 OH





- PAR + XO2

+ HO2




O3 + OLE = 0.5 ALD2 + 0.74 FORM

+ 0.22 XO2 + 0.1 OH

+ 0.33 CO + 0.44 HO2





NO3 + OLE = 0.91 XO2 + FORM

+ 0.09 XO2N + ALD2

+ NO2 - PAR




O + ETH = FORM + 1.7 HO2

+ CO + 0.7 XO2

+ 0.3 OH




OH + ETH = XO2 + 1.5 FORM

+ 0.22 ALD2 + HO2




O3 + ETH = FORM + 0.42 CO

+ 0.12 HO2




TOL + OH = 0.44 HO2 + 0.08 XO2

+ 0.36 CRES + 0.56 TO2




TO2 + NO = 0.9 NO2 + 0.9 HO2

+ 0.9 OPEN




TO2 = CRES + HO2




OH + CRES = 0.4 CRO + 0.6 XO2

+ 0.6 HO2 + 0.3 OPEN








CRO + NO2 = products




OPEN = C2O3 + HO2

+ CO



OPEN + OH = XO2 + 2 CO

+ 2 HO2 + C2O3 + FORM




OPEN + O3 = 0.03 ALD2 + 0.62 C2O3

+ 0.7 FORM + 0.03 XO2

+ 0.69 CO + 0.08 OH

+ 0.76 HO2 + 0.2 MGLY




OH + XYL = 0.7 HO2 + 0.5 XO2

+ 0.2 CRES + 0.8 MGLY

+ 1.1 PAR + 0.3 TO2




OH + MGLY = XO2 + C2O3




MGLY = C203 + HO2

+ CO



O + ISOP = 0.6 HO2 + 0.8 ALD2

+ 0.55 OLE + 0.5 XO2

+ 0.5 CO + 0.45 ETH

+ 0.9 PAR





+ 0.67 HO2 + 0.13 XO2N

+ ETH + 0.4 MGLY

+ 0.2 C203 + 0.2 ALD2




O3 + ISOP = FORM + 0.4 ALD2

+ 0.55 ETH + 0.2 MGLY

+ 0.1 PAR + 0.06 CO

+ 0.44 HO2 + 0.1 OH








XO2 + NO = NO2




XO2N + NO = products




XO2 + XO2 = products




SO2 + OH = HO2 + SULF











The temperature dependent reaction rate is given by

where Ea is the activation energy, R is the gas constant, T is the temperature in degrees Kelvin and kr is the rate at this temperature given in Table 1.

4.2.2 Treatment of surface removal processes

Dry deposition is an important removal process of a lot of substances of interest for air pollution computation. The computation of the surface uptake processes follows in principle the work of Wesely [Wesely 1988].

4.3 The Forecast system
4.3.1 The Ozone model REGOZON

The principle design of the model REGOZON [Mieth 1993] is similar to the MEMO/DYMOS system. A wind field model is coupled with a dispersion model and a chemical module. For the computation of meteorological values and the dispersion, an enhanced version of the Eulerian grid model REWIMET [Heimann 1985] is used. It is based on the conservation laws for impulse, mass, energy and passive constituents comparable to MEMO/DYMOS but with a higher level of approximations. A hydrostatically stratified atmosphere is assumed, which is dry and incompressible. The model equations are expressed in three vertical layers. The first (surface layer) follows the ground level and has a fixed vertical thickness of 50 m above ground. It is turbulent mixed and its physical behaviour is strongly coupled with the surface characteristics. Emissions from traffic, from households and from industrial sources with low emission heights are introduced into the surface layer. The second layer (mixed layer) reaches from the upper level of the surface layer to the upper level of the atmospheric boundary layer, up to the mixing height. This layer is also turbulently mixed and shows the characteristic diurnal variation of the thickness of the atmospheric boundary layer. Emissions from higher emission sources, for example high stacks from power stations, enter the mixed layer. The third layer (temporary layer) is located above the mixed layer. It is assumed to be free of turbulence. Since the atmospheric boundary layer can expand to the suprascale inversion, it is possible for the temporary layer to disappear. It will be recreated when the atmospheric boundary layer sinks. No substances are emitted in this layer but it transports the suprascale background concentrations of ozone and ozone precursor substances above the atmospheric boundary layer. The vertical model structure is shown in Figure 3. For the computation of the temperature regime a surface energy budget routine has been added.

The computation of the dispersion is carried out immediately after the determination of the meteorological values. The transport model uses the same vertical structure. The dispersion equation is solved in two dimensions but with an allowed vertical exchange according to the determined stability. The chemical changes form a source or sink term in the dispersion equation. The photochemical scheme of CBM IV [Gery 1988] is applied; the chemical module is nearly identical with the DYMOS module. As characteristic time steps for the solution of the chemical system are much smaller, the chemical computations have been decoupled from the determination of transport. Usually, chemistry is not computed at every transport time step [Mieth 1996].

Vertical structure of REGOZON

Dry deposition velocities were computed using a resistance approach as a function of land utilization as described in [Wesely 1988] , the computation of biogenic emissions uses the methods from [Pierce 1990].

The forecast system REGOZON is applicable under the following restrictions:

    • relatively flat terrain without strong orographic features and land/see breeze effects
    • fair weather conditions with moderate wind speeds
    • stable, barotropic synoptic situations.

Because of the strong constraints and the rough vertical model resolution, the computational time for the REGOZON model is small. A comparison between this relatively simple model and a nonhydrostatic 3d model with 35 vertical layers [Kapitza 1992] show a similar ozone production for the selected episodes which conform with the above restrictions. The extent of computation time requiered by the advanced MEMO/DYMOS system is about an order of magnitude higher.

4.3.2 The Multilayer dispersion model MUSE

Whereas REGOZON is a combined mesoscale meteorological and dispersion model, MUSE [Sahm 1995] is a dispersion model only, requiring the meteorological quantities usually computed by MEMO. On one hand it enhances the computation time in comparison to REGOZON, on the other hand it could be applied to a lot more situations, because of the lower constrains of the meteorological parts. For example, whereas REGOZON seems to be unable to describe the transport processes related to sea breeze, MUSE is able to perform this task.

The model equations solved by MUSE are those of the fully 3D version of the dispersion model DYMOS. However, the multilayer model MUSE differs from the ECOSIM main dispersion model in the following points:

    • Instead of discretizing the model domain into 10-30 non-equidistantly distributed vertical layers, three layers are being used in the vertical direction. This modelling concept is very similar to the REGOZON model structure. The depth of the lower layer (which practically corresponds to the surface layer) is prescribed and kept constant in time. The depth of the middle layer is permitted to vary following the diurnal variation of the mixing height. The latter is described by Deardorff's prognostic equation [Deardorff 1974] which was shown to lead to realistic results when the variation of the mixing height with time is strongly influenced by surface heating [Pielke 1984]. For periods of stability, an algebraic parameterization based on the friction velocity and the Monin-Obukhov length is applied [Moussiopoulos 1991]. The limit of the upper layer coincides with the fixed domain top. The assumption of a time dependent depth of the middle and upper layer leads to the need of an entrainment term [Ngyen 1989].
    • The description of vertical transport due to turbulent diffusion at the lower limit of the middle layer is based on the turbulent kinetic energy defined at the interface of the lower and middle layers. To avoid unrealistic vertical diffusion rates associated with the growth of the mixing height, one-sided concentration gradients are calculated at this interface. Turbulent diffusion in the upper limit of the middle layer is neglected.
    • Advective transport is described using the scheme designed by Smolarkiewicz [Smolarkiewicz 1984]. The average wind speed in the middle and upper layers is calculated by integrating the fluxes of the corresponding layers of the fully 3D wind field model.
    • The differential equation system in MUSE is solved with a backward difference (BDF) solution procedure (by applying the Gauss-Seidel iteration scheme) [Kessler 1995]. Because of the nature of this semi-implicit algorithm, vertical diffusive transport and chemical transformation of pollutants have to be treated separately.
    • Aspects related to the formation of photochemical oxidants could be analyzed with the multilayer model MUSE in conjunction with the chemical reaction mechanism KOREM (a modified version of the Bottenheim-Strausz mechanism [Bottenheim 1982]). This is a more compressed chemical reaction scheme.

With the above simplifications, the MUSE code becomes approximately 5 times faster than the full 3D version. At the same time, the memory requirements are reduced by 90 %. The model configuration MEMO/MUSE can be used for forecast purposes in more complex terrain than REGOZON under defined conditions.

4.4 The Groundwater model

The arithmetical models that were used in the framework of the ECOSIM project, simulate the basic hydrological and transport processes in the studying area in a deterministic way. The simulation with arithmetical models provides the possibility of the estimation of the contaminants' concentrations at every point of the particular area and at any time in the past, the present and the future.

4.4.1 Simulation of the hydraulic conditions

For the estimation of the water movement module, the MODFLOW model was implemented. MODFLOW is a three-dimensional finite-difference ground-water flow model, developed by the US Geological Survey Department.

Model description. MODFLOW is a three-dimensional finite-difference ground-water flow model. It has a modular structure that allows it to be easily modified to adapt the code for a particular application. MODFLOW simulates steady and nonsteady flow in an irregularly shaped flow system in which aquifer layers can be confined, unconfined, or a combination of confined and unconfined. Flow from external stresses, such as flow to wells, areal recharge, evapotranspiration, flow to drains, and flow through river beds, can be simulated. Hydraulic conductivities or transmissivities for any layer may differ spatially and be anisotropic (restricted to having the principal direction aligned with the grid axes and the anisotropy ratio between horizontal coordinate directions is fixed in any one layer), and the storage coefficient may be heterogeneous. The model requires input of the ratio of vertical hydraulic conductivity to distance between vertically adjacent block centres. Specified head and specified flux boundaries can be simulated as can a head dependent flux across the model's outer boundary that allows water to be supplied to a boundary block in the modelled area at a rate proportional to the current head difference between a "source" of water outside the modelled area and the boundary block. MODFLOW is currently the most used numerical model in the US Geological Survey for ground-water flow problems.

MODFLOW is divided into a series of components called "packages." Each package performs a specific task. Some of the packages are always required for a simulation and some are optional. The input for each package is contained in a separate text file.

Basic Package (BAS1). The basic package is always required. The input to the basic package includes the grid dimensions, the computational time steps, and an array identifying which packages are to be used.

Output Control (OUT1). The output control file controls what information is to be output from MODFLOW and when it is to be output.

Block Centered Flow Package, Version 1 (BCF1). The block centered flow package performs the cell by cell flow calculations. The input to this package includes layer types and cell attributes such as storage coefficients and transmissivity.

River Package (RIV1). The river package is used to simulate river type boundary conditions.

Recharge Package (RCH1). The recharge package is used to simulate recharge to the groundwater from precipitation.

Well Package (WEL1). The well package is used to simulate well type boundary conditions.

Drain Package (DRN1). The drain package is used to simulate drain type boundary conditions.

Evapotranspiration Package (EVT1). The evapotranspiration package is used to simulate the effect of evapotranspiration in the vadose zone.

General Head Package (GHB1). The general head boundary package is used to simulate head type boundary conditions where the value of the boundary condition depends on the head in the cell. Commonly used to simulate lakes.

Strongly Implicit Procedure Package (SIP1). The SIP package is an iterative solver based on the strongly implicit procedure.

Pre-conditioned Conjugate Gradient Package, (PCG2). The PCG2 package is an iterative solver based on the preconditioned conjugate gradient technique.

Slice Successive Over-relaxation Package (SOR1). The SSOR package is an iterative solver based on the slice successive overrelaxation technique.

Method. The ground-water flow equation is solved using the finite-difference approximation. The flow region is considered to be subdivided into blocks in which the medium properties are assumed to be uniform. The plan view rectangular discretization results from a grid of mutually perpendicular lines that may be variably spaced. The vertical direction zones of varying thickness are transformed into a set of parallel "layers". Several solvers are provided for solving the associated matrix problem; the user can choose the best solver for the particular problem. Mass balances are computed for each time step and as a cumulative volume from each source and type of discharge.

In order to use MODFLOW, initial conditions, hydraulic properties, and stresses must be specified for every model cell in the finite-difference grid.

Output options. Primary output is head, which can be written to the listing file or into a separate file. Other output includes the complete listing of all input data, draw down, and budget data. Budget data are printed as a summary in the listing file, and detailed budget data for all model cells can be written into a separate file.

4.4.2 Simulation of the pollutants transport

Model oescription. MT3D model, is a model for simulation of advection, dispersion and chemical reactions of contaminants in groundwater flow systems in either two or three dimensions. The model uses a mixed Eulerian-Lagrangian approach to the solution of the advective-dispersive-reactive equation, based on combination of the method of characteristics and the modified method of characteristics. This approach combines the strength of the method of characteristics for eliminating numerical dispersion and the computational efficiency of the modified method of characteristics. The model program uses a modular structure similar to that implemented in the US Geologic Survey modular three-dimensional finite-difference groundwater flow model, referred to as MODFLOW, (McDonald and Harbaugh, 1988). The modular structure of the transport mode makes it possible to simulate advection, dispersion, source/sink mixing, or chemical reactions independently without reserving computer memory space for unused options; new packages involving other transport processes can be added to the model readily without having to the existing code.

The MT3D transport model was developed for use with any block-centered finite-difference flow model such as MODFLOW and is based on the assumption that changes in concentration field will not affect the flow field significantly. After a flow model is developed and calibrated, the information needed by the transport model can be saved in disk files which are then retrieved by the transport model.

The MT3D transport model can be used to simulate changes in concentration of single-species miscible contaminants in groundwater considering advection, dispersion and some simple chemical reactions, with various types of boundary conditions and external sources or sinks. The chemical reactions included in the model are equilibrium-controlled linear or non-linear sorption and first-order irreversible decay or biodegradation. Currently, MT3D accommodates the following spatial discretization capabilities and transport boundary conditions: (1) confined, unconfined or variably confined/unconfined aquifer layers; (2) inclined model layers and variable cell thickness within the same layer; (3) specified concentration or mass flux boundaries; and (4) the solute transport effects of external sources and sinks such as wells, drains, rivers, areal recharge and evapotranspiration. The primary modules of MT3D are:

    • Basic BTN handles basic tasks that are required by the entire transport model. Among these tasks are definition of the problem, specification of the boundary and initial conditions, determination of the step size, preparation of mass balance information, and printout of the simulation results.
    • Flow Model FMI Interfaces with a flow model. Currently, the interfacing is done through an unformatted disk file containing heads and flow terms. The FMI Package reads the contents of this file and prepares heads and flow terms in the form needed by the transport model.
    • Advection ADV Solves the concentration change due to advection with one of the three mixed Eulerian-Lagrangian schemes included in the package: MOC, MMOC, or HMOC.
    • Dispersion DSP Solves the concentration change due to dispersion with the explicit finite difference method.
    • Sink & Source SSM Solves the concentration change due to fluid Mixing sink/source mixing with the explicit finite difference method. Sink/source terms may include wells, drains, rivers, recharge and evapotranspiration. The constant-head boundary and general-head-dependent boundary are also handled as sink/source terms in the transport model.
    • Chemical RCT Solves the concentration change due to chemical reactions. Currently, the chemical reactions include linear or nonlinear sorption isotherms and first-order irreversible rate reactions (radioactive decay or biodegradation).
    • Utility UTL Contains a number of utility modules that are called upon by primary modules to perform such general-purposed tasks as input/output of data arrays.

Method. The advective-dispersive-reactive equation describes the transport of miscible contaminants in groundwater flow systems. Most numerical methods for solving the advective-dispersive-reactive equation can be classified as Eulerian, Lagrangian or mixed Eulerian-Lagrangian. In the Eulerian approach, the transport equation is solved with a fixed grid method such as the finite-difference or finite-element method. The Eulerian approach offers the advantage and convenience of a fixed grid, and handles dispersion/reaction dominated problems effectively. For advection-dominated problems which exist in many field conditions, however, an Eulerian method is susceptible to excessive numerical dispersion or oscillation, and limited by small grid spacing and time steps. In the Lagrangian approach, the transport equation is solved in either a deforming grid or deforming coordinate in a fixed grid. The Lagrangian approach provides an accurate and efficient solution to advection dominated problems with sharp concentration fronts. However, without a fixed grid or coordinate, a Lagrangian method can lead to numerical instability and computational difficulties in non-uniform media with multiple sinks/sources and complex boundary conditions. The mixed Eulerian-Lagrangian approach attempts to combine the advantages of both the Eulerian and the Lagrangian approaches by solving the advection term with a Lagrangian method and the dispersion and reaction terms with an Eulerian method.

The numerical solution implemented in MT3D is a mixed Eulerian-Lagrangian method. The Lagrangian part of the method, used for solving the advection term, employs the forward tracking method of characteristics (MOC), the backward-tracking modified method of characteristics (MMOC), or a hybrid of these two methods. The Eulerian part of the method, used for solving the dispersion and chemical reaction terms, utilises a conventional block-centered finite-difference method.

The MT3D transport model uses an explicit version of the block-centered finite-difference method to solve the dispersion and chemical reaction terms. The limitation of an explicit scheme is that there is a certain stability criterion associated with it, so that the size of time steps cannot exceed a certain value. However, the use of an explicit scheme is justified by the fact that it saves a large amount of computer memory which would be required by a matrix solver used in an implicit scheme. In addition, for many advection-dominated problems, the size of transport steps is dictated by the advection process, so that the stability criterion associated with the scheme for the dispersion and reaction processes is not a factor. It should be noted that a solution package based on implicit schemes for solving dispersion and reactions could easily be developed and added to the model as an alternative solver for mainframes, more powerful personal computers, or workstations with less restrictive memory constraints.

The implementation of the model consists of the following actions: definition of the simulation problem, specification of the initial and boundary conditions, determination of appropriate transport step size, preparation of global mass balance information and output of simulation results. The Flow Model Interface Package interfaces with a flow model to obtain the flow solution from the flow model.

Output options. The program generates a standard output file and several optional output files. The standard output file is generated every time the model is run. The optional output files are generated only if they are requested. The amount, type, and frequency of information to be written on the output files are controlled by the user-specified options in the input file to the Basic Transport Package. The functions of these output files are listed below:

    • The Standard Output File: contains echo of input data to ensure they are read in properly, printout of simulated concentrations, and some other essential information such as model-calculated parameters and mass budget at user-specified times and at the end of the simulation.
    • The Unformatted Concentration File (default name: MT3D.UCN): contains concentrations saved at user-selected times in the unformatted form for a continuation run or for post-processing purposes.
    • The Observation Point File (default name: NMD.OBS): contains concentrations versus total elapsed time at user-specified observation points at every transport step.
    • The Mass Balance Summary File (default name: MT3D.MAS): contains a one-line summary of mass budget at every transport step for checking purposes.
4.5 The Coastal water model (Princeton Ocean Model - POM)

The Princeton Ocean Model (POM) is an ocean circulation numerical model designed by A.Blumberg and G.Mellor (1987) for both coastal and open ocean studies. It is a public domain model that is being used by a large number of research and academic institutes all over the world. Due to its ability to simulate both shallow water and deep ocean dynamics, it has been used for a variety of applications ranging from small scale coastal management problems to general circulation studies of the Atlantic Ocean. A full description of the numerical scheme can be found in Blumberg and Mellor (1987), whereas a web site dedicated to the model is available on POM home page (

4.5.1 Main features

POM is a three-dimensional, primitive equation numerical model. The prognostic variables are the three components of the velocity U-V-W, the temperature T and the salinity S fields. The equation of state is used for the computation of potential density. Two more prognostic equations are used to calculate turbulent kinetic energy and turbulent macroscale. These equations are part of the Mellor - Yamada 2.5 turbulence closure scheme used for the calculation of vertical diffusivity. Horizontal diffusivities are calculated according to the Smagorinsky formula. A set of vertical integrated equations of continuity and motion are also solved to provide free surface variations.

The principal attributes of the model are as follows:

    • It contains an embedded second moment turbulence closure sub-model to provide vertical mixing coefficients
    • Horizontal diffusion terms are computed using the Smagorinsky (1963) formula that relates these terms to the scales of motion being resolved in the model and to the local deformation field.
    • It is a sigma coordinate model in that the vertical coordinate is scaled on the water column depth
    • The horizontal grid uses curvilinear orthogonal coordinates and a "Arakawa C" differencing scheme
    • The horizontal time differencing is explicit whereas the vertical differencing is implicit. The latter eliminates time constraints for the vertical coordinate and permits the use of fine vertical resolution in the surface and bottom boundary layers.
    • The model has a free surface and split time step. The external model portion of the model is two-dimensional and uses a short time step based on the Courant-Friedrichs-Levy (CFL) stability condition and the external wave speed. The internal mode is a three-dimensional and uses a long time step based on CFL condition and the internal wave speed.
    • Complete thermodynamics have been implemented.

By and large, the 2.5 turbulence closure scheme (Mellor and Yamada,1982) seems to do a fair job simulating mixed layer dynamics although there have been indications that calculated mixed layer depths are a bit too shallow (Martin, 1985). At least a part of the deficit is due to internal wave velocity shears (Mellor, 1989). Also, wind forcing may be spatially smoothed and temporally smoothed. It is known that the latter process will reduce mixed layer thickness (Klein, 1980). Further study is required to quantify these effects.

The sigma coordinate system is probably a necessary attribute in dealing with significant topographic variability such as that encountered at estuaries or over continental shelf breaks and slopes. Together with the turbulence sub-model, the model produces realistic bottom boundary layers which are important in coastal waters (Mellor and Blumberg, 1985) and in tidally driven estuaries (Oey et. al., 1985a,b) which the model can simulate since it does have a free surface. More recently, it was found that bottom topography layers are important for deep water formation processes (Zavatarelli and Mellor, 1995, Jungclaus and Mellor, 1996) and possibly for the maintenance of the baroclinicity of oceans basins (Mellor and Wang, 1996).

One well-known problem with the sigma coordinate models is the error introduced by the pressure gradient in areas with steep topography. In POM the problem is handled by subtracting a mean density profile before computing the density gradients. It has been shown (Mellor et. al., 1994) that using this technique and by adjusting the topography so that dxH/H is small, the horizontal pressure gradient has vanished.

The horizontal finite difference scheme is staggered and, in the literature, has been called an Arakawa C-grid. The horizontal grid is a curvilinear coordinate system, or as a special case, a rectilinear coordinate system may be easily implemented.

4.5.2 Boundary conditions

The surface boundary conditions are surface heat flux, surface salinity flux and wind stress. These fluxes can be either explicitly prescribed or interactively computed by the model using bulk formulas. The bulk formulas are used to compute these fluxes using the model's sea surface temperature (SST) and atmospheric parameters, namely, the wind field w, air temperature Ta, relative humidity Rh, cloud cover C and precipitation P. This later approach is a modern technique developed during the last decade when the fast development of computer technology made feasible the production of reliable prognostic or diagnostic atmospheric model results. This interactive way of flux computation allows simulation of feedback mechanisms that are known to play a crucial role on air-sea interaction

Finally, the model can handle open boundaries through appropriate user defined boundary conditions. A Sommerfeld radiation condition is used for the normal to the open boundary internal velocity. For the tangential component of velocity and the scalar quantities, salinity and temperature, an upstream advection equation is used. When there is inflow through the open boundary, temperature and salinity are prescribed form the initial data. For the external mode a free radiation condition is used for the normal mean velocities and a zero-gradient condition for the free surface elevation.

5 Model input requirements (air domain)

The quality of the simulation results depend heavily on the accuracy of the input data. The preparation of these data has to be carried out with great care. There are three main classes of data:

    • static data (digital elevation map)
    • slowly varying data (land use)
    • rapid varying data (nearly all emissions, boundary concentration values, synoptic meteorological conditions).

Another important task is the proper selection of the model domain. Even if the area of interest is relatively small it has to be ensured that all orographic features affecting the mesoscale wind field are included in the model domain. For example, isolated summits should not lay directly at the border of the model domain. If the region of interest is affected by sea breeze, a considerable amount of water has to be added to the model domain. If there are strong emission sources close to the planned model area the model domain should be enlarged to include them.

5.1 Digital elevation map and land use data

The elevation map is necessary to compute the wind flow. As the model cannot resolve the subgrid scale features, the orography heights above sea level must be provided in the chosen computation grid as grid averages.

The land use data influences the computation of the thermal properties of the surface and of ground resistance. Additionally, the surface uptake processes and the biogenic emissions are determined by land use parameters. Because the land use or at least the agricultural use may vary weekly, seasonally or yearly, a regular update may be necessary. The most accurate means of data preparation is the use of satellite images. All land use types significant for the model domain should be included in the inventory. The following types are distinguished and must be provided in percent per grid:

    • sea water
    • fresh water
    • arid land
    • low vegetation
    • grassland
    • farmland
    • coniferous forest
    • mixed forest
    • deciduous forest
    • sealed area
    • suburban area
    • urban area
5.2 Emission inventory

The analysis of the air quality of an urban region by means of simulation requires a detailed emission inventory. Especially in order to calculate ozone, which is really a secondary air pollutant, as much as possible information about precursor emissions should be collected. For an ozone calculation, the inventory should contain at least emission rates for nitrogen oxides and volatile organic compounds. For investigating nonreactive transport processes, the user can select his substance of interest.

The anthropogenic emissions can be given as point sources with a concrete coordinate or they can be aggregated on any grid size smaller than the selected computational grid. Usually, the environmental agencies collect the data as yearly emission rates. So the typical unit for an emission source in the ECOSIM emission data base is kg per year. Every emitter can have its own emission height or the emission height of its emission class which can be a concrete height or might be distributed over a couple of levels. The emissions sources are grouped together in order to simplify the accessibility for scenario analysis. At the recent stage, the ECOSIM emission inventory distinguishes 7 different emission classes:

    • power stations and large emitters (emissions from high stacks)
    • small industries and crafts
    • households
    • moving" traffic (except of city centre traffic)
    • city centre traffic
    • traffic evaporation
    • air port and air craft emissions

The ECOSIM emission inventory must include at least one of these emission classes and one substance. The inventory groups NOx (sum of all nitrogen oxide emissions) and VOC (volatile organic compounds without methane) into one database, respectively. The user can select the appropriate relation between the species (e.g. the NO2/NO split) for every emission class or can simply use the default values reflecting typical conditions.

For episode simulations, and especially for forecast purposes, the dynamic behaviour of emissions must be considered. Every emission class for every substance can have its own diurnal variation pattern for

    • normal working days
    • saturday
    • sunday.

The system automatically searches for the day of the week and uses its dynamics. Ozone precursors have both anthropogenic and biogenic sources. While anthropogenic emissions are given by the emission inventory the biogenic emissions must be computed on-line. Especially for hot summer episodes, biogenic emissions of volatile organic compounds can easily reach the same total amount like anthropogenic emissions. These emissions are a function of vegetation type, temperature and can be insulation dependent. The model system uses the approach given by Pierce [Pierce 1990].

5.3 Initial and background data
5.3.1 Meteorological data

The prognostic model MEMO is a set of partial differential equations in 3 spatial directions and in time. To solve these equations, information on the initial state in the whole domain and on the development of all relevant quantities at the lateral boundaries is required.

1) Initial state. To generate an initial state for the prognostic model, a diagnostic model [Kunz 1991] is applied using measured temperature and wind data. Both temperature and wind data can be provided as:

    • surface measurements i.e. single measurement directly above the surface (not necessary)
    • upper air soundings i.e. soundings that consist of several (at least two) measurements at different height levels at a constant geographical location (at least one sounding for temperature and wind velocity is necessary)

2) Time dependent boundary conditions. Information on quantities at the lateral boundaries can be taken into account as surface measurements and upper air soundings.

The forecast model REGOZON uses the geostrophic wind, a vertical temperature profile obtained by an upper air sounding, cloud cover fraction, water temperature and the soil temperature in 1 m depth. The implemented model in the region of Berlin uses predefined interfaces to observations given in WMO (World Meteorological Organisation) format. It can take into account time dependent synoptic conditions. Interfaces have been established to incooperate weather forecast results in 6h time steps.

5.3.2 Background concentrations

Air quality models in the urban scale should include mean concentration values from outside the model domain as background data. For example, ozone could exhibit a relatively high large scale background level with a strong diurnal variation near the ground. If these background levels would be neglected, even the computation of additional ozone production would fail due to the nonlinear character of the ozone chemistry. For DYMOS, the background values of all relevant substances can be defined separately for every layer. Also, these background values can vary in any desired time interval.

For a number of simulation purposes, there is no knowledge of the time dependencies of the background concentrations. But especially surface ozone, and also its background, shows a rather strong diurnal variation. To allow a more realistic simulation, the model can run on request in a box model background mode. In this mode, the model uses the background values only for initialising the concentrations. A box model with a assumed infinite horizontal extension performs the determination of the time dependent background values. The vertical structure and the formulation of the vertical processes in the box is similar to the conditions in the model domain; the box model uses averaged meteorological values from the model domain as time depending input parameters.

6 Model output (air domain)

The MEMO output contains the fields of all relevant meteorological variables such as wind, pressure, diffusion coefficients and temperature. These data serve as input for the DYMOS model which calculates the transport, chemical changes, deposition and biogenic emissions on the basis of the meteorological data. The meteorological fields from MEMO are usually provided every hour. If the data transfer between MEMO and DYMOS is no longer a critical problem, the MEMO output will be used in shorter time intervals.

The final output from DYMOS and REGOZON with relevant information for the user consists of the desired concentration field. The user can specify the substance(s) of interest, the output interval and the number of levels from the ground he wants to see. The graphic representation of the results uses all common features like isoline plots, contour plots or grid representation in selected sensitivity band. All simulation results can be stored together with the case descriptor and the input data set.

7 Computer user notes and installation guide

For the hardware and software requirements of the demonstrator the user is referred to Section 2.4. MEMO and REGOZON require at least a fast workstation (e.g. SUN, Hewlett Packard, IBM, DEC) with sufficient RAM depending on the numerical grid.

For an efficient and reliable execution of the DYMOS system a parallel machine running PVM [Geist 1994] or a workstation cluster of at least 4 CPUs is strongly recommended. The parallel program version is able to guarantee an execution in a fraction of real time.

The ECOSIM Demonstrator is distributed as a compressed (gzip) tar archive on tape, or can be downloaded from the ESS ftp server. General instruction for the installation of an ESS application can be found in the AirWare Reference Manual

Installation of the model servers. In the special case that the model and database servers are implemented on the same machine as the ECOSIM Demonstrator (e.g., in the case of the Berlin Demonstrator), the corresponding files are available in a separate tar file, www.tar. The application server environment consists of the http daemon and individual models or data servers (e.g., regozon and blume in the case of Berlin). To install the server environment,

    • move the file www.tar to a appropriate location on the target machine;
    • extract the tar archive with: tar xvf www.tar

The http daemon is driven by a number of configuration files (www/config/*.conf). The following variables are to be redefined in the configuration files (according to the new location or path):

    • ServerRoot
    • CacheRoot
    • Pass
    • Exec

To start http daemon use: httpd_3.0A -v -r config/httpd.conf -p 8080. The individual servers, models or data servers (regozon, blume.exe) are in the directory www/cgi-bin. Blume data files are in the directory www/blume. Regozon is driven by its configuration file www/fort.94.

8 The Demonstrator implementation

8.1 REGOZON forecasting model and the BLUME network

The REGOZON ozone forecasting model is implemented as a remote server; the executable regozon is located in the directory ./www/bin

For details of the model operation, assumptions, coefficients etc., please consult the REGOZON description on the ECOSIM web pages. In the Berlin Demonstrator, REGOZON is implemented with the following constraints:

    • The models run with a fixed 30 minutes output time step to allow direct comparison with the BLUME monitoring network observation data.
    • The first model output for a daily (22 hours) run is provided for 01:30.
    • The models run for a maximum of 22 hours.
    • The model uses a fixed 2 km grid, for a 100 km by 100 km domain.

Please note that all static input data that are NOT edited by the ECOSIM interface (see below) are provided at the server side; this includes meteorological forecasts in the case of running REGOZON in forecasting mode, i.e., into the future (today or tomorrow).

REGOZON functionality. The basic functionality includes the following options:

    • editing and submitting new simulation requests;
    • display and analysis of a simulation scenario result;
    • comparison of a simulation scenario with observation data.
8.1.1 REGOZON forecasting model: scenario editing

The CREATE button in the main REGOZON interface leads to the scenario editing functions. This will prepare a new scenario for simulation by the REGOZON server. The current date is used as a default, and the previous scenario as a TEMPLATE to generate a new one. Please note that a new name has to be specified for the new scenario, as the system will not overwrite an existing one.

Changing the scenario name. To change the scenario name, position the mouse pointer in the corresponding window, and use the keyboard to enter a new name; use the Delete or Backspace keys to delete the old string.

Editing scenario parameters. The following scenario parameters can be set or selected by the user:

    • set of output variables to be computed (toggle on/off);
    • set of output layers to be displayed (toggle on/off);
    • date of the simulation or forecast (the default date for any new scenario is set to the current date);
    • duration of the simulation
    • name of the new scenario;
    • scaling of the emission matrices.

To change a parameter, position the mouse pointer over its window and click with the left mouse button; this will either toggle the value (on/off) or create a pop-up window with a set of alternatives to choose from. Text strings like the scenario name can be edited with the keyboard.

Emission matrices. To characterize the emission scenario for a simulation, the system offers a set of emission matrices that can be scaled. These matrices cover:

    • the basic substances: NO2, VOC, CO;
    • the basic source categories, i.e.:
      • traffic;
      • inner city traffic;
      • high industrial stacks;
      • low industrial stacks;
      • households;
      • airports.

The user can select and display in one of the three smaller windows any one of these matrices, and specify a scaling factor for them; the fourth window shows their sum total, i.e., substance specific total annual emissions (in tons/year).

To select a substance, click its name and total emission value in the lower right selector window labelled Total Emissions of the main editor pop-up tool. To edit the multiplier, the standard editor tool of the embedded expert system is used.

8.1.2 REGOZON forecasting model: scenario display

When entering the REGOZON forecasting model page, a default scenario specified in the ./defaults data file neco is loaded.

The LOAD button triggers a file selector tool, that offers all available scenarios, including those that are pending (still being computed) or terminated with an error message. The display control options for a given scenario include:

    • selection of the time step with the tape deck, including the possibility of continuous animation;
    • the selection of the variable to be displayed (O3 or NO2 for REGOZON);
    • the selection of a vertical layer. The latter is selected by clicking on the appropriate symbol in the layer selector.

The display of the model results uses a color coded map overlay; the color coding of the display can be controlled with the cell-grid controller tool. The controller supports:

    • selection of the color ramp;
    • selection of display style (opaque, translucent, isolines, interpolated);
    • selection of upper and lower display bounds.
    • In addition, the main icon menu supports:
    • zooming in and out of the map display;
    • connection to the BLUME monitoring data server to load monitoring data for comparison;
    • pseudo 3D display of the concentration field;
    • monitoring data comparison;
    • GIS link for the selection of background maps layers;
    • standard info and exit icons.

3D display. The 3D icon triggers the 3D display of the current concentration overlay. The concentration value will be shown as the vertical dimension of the model result grid, using the same color coding as in the 2D display. Vertical positioning, exaggeration, and 3D rotation of the image can be controlled interactively.

8.1.3 BLUME monitoring data acquisition

Please note that this option is only available for the BERLIN implementation. The item: BLUME Monitoring Network in the case-study and function selector leads to the Monitoring network interface. The BLUME interface is designed to provide the following functions:

    • download a new data set for O3 and NO2 (a day at a time) from the BLUME database server;
    • pre-view that daily data set for quality control;
    • integrate the data set with the overall monitoring time-series data base of ECOSIM.

Scenario information. The current date and variable are shown in the upper right corner of the screen, corresponding to the spatially interpolated concentration data shown in the map window. By default, the data set of the last (successful) download is shown.

Animation control. The time series data can be displayed one observation set at time, or in continuous animation. This is controlled by the buttons of the tape deck. In parallel, the analog clock indicates the time of the measurement.

Station and variable selection. For the display of the data, the user can select:

    • the variable (currently: ozone or nitrogen-dioxide)
    • an individual station for selective display, see below).

For this, two parallel selector boxes are offered on the screen.

Station current value display list. For the current time-step, as indicated by the clock above, the concentration values of the current substance for the individual stations are displayed in a list of stations.

Individual station time-series display. For the individual stations selected (see above), the time series for the current day is shown as 48 half-hourly values.

Concentration histogram (cellgrid-controller). The histogram tool displays the (spatial) frequency distribution of concentration values (classes) corresponding to the large map window with the interpolated measurement data.

Main icon menu. The main icon menu of the BLUME interface offers the following functions:

    • zoom in, zoom out
    • connect to the BLUME database server
    • object-data base and time-series display
    • data import: attach data to the time-series database
    • embedded GIS link
    • help/info hypertext
    • return/EXIT
8.2 The MEMO/MUSE meteorological and air quality models

The MEMO meteorological forecasting model is implemented as a remote server; the corresponding URL is specified in the default file ./defaults/neco in the Athens Demonstrator version. For details of the model operation, assumptionds, coefficients etc., please consult the MEMO description on the ECOSIM web pages. In the Athens Demonstrator, MEMO is implemented with the following features and constraints:

    • the model is implemented for a domain of 72 by 72 kilometres with a 2 km grid (36 by 36 cells of 2 by 2 km);
    • meteorological conditions are limited to a set of predefined, typical scenarios.

Please note that all static input data that are NOT edited by the MEMO interface (see below) are provided at the server side; this includes initial and boundary conditions of the model, and any dynamic forcing applied.

MEMO functionality. The basic functionality includes the following options:

    • editing and submitting new simulation requests;
    • display and analysis of a simulation scenario result (same as for REGOZON).
8.2.1 MEMO meteorological model: scenario editing

The CREATE button in the main MEMO interface leads to the scenario editing functions. This will prepare a new scenario for simulation by the MEMO server. The current date is used as a default, and the previous scenario as a TEMPLATE to generate a new one. Please note that a new name has to be specified for the new scenario, as the system will not overwrite an existing one.

Changing the scenario name. To change the scenario name, position the mouse pointer in the corresponding window, and use the keyboard to enter a new name; use the Delete or Backspace keys to delete the old string.

Editing scenario parameters. The following scenario parameters can be set or selected by the user:

    • set of output variables to be computed (toggle on/off);
    • set of output layers to be displayed (toggle on/off);
    • a symbolic meteorological scenario;
    • duration of the simulation;
    • time step of the simulation;
    • name of the new scenario;
    • scaling of the emission matrices.

To change a parameter, position the mouse pointer over its window and click with the left mouse button; this will either toggle the value (on/off) or create a pop-up window with a set of alternatives to choose from. Text strings like the scenario name can be edited with the keyboard.

Emission matrices. To characterize the emission scenario for a simulation, the system offers a set of emission matrices that can be scaled. These matrices cover:

    • the basic substances: NO2, VOC, CO;
    • the basic source categories, i.e.:
      • traffic;
      • inner city traffic;
      • high industrial stacks;
      • low industrial stacks;
      • households;
      • airports.

To select a substance, click its name and total emission value in the lower right selector window of the main editor pop-up tool.

The user can select and display in one of the three smaller windows any one of these matrices, and specify a scaling factor for them; the fourth window shows their sum total, i.e., substance specific total annual emissions (in tons/year). To edit the multiplier, the standard editor tool of the embedded expert system is used.

8.2.2 MEMO/MUSE: scenario display

Please refer to chapter 8.1.2. REGOZON forecasting model: scenario display.

8.3 The POM (Princeton Ocean Model) coastal water quality model

For details of the model operation, assumption, coefficients etc., please consult the POM description on the ECOSIM web pages. In the Athens Demonstrator, POM is implemented with the following features and constraints:

    • the model is implemented for a domain of 100 by 100 kilometres with a 2 km grid (50 by 50 cells of 2 by 2 km);
    • meteorological and oceanographical conditions are limited to a set of predefined monthly scenarios;
    • the model runs are limited to a set of predefined length, ranging from 1 to 28 days;
    • only a conservative tracer is supported;
    • there is an option to relocate sources, but they have to be predefined in TEMPLATES (or scenario examples) since the interface does NOT allow to add a new source.

Please note that all static input data that are NOT edited by the POM interface, see below, are provided at the server side; this includes initial and boundary conditions of the model, and any dynamic forcing applied.

POM can be implemented as a remote server or as a local executable, with the same client-server protocol, which requires the installation of a local www server (httpd), see: ECOSIM Installation Guide.

The corresponding server URL is specified in the default file ./defaults/neco in the Athens Demonstrator version:

    • .cliser.sernames.n: 3
    • Athens MEMO
    • .cliser.sernames.node.1:
    • MP1 POM
    • .cliser.sernames.node.2: mp1:8080/cgi-bin/sarpom
    • Athens GW
    • .cliser.sernames.node.3:

To modify the list of servers and potential alternates and backups, modify the node name in the default file; when adding a server, please make sure to also update .cliser.sernames.n accordingly.

POM functionality. The basic functionality includes the following options:

    • editing and submitting new simulation requests;
    • display and analysis of a simulation scenario result (same as for REGOZON).
8.4 The ECOSIM groundwater quality model

For details of the model operation, assumption, coefficients etc., please consult the Groundwater Models description on the ECOSIM web pages. In the Athens Demonstrator, the Groundwater Model is implemented with the following features and constraints:

    • Steady state simulation
    • One unconfined aquifer (limestone - dolomites)
    • All aquifer properties are isotropic and homogeneous

Please note that all static input data that are NOT edited by the ECOSIM interface, see below, are provided at the server side; this includes initial and boundary conditions of the model, and any dynamic forcing applied. The groundwater model is implemented as a remote server through the http protocol. The corresponding server URL is specified in the default file ./defaults/neco in the Athens Demonstrator version:

    • Athens GW
    • .cliser.sernames.node.3:

To modify the list of servers and potential alternates and backups, modify the node name in the default file; when adding a server, please make sure to also update .cliser.sernames.n accordingly.

Groundwater Model functionality. The basic functionality includes the following options:

    • editing and submitting new simulation requests; (the parameters that can be set are an initial concentration and the duration of the model run);
    • display of the scenario result (a single static image representing the situation at the end of the specified run-time).

9 Bibliography

Blackadar, A. K. (1962) The Vertical Distribution of Wind and Turbulent Exchange in a Neutral Atmosphere, J. of Geophys. Res 67, 3095-3102.

Blumberg, A.F. and Mellor G.L. (1987) A description of a three dimensional coastal ocean circulation model, in Three Dimensional Coastal Ocean Circulation Models, Coastal Estuarine Sci., vol. 4, edited by N.S. Heaps, pp. 1-16, AGU, Washinghton D.C.

Bottenheim, J.W. & Strausz, O.P. Modelling study of a chemically reactive power plant plume, Atmos. Environ., 1982, 16, 85-97.

Carpenter, K.M. (1982) Note on the paper 'Radiational condition for the lateral boundaries of limited-area numerical models' by Miller, M.J. and Thorpe, A.J. (Q.J. 107, 615-628), Quart. J. R. Met. Soc.108, 717-719.

Deardorff, J. W. Three-dimensional numerical study of the height and the mean structure of a heated planetary boundary layer, Bound. Layer Meteorol., 1974, 7, 81-106.

Eppel D. P., Kapitza H., Claussen M., Jacob D., Koch W., Levkov L., Mengelkamp H-T., Werrmann N. (1995) The Non-Hydrostatic Mesoscale Model GESIMA. Part II: Parametrizations and Applications. Beitr. Phys. Atmosph. 68, 15-41.

EUROTRAC (1992) Annual Report 1991, Part 5.

Flassak, Th. Ein nicht-hydrostatisches mesoskaliges Modell zur Beschreibung der Dynamik der planetaren Grenzschicht, Fortschr.-Ber. VDI Reihe 15, Nr. 74, VDI-Verlag, Düsseldorf, 1990.

Flassak, Th. and Moussiopoulos, N. (1988) Direct solution of the Helmholtz equation using Fourier analysis on the CYBER 205, Environmental Software 3, 12-16.

Geist Al, Beguelin A., Dongarra J., Jiang W., Manchek R., Sunderam V. (1994) PVM: Parallel Virtual Machine, A User`s Guide and Tutorial for Networked Parallel Computing. The MIT Press, Cambridge, Massachusetts

Gerlach J., Schmidt M. Parallele Implementierung von Modellen für die Luftschadstoffanalyse. Proc. 8. Symposium Simulationstechnik Berlin, 1993.

Gery M. W., Whitten G. Z. and Killus J. P. (1988) Development and Testing of the CBM-4 for Urban and Regional Modeling. EPA-600/3-88-012, U.S. EPA.

Graf, J. & Moussiopoulos, N. Intercomparison of two models for the dispersion of chemically reacting pollutants, Contr. Phys. Atmos., 1991, 64, 13-25.

Harten, A. (1986) On a large time-step high resolution scheme, Math. Comp. 46, 379-399.

Heimann D. (1985) Ein Dreischichten-Modell zur Berechnung mesoskaliger Wind- und Immissionsfelder über komplexem Gelände. Dissertation, University of Munich.

Helmis, K. G. Research on the impact on air quality of the Greater Athens Area due separately to the operation of the new airport at Spata and the airport at Hellenikon, Technical Report, University of Athens, 1994.

Jungclaus, J. H. and G. L. Mellor (1996)A three dimensional model study of the Mediterranean outflow, Journal of Marine Systems, Submitted

Kapitza H. and Eppel D. P. (1992) The Non-Hydrostatic Mesoscale Model GESIMA. Part I: Dynamical Equations and Tests. Beitr. Phys. Atmosph. 65, 129-146.

Kessler, Chr. Entwicklung eines effizienten Lösungsverfahrens zur modellmäßigen Beschreibung der Ausbreitung und Umwandlung reaktiver Luftschadstoffe, PhD Dissertation, Universität Karlsruhe, 1995.

Klein P.A. (1980)A simulation of the effects of air-sea transfer variability on the structure of the marine upper layers, Journal of Physical Oceanography, 10, 1824-1841.

Klemp, J.B. and Durran, D.R. (1983) An upper boundary condition permitting internal gravity wave radiation in numerical mesoscale models, Mon. Weather Rev.111, 430-444.

Kunz, R. (1991) Enwicklung eines diagnostischen Windmodells zur Berechnung des Anfangszustandes für das dynamische Grenzschichtmodell MEMO, Diplomarbeit Universität Karlsruhe.

Liu, Peek, Jones, Buus, and Nye (1994) Managing Internet Information Services. 630 pp., O`Reilly & Associates, Inc. Sebastopol, CA. ISBN 1-56592-062-7.

Lux, Th.; Mieth, P.; Schäfer, R.-P.; Schmidt, M.; Sydow, A.; Unger, S.: "Parallel Computer Application: Integrated Air Pollution Analysis", Proc. IMACS/IFAC 2nd International Symposium on Mathematical and Intelligent Models in System Simulation, Brussels,, 1993

Martin P.J. (1985) Simulation of the mixed layer at OWS November and Papa with several models, Journal of Geophysical Research, 90, 903-916

Mellor G.L., (1989) Retrospect on oceanic boundary layer modeling an second moment closure, Hawaiian Winter workshop on "Parameterization of Small-Scale Processes", January 1989, University of Hawaii, Honolulu.

Mellor G.L. and Blumberg A.F., (1985) Modeling vertical and horizontal diffusivities with the sigma coordinate system, Mon. Wea. Rev., 113, 1380-1383

Mellor G.L., Ezer T., Oey L.Y., (1994) The pressure gradient conundrum of sigma coordinate models, J. Atmos. Oceanic Technol., 11(4), 1126-1134.

Mellor, G. L. and Wang X. H., (1996) Pressure compensation and the bottom boundary layer, Journal of Physical Oceanography, 26(10), 2214-2222

Mellor G.L and Yamada T., (1982) : Development of a turbulence closure model for geophysical fluid problems, Rev. Geophys., 20, 851-875.

Mieth P. and Unger S. (1993) Szenarienrechnungen und Fallstudien zur Ozonanalyse im Berliner Raum mittels eines mesoskaligen Ausbreitungsmodells. GMD FIRST im Auftrag der Senatsverwaltung für Stadtentwicklung und Umweltschutz Berlin und des Ministeriums für Umwelt, Naturschutz und Raumordnung des Landes Brandenburg, Berlin.

Mieth P. and Unger S. (1996) Comparison of identical emission reduction measures under different meteorological conditions, in Air Pollution IV (ed B. Caussade, H. Power, C.A. Brebbia), Proceedings of the Air Pollution `96 Conference, Toulouse, France, 1996, Computational Mechanics Publications, Southampton, 1996

Moussiopoulos, N. The EUMAC Zooming Model, a tool for local-to-regional air quality studies, Meteorology and Applied Physics, in press.

Moussiopoulos, N. Mathematische Modellierung mesoskaliger Ausbreitung in der Atmosphäre, Fortschr.-Ber. VDI Reihe 15, Nr. 64, VDI-Verlag, Düsseldorf, 1989.

Moussiopoulos, N. Atmospheric Transport Phenomena, Giahoudi - Giapouli O.E, Thessaloniki, 1991, in Greek.

Moussiopoulos, N., Sahm, P., Gikas, A., Karagiannidis, K., Karatzas, K. & Papalexiou, S. Analysis of air pollutant transport in the Athens basin and in the Spata area with a three-dimensional dispersion model, in Air Pollution `95 (ed C.A. Brebbia, N. Moussiopoulos & H. Power), Proceedings of the Air Pollution `95 Conference, Porto Carras, Greece, 1995, Computational Mechanics Publications, Southampton, 1995.

Moussiopoulos, N. et al. Analysis of air mass transport in the Athens basin and in the Spata area with the non hydrostatic mesoscale model MEMO, Technical Report of the LHTEE, No 94092, Aristotle University Thessaloniki, 1994.

Mousiopoulos, N. et al. Impact on the air quality due to the transfer of Athens airport from Hellenikon to Spata, Technical Report of the LHTEE, No 94091, Aristotle University Thessaloniki, 1994.

Mousiopoulos, N. et al. Analysis of the impact of the operation of the Athens airport on the photochemical smog formation in the Athens basin, Technical Report of the LHTEE, No 94113, Aristotle University Thessaloniki, 1994.

Moussiopoulos, N., Flassak, Th., Berlowitz, D., Sahm, P. (1993) Simulations of the Wind Field in Athens With the Nonhydrostatic Mesoscale Model MEMO, Environmental Software 8, 29-42.

Moussiopoulos, N., Sahm, P., Kessler, Chr., Kunz, R. (1994) Simulation of Mesoscale Wind Flow and Photosmog Formation in the Greater Athens Area, in The EUMAC Zooming Model (EZM) (Moussiopoulos, N. ed.), EUMAC reports, in print.

Moussiopoulos, N. (1987) An efficient scheme to calculate radiative transfer in mesoscale models, Environmental Software 2, 172-191.

Moussiopoulos, N. (1989) Mathematische Modellierung mesoskaliger Ausbreitung in der Atmosphäre, Fortschr.-Ber. VDI, Reihe 15, Nr. 64, pp. 307.

Moussiopoulos, N. and Flassak, Th. (1989) A fully vectorized fast direct solver of the Helmholtz equation, in: Applications of supercomputers in engineering: Algorithms, computer systems and user experience (Brebbia, C.A. and Peters, A., eds.), Elsevier, Amsterdam 67-77.

Ngyuen, Th. H. Ein Vier-Schichten-Modell zur Berechnung der Aus­breitung und der chemischen Umwandlung von Schadstoffen in der Atmos­phäre, Fortschr.-Ber. VDI Reihe 15, Nr. 66, VDI-Verlag, Düsseldorf, 1989

Oey L.Y., Mellor G.L. and Hires R.I. (1985a) A three dimensional simulation of the Hudson-Raritan estuary, Part I: Description of the model and model simulations, Journal of Physical Oceanography, 15, 1676-1692

Oey L.Y., Mellor G.L. and Hires R.I. (1985b) A three dimensional simulation of the Hudson-Raritan estuary, Part II: Comparison with observation, Journal of Physical Oceanography, 15, 1693-1709

Orlanski, J. (1976) A simple boundary condition for unbounded hyperbolic flows, J. Comput. Phys. 21, 251-269.

Pielke, R. A. Mesoscale Meteorological Modelling, Academic Press, 1984.

Pierce T. E., Lamb B. K. and Van Meter A. R. (1990) Development of a Biogenic Emissions Inventory System for Regional Scale Air Pollution Models. Air & Waste Management Association, Presentation at the 83rd Annual Meeting & Exhibition, Pittsburgh.

Sahm P. and Moussiopoulos N. (1995), MUSE - a multilayer dispersion model for reactive pollutants, in Air Pollution III (H. Power, N. Moussiopoulos and C.A. Brebbia, eds), Computational Mechanics Publications, Southampton, Vol. 1, 359-368.

Schumann, U. and Volkert, H. (1984) Three-dimensional mass- and momentum- consistent Helmholtz-equation in terrain-following coordinates, in Notes on Numerical Fluid Mechanics, Vol. 10 (Hackbusch, W., ed.), Vieweg, Braunschweig, 109-131.

Smagorinsky J. (1963) General circulation experiments with the primitive equations, I , The basic experiment, Mon. Weather Rev., 91, 99-164

Smolarkiewicz, P. K. A fully multidimensional positive definite advection transport algorithm with small implicit diffusion, J. Comput. Phys., 1984, 54, 325-362.

Sydow, A.; Schmidt, M.; Lux, Th.; Schäfer, R.-P.; Mieth, P.: "Air Pollution Modelling and Simulation (an Approach of Parallel Simulation)", Proc. EUROSIM Simulation Congress, Capri, Italy, 1993

Sydow, A.; Schmidt, M.; Unger, S.; Lux, Th.; Mieth, P.; Schäfer, R.-P.: "Parallel Simulation of Air Pollution Transport for Urban Decision Making", Proc. 12th IFAC World Congress, July 19 - 23, 1993, Sydney, Australia, 1993

Sydow, A.; Schmidt, M.; Lux ,Th. ;Schäfer,R.-P. ; Mieth, P. : "Simulation of Air Pollutant Dispersion on Parallel Hardware. Simulation Practice and Theory 1", 57-64, 1993

Sydow, A.; Lux, Th.; Mieth, P.; Schmidt, M.; Unger, S.: "Computational Aspects of Dynamic Models for Simulating Transport and Chemical Conversion of Air Pollutants", Proc. 14th IMACS World Congress on Computation and Applied Mathematics, Atlanta, Georgia, USA, 1994

Toll Abello, I., Moussiopoulos, N., Sahm, P., Samaras, Z., Zachariadis, T. & Aslanoglou, M. A new emissions inventory for Athens and its use to improve quality predictions, in Air Pollution `95 (ed C.A. Brebbia, N. Moussiopoulos & H. Power), Proceedings of the Air Pollution `95 Conference, Porto Carras, Greece, 1995, Computational Mechanics Publications, Southampton, 1995, in press.

Wesely M. L. (1988) Improved Parameterizations for Surface Resistance to Gaseous Dry Deposition in Regional-Scale, Numerical Models. EPA-600/3-88-025. U.S. EPA.

Zavatarelli M. and Mellor G.L. (1995) A Numerical Study of the Mediterrannean Sea Circulation, Journal of Physical Oceanography, 25, 1384-1414.

Copyright 1995-2002 by:   ESS   Environmental Software and Services GmbH AUSTRIA