A large-scale global air quality model, running efficiently on a single vector processor, is enhanced to make more realistic and more long-term simulations feasible. Two strategies are combined: non-uniform grids and parallel processing. The communication through the hierarchy of non-uniform grids interferes with the inter-processor communication. We discuss load balance in the decomposition of the domain, I/O, and inter-processor communication. A model shows that the communication overhead for both techniques is very low, whence non-uniform grids allow for large speed-ups and high speed-up can be expected from parallelization. The implementation is in progress, and results of experiments will be reported elsewhere

Berkvens, P.J.F.

Botchev, M.A.

English

Bochev, Mikhail A.

NARCIS 

Parallel processing and non-uniform grids in global air quality modeling

Bochev, M.A.

University of Twente Research Information

Parallel Processing and Non-Uniform Grids in Global AirQuality Modelling1P.J.F. Berkvens2 and M.A. Botchev3AbstractA large-scale global air quality model, running efficiently on a single vec-tor processor, is enhanced to make more realistic and more long-term simu-lations feasible. Two strategies are combined: non-uniform grids and parallelprocessing. The communication through the hierarchy of non-uniform gridsinterferes with the inter-processor communication. We discuss load balance inthe decomposition of the domain, I/O, and inter-processor communication. Amodel shows that the communication overhead for both techniques is very low,whence non-uniform grids allow for large speed-ups and high speed-up can beexpected from parallelization. The implementation is in progress, and resultsof experiments will be reported elsewhere.1 IntroductionWe consider computational aspects of global air quality models (AQMs). We de-scribe two methods for enhancing the efficiency of AQMs, namely two-way nesting(TWN) and domain decomposition with parallelization (DDP). We model the com-puter time and the wall clock time used by the computations and the communica-tions. We estimate the speed-ups that can be achieved by using TWN and DDP.Numerical methods for AQMs are typically based on a grid-like discretization ofthe atmosphere, but Lagrange-type methods exist. The species concentrations evo-lution on the grid is computed by a time-steppingmethod. In off-line models, whichwe consider here, the data for the meteorological and other mechanisms that influ-ence the concentrations are obtained during execution of the model program fromstored data sets or from a concurrently running dynamical meteorological model.Results like time series of concentrations, averages, and budgets (contribution toevolution per mechanism) are stored while the program is running.The costs of simulations with a comprehensive AQM on a fine grid are highcompared to the hardware available today. The number of grid cells is then at leastbased on cells of (longitude latitude) and layers. The numberof timesteps is based on a -year run with -hour timesteps. The numberof species is . The I/O is byte per real time second.These high costs call for very efficient computations. Apart from optimization ofthe numerical methods and the code, there are two additional techniques to enhance1This work has been done within the NWO project no. 613-302-0402CWI, P.O. Box 94079, 1090 GB Amsterdam, The Netherlands, e-mail: berkvens@cwi.nl3FOM-Institute for Plasma Physics, P.O. Box 1207, 3430 BE Nieuwegein, The Netherlands, e-mail:botchev@rijnh.nlBerkvens, P.J.F. and Botchev, M.A. (2002) Parallel Processing and Non-Uniform Grids in Global Air Quality Modeling. In: Proceedings of the Second Conference on Air Pollution Modelling and Simulation, APMS'01, April 9-12, 2001, Champs-sur-Marne, France. pp. 215-224. Springer Verlag. ISBN 3-540-42515-2 the efficiency of a method: two-way nesting (TWN) is meant to assign the availablecomputer time more efficiently, domain decomposition with parallelization (DDP)is meant to assign the available wall clock time more efficiently.A TWN technique has been implemented and tested [2], and combined with theupgrade TM5 of a particular AQM called TM3 [5, 9]. It has been applied to a casestudy [8]. The authors are currently implementing and testing DDP with TM3 andTM5. We present a model for the computer and wall-clock times for the calculationsand communications when running an AQM code equipped with TWN and DDP.Performance measurements are not available yet.The contents are as follows. Sect. 2 describes computational aspects of the nu-merical methods in an AQM, biased towards TM5. Sect. 3 discusses computationalaspects of our TWNmethod. Sect. 4 does the same for DDP. In Sect. 5 the computerand wall clock time speed-ups of a method are defined, and their enhancement byTWN and DDP is demonstrated theoretically. In Sect. 6 we give our conclusions.2 Numerical Methods for AQMThree-dimensional AQMs covering the whole atmosphere typically contain the fol-lowing mechanisms (with unit amounts of work):advection ( )turbulent mixing, cumulus convection ( )dry deposition, wet deposition ( )volume and surface emission ( )kinetic chemistry, photochemistry ( )For the unit amounts of work we refer to Table 1. We restrict ourselves to gas phasechemistry, but the computational aspects of aqeous and solid phase chemistry canbe modelled in a similar vein.Parameters quantifying these mechanisms have to be provided. They include(proportionality with numbers of grid cells) [number of inputs per month real time]:wind velocity ( ), surface pressure ( ), large scale precipitation rates( ), humidity ( ), temperature ( ) [ ];turbulent diffusivity ( ), up/downdraft en/detrainment rates ( ),bottom/top of convection layer ( ), convective precipitation rates ( ),cloud cover, cloud bottom/top ( ) [ ];dry deposition rates ( ), emission rates ( ), land use ( ) [ ].Results that are stored may include full running model status ( ) ( timesper month) for restarts, zonal ( ) and time ( ) averages of concentrations,and budgets ( ), all times per month real time. The above properties ofAQMs are partly taken from [7] and the current TM3 and TM5 code [9].Operator splitting makes the computations tractable by calculating the changesdue to the various mechanisms consecutively rather than simultaneously. Reductionof the degree of operator splitting is studied [1, 4] and gives better accuracy at thesame computational costs. Explicit time-stepping is used for advection, deposition,and emission, implicit time-stepping for the remaining mechanisms.We consider the case of a constant number of cells per grid column and of aconstant number of split steps (for the mechanisms) per timestep. Then the com-putational work load per grid column and per time-step is nearly constant. Themain differences are caused by the photochemistry which is active only in the sun-lit part of the atmosphere, and by differences in the vertical extent of the cumulusconvection. We treat this amount of work related to the time-stepping of themechanisms per column as a constant. In the nesting algorithm, data needs to becopied from parent to child and vice versa. The unit communication volume for acolumn to be updated by copying tracer data from or to a parent region is . Thecommunication volume (i.e. amount of data in bytes) of the I/O per column ofthe basic region is not constant per timestep but rather per real time second. Becausewe need their values later on, we will give rough estimates:flop, (1)flop, (2)byte. (3)We have used the data from Table 1, where the symbols are explained and rough es-timates for their values are given. The (insignificant) costs for coarsening input datafrom fine regions to coarse regions (to have the data on overlapping regions corre-spond) are not modelled explicitly. The cell numbers are for a regular grid withoutnesting. The machine characteristics are estimated from specifications provided bySARA, owner of the TERAS platform TM5 is implemented on. To simplify the cal-culations model below, , , and are assigned the same values. Thisdoes not affect the conclusions. For later use we define:(4)3 Two-Way NestingStarting from a regular three-dimensional grid with a basic resolution, TWN usespredefined fixed block-shaped nesting regions within the basic grid. They may beTable 1: Rough estimates for model quantitiesquantity value meaning360 cells in west-east direction180 cells in south-north direction60 cells in vertical direction30 cells in troposphere part of column64800 cells in the horizontal directions3888000 total number of cells50 transported species20 dry-deposited species20 wet-deposited species10 volumewise emitted species (e.g. aircraft)10 areawise emitted species (e.g. car traffic)100 chemically active species100 total number of species15 average number of reactions+coreactants per species8 number of mechanisms5 land use characteristics25 flop work for advection step15 flop work for vertical mixing proportional to5 flop work for vertical mixing proportional to1 flop work for vertical mixing proportional to5 flop work for dry or dry deposition5 flop work for emission100 flop work for (photo)chemistry1 flop work for budgets and averages8 byte number of bytes per real120/month large-scale meteorological inputs per month240/month subgrid scale meteorological inputs per month1/month surface and emission inputs per month1/month outputs per month360/month time-steps in basic region per month2592000 s seconds per 30-day month250Mflop/s average calculation speed250Mbyte/s communication speed among processors250Mbyte/s I/O speed15 s latency time in interprocessor communicationrecursively nested. Each nesting child region is spatially and/or temporally refinedwith respect to its parent region, which contains it. Because local refinement in thevertical direction is less important, we only consider horizontal spatial refinement.A parent and a child are coupled by advection through their interface. The concen-trations in the child’s interior evolve autonomously. On the interfaces the advectionflux is determined from data on both sides, so there is a two-way coupling.We give an example for a time refinement factor of to illustrate how the nestingalgorithm works, but it works for any refinement factor. Let , , , denote-advection steps and a step of the remaining mechanisms in the parent region,and similarly , , , in the child region. Then a symmetric algorithm reads:(5)corresponding to one time-step in the parent region and two in the child region.Before each and step data is copied from the parent to the interface of its child.After a time-step in the child has been finished, its data is copied back to its parent.Our nesting algorithm conserves mass, positivity, and homogeneity (in terms of themixing ratio) of the solution if the underlying advection scheme does so [3, 2], whileallowing symmetric operator splitting in each region (except for the interface cells).When and consist of substeps for various mechanisms, their order should bereversed in their second occurrence in a timestep to maintain symmetric splitting.We now model the computational work for the TWN method. Letbe the number of nesting regions ( indicating the basic region). Let, , , be the spatial dimensions and the time interval, per nestingregion, and define . ( and are the same for each region.)Let , , , be the spatial and temporal resolutions per nest, and define. ( is the same in each region.) We model as:(6)For the I/O volume we model:(7)Assuming region has the largest refinement, we model the computational workwith the original method without TWN, that is needed to get at least the sameresolution everywhere as with TWN, as:(8)4 Domain Decomposition with ParallelizationIn DDP all (nesting) regions are divided into as many subdomains as available paral-lel processors. The computational work is distributed over the processors as evenlyas possible. For the parallel processing we use an SPMD approach with MPI. Extraoverhead is caused by the communication among the processors, and by the fact thatin MPI only one processor is guaranteed to be available for I/O. For long enoughcomputations, we expect that I/O is the largest intrinsically serial part of the compu-tations. Other serial parts include initializations (‘bookkeeping’ of data structures)and ‘finalizations’.If processors are assigned, the basic region and each of the zoomregions are decomposed as evenly as possible into parts in the west-east directionand into parts in the south-north direction, such that a Cartesian grid of rectangu-lar subdomains results in each region. Assuming a corresponding virtual Cartesiantopology of processors, each processor is assigned the corresponding sub-domain of each region. Since the amount of work per grid column and per timestep is fixed (Sect. 2) and because no self-adjusting time-steps are used, we have theadvantage that only this static load balancing is needed.Parallel computations need extra computer time due to inter-processor commu-nication. It consists of three parts, namely intra-region communication, inter-regioncommunication, and scatter and gather of I/O. Intra-region communication is asso-ciated with the advection calculations in a single region, where each processor needsdata from neighbouring subdomains, residing on other processors, to determine thefluxes on and near its boundaries. Inter-region communication is needed for advec-tion calculations in multiple regions, where for a particular region a processor mayneed data from its parent or child to compute fluxes at or near the interfaces. I/O-related communication consists of the scattering of input data from processor tothe other processors and of the gathering of output data back to processor .We approximate the communication volumes as follows. The intra-region com-munication volume for exchanges between neighbouring subdomains, the inter-region communication volume for updates, and the communication volumefor the scatter and gather of the I/O are modelled as:(9)(10)(11)where we have neglected the fact that on the boundaries of a (nesting) region, intra-region communication may be absent. Such communication procedures have to beset up , , and times respectively. Let the factor account for the factthat for any region its overlap in the parent may be distributed over more thanprocessor. Each child has only one parent, which is larger than or at most equal insize to the child, so each particular child subdomain has , , or ‘parent’ processorsto communicate with. We then have:(12)(13)(14)The total parallelization communication volume and the total number ofcommunication set-ups are:(15)5 Speed-UpWe determine the speed up due to nesting ( ) and due to parallelization ( ).5.1 NestingWe define the speed-up of the nesting method as the ratio of computer timeconsumed for a particular choice of nesting regions and the amount of computertime needed when no nesting is used. Writing for the calculation speed (floprate) and for I/O speed (byte rate) of the actual computer, we express:(16)In actual computations there is overhead, which we have modelled with an overheadfactor depending on the coding, the compiler, and the computer.Define and . Substituting equations (6) to(8) and flop byte byte in equation (16) yields:(17)where . We use the estimates of , , and to concludeand . If then:(18)5.2 ParallelizationWe define the speed-up of the parallelized method as the ratio of wall clock timesspent when multiple processors (and subdomains) are used and when a single pro-cessor is used. Let be the scalar fraction of the calculation work withTWN on processor, see equation (16). Introducing for the communicationspeed (byte rate) of the actual computer, we define as:(19)Here is the latency time for a communication set-up, and is an organi-zational overhead factor related to the division of work into smaller amounts.Substituting from the equations in Sects. 3 and 4, and dividing the numeratorand the denominator by , we find for :(20)where use was made of flop byte byte and where(21)(22)(23)(24)(25)(26)Since and (a region is not divided into nearly asmany subdomains as cells), we have . The terms in equation (20) which‘threaten’ the scalability most, are and, apart from the factor . ,being the numbers of cells in the horizontal directions in region , we have:(27)The first quotient in the above equation is no larger than , so (e.g.)is certainly guaranteed if:(28)For realistic values (e.g. , , and ) this inequality holdsfor realistic numbers of processors in region . The second ‘threatening’ term canonly become near or larger than (and therefore dangerous) if not, which only happens for unrealistically coarse resolutions in the basic region.Using the smallness of many terms we approximate equation (20) as follows:(29)Defining as the scalar part of the computer time spent on actual computationsby processors [6], then , and(30)We expect to increase with the number of processors , but we assumefor large enough problems. Then the speed-up is almost proportionalto , depending on the smallness of .5.3 ExampleWe give a realistic example of computations with an AQM and the expected speed-ups from TWN and DDP. First consider TWN with three regions: region containsregion , which in turn fully contains region . For the dimensions of the regionswe specify , , and , whereas, , and .With these numbers, we find , which is a considerable speed-up factor.Next we consider DDP with , , and . Since the aboveresolutions and nesting are nowhere near the danger zone for (almost) perfect speed-up, we find which clearly favours DDP.6 ConclusionWe propose two techniques which make atmospheric chemistry simulations farcheaper. The computer time can be strongly reduced by Two-Way Nesting. Thewall clock time is expected to be reduced according to almost ideal speed-up forlarge problems with Domain Decomposition with Parallelization. This brings morerealistic simulations within reach. Developments towards this goal are underway.References[1] P. J. F. Berkvens, M. A. Botchev, M. C. Krol, W. Peters, and J. G. Verwer. Solv-ing vertical transport and chemistry in air pollution models. Technical ReportMAS-R0023, CWI (Centre for Mathematics and Computer Science), Amster-dam, The Netherlands, Aug 2000.[2] P. J. F. Berkvens, M. A. Botchev, W. M. Lioen, and J. G. Verwer. A zoom-ing technique for wind transport of air pollution. In R. Vilsmeier, D. Hänel,and F. Benkhaldoun, editors, Finite Volumes for Complex Applications. HermesScience Publications, 1999. (Abridged from [3]).[3] P. J. F. Berkvens, M. A. Botchev, W. M. Lioen, and J. G. Verwer. A zoom-ing technique for wind transport of air pollution. Technical Report MAS-R9921, CWI (Centre for Mathematics and Computer Science), Amsterdam, TheNetherlands, Aug 1999.[4] P. J. F. Berkvens, M. A. Botchev, and J. G. Verwer. On the efficient treatmentof vertical mixing and chemistry in air pollution modelling. In M. Deville andR. Owens, editors, CD-ROM Proceedings of the 16th IMACS World CongressOn Scientific Computation, Applied Mathematics and Simulation, 2000. EcolePolytechnique Fédérale de Lausanne, Lausanne, Switzerland, Aug 21-25, 2000.[5] F. J. Dentener, J. Feichter, and A. Jeuken. Simulation of transport of Rnusing on-line and off-line global models at different horizontal resolutions: Adetailed comparison with measurements. Tel. Ser. B, 51:573–602, 1999.[6] D. R. Emerson. Introduction to parallel computers: architecture and algorithms.In P. Wesseling, editor, High Performance Computing in Fluid Dynamics, vol-ume 3 of ERCOFTAC Series, chapter 1, pages 1–42. Kluwer Academic Pub-lishers, Dordrecht, Boston, London, 1996.[7] M. Heimann. The global atmospheric tracer model TM2. Technical report,DKRZ (Deutsches Klimarechenzentrum), Hamburg, Germany, 1995.[8] M. C. Krol, W. Peters, M. A. Botchev, and P. J. F. Berkvens. A new algorithm fortwo-way nesting in nlobal models: principles and applications. In B. Sportisse,editor, Proceedings of the Second Conference on Air Pollution Modelling andSimulation 2001, Champs-sur-Marne, France. Springer Verlag, to appear in2001. submitted.[9] Url http://www.phys.uu.nl/ peters/TM3/TM3S.html.

Parallel Processing and Non-uniform Grids in Global Air Quality Modelling

P. J. F. Berkvens

M. A. Botchev

Crossref

A zooming technique for wind transport of air pollution.

Anewalgorithmfor two-way nesting in nlobal models: principles and applications.

Introductionto parallelcomputers: architectureand algorithms.

On the efﬁcient treatment of vertical mixing and chemistry in air pollution modelling.

Simulation of transport of Rn using on-line and off-line global models at different horizontal resolutions: A detailed comparison with measurements.

Solving vertical transport and chemistry in air pollution models.

The global atmospheric tracer model TM2.

Parallel processing and non-uniform grids in global air quality modeling

Abstract

Similar works

Full text

Available Versions

NARCIS

University of Twente Research Information

Crossref