Geoscientific Model Development the 1-way On-line Coupled Atmospheric Chemistry Model System Meco(n) – Part 1: Description of the Limited-area Atmospheric Chemistry Model Cosmo/messy

The numerical weather prediction model of the Consortium for Small Scale Modelling (COSMO), maintained by the German weather service (DWD), is connected with the Modular Earth Submodel System (MESSy). This effort is undertaken in preparation of a new, limited-area atmospheric chemistry model. Limited-area models require lateral boundary conditions for all prognostic variables. Therefore the quality of a regional chemistry model is expected to improve, if boundary conditions for the chemical constituents are provided by the driving model in consistence with the meteorological boundary conditions. The new developed model is as consistent as possible, with respect to atmospheric chemistry and related processes, with a previously developed global atmospheric chemistry general circulation model: the ECHAM/MESSy Atmospheric Chemistry (EMAC) model. The combined system constitutes a new research tool, bridging the global to the meso-γ scale for atmospheric chemistry research. MESSy provides the infrastructure and includes, among others, the process and diagnostic submodels for atmospheric chemistry simulations. Furthermore , MESSy is highly flexible allowing model setups with tailor made complexity, depending on the scientific question. Here, the connection of the MESSy infrastructure to the COSMO model is documented and also the code changes required for the generalisation of regular MESSy submod-els. Moreover, previously published prototype submodels for simplified tracer studies are generalised to be plugged-in and used in the global and the limited-area model. They are used to evaluate the TRACER interface implementation in the new COSMO/MESSy model system and the tracer transport characteristics, an important prerequisite for future atmospheric chemistry applications. A supplementary document with further details on the technical implementation of the MESSy interface into COSMO with a complete list of modifications to the COSMO code is provided.


Introduction
Since atmospheric chemistry related processes are often governed by local features, e.g.emissions of a power plant or fire emissions, which are not sufficiently resolved in global models, smaller scale models for atmospheric chemistry are required.In order to investigate those smaller scale effects, a regional model including a chemistry setup consistent with that of the global model is desirable.This goal is achieved by connecting MESSy to the regional weather prediction and climate model of the Consortium for Small Scale Modelling (COSMO model1 , Steppeler et al., 2003;Doms and Schättler, 1999), resulting in the model system COSMO/MESSy.
Here, an overview of the functionality of COSMO/MESSy is provided.First, the COSMO model and the MESSy interface structure are briefly recapitulated in Sects.2.1 and 2.2, respectively.Section 3 describes the connection of the MESSy infrastructure submodels to COSMO and Sect. 4 depicts the required generalisation of the regular submodels.
A first application of some previously published submodels3 focusing on the evaluation of the functionality of the TRACER interface in COSMO/MESSy and the tracer transport characteristics is presented in Sect.5, before the potential of the newly built model system is shortly summarised in Sect.6.This is the first part of three companion articles.It documents the implementation of the MESSy infrastructure (plus a few diagnostic submodels to test mainly the TRACER infrastructure) into the COSMO model.The resulting model is a regional model of the atmosphere, enabled for (now standardised) further extensions into a regional atmospheric chemistry model.The second part (Kerkweg and Jöckel, 2012) is about a different issue, however, requires the COSMO/MESSy model of the first part: there, we present a newly developed coupling technique that provides the boundary data required by the regional model directly (on-line) from ECHAM5/MESSy via the newly developed Multi-Model-Driver (MMD) library and two corresponding submodels.This on-line coupled system is called MECO(n), i.e.MESSy-fied ECHAM and COSMO models nested n times.The third article (Hofmann et al., 2012) provides a meteorological evaluation of the nested system, described technically in Part 2. This evaluation, focusing on distinct meteorological events on synoptic scale, and on the question if and how they can be reproduced by MECO(n), is a prerequisite for further applications with chemistry, like chemical weather (air pollution) forecasts, measurement campaign analyses etc.

Model description
The model system COSMO/MESSy consists of two model components: the numerical weather prediction model of the Consortium for Small Scale Modelling (COSMO model, Doms and Schättler, 1999) and the Modular Earth Submodel System (MESSy, Jöckel et al., 2005).Both components are briefly introduced.Currently, the model system is based on the COSMO model version cosmo_4.8_clm12 and the MESSy version 2.41.

The COSMO model
The COSMO model is a non-hydrostatic limited-area weather prediction model (Steppeler et al., 2003;Doms and Schättler, 1999).It has been designed for both operational numerical weather prediction (NWP) and for research on the meso-β and meso-γ scale.Additionally to the short range NWP application, the COSMO model was set up for regional climate studies by the CLM-community4 (Rockel et al., 2008).The COSMO model is based on the primitive thermo-hydrodynamical equations for compressible flow in a moist atmosphere and uses rotated geographical coordinates and a generalised terrain following height coordinate.The Arakawa C-grid with Lorenz vertical grid staggering is used.The COSMO model treats the following prognostic variables: horizontal and vertical Cartesian wind components, pressure perturbation, temperature, specific humidity and cloud water content.There are additional optional prognostic variables such as the cloud ice content, the turbulent kinetic energy, the specific water contents of rain, snow and graupel.Different time integration schemes can be selected: -a second-order leapfrog, horizontally explicit, vertically implicit, time-split integration scheme, including extensions proposed by Skamarock and Klemp (1992), -a three time-level 3-D semi-implicit scheme (Thomas et al., 2000), -several options for two time-level 2nd and 3rd order Runge-Kutta split-explicit schemes after Wicker and Skamarock (2002) and -a TVD-variant (Total Variation Diminishing) of a 3rd order Runge-Kutta split-explicit scheme.
A variety of physical processes are taken into account through parameterisations, e.g.grid-scale and subgrid-scale clouds, precipitation, moist and shallow convection, radiation, a soil model and so forth.
The COSMO model is available to the scientific community 5 .To make the COSMO model applicable on the climate scale, a number of extensions have been implemented by the CLM-community, for instance, netCDF6 output and a restart7 facility.Hereafter, the COSMO model is named COSMO-CLM when features are discussed which are only used by the climate community.

MESSy
The Modular Earth Submodel System (MESSy, Jöckel et al., 2005Jöckel et al., , 2010) is a multi-institutional project providing the infrastructure to expand state-of-the-art geoscientific domain models (e.g. of the atmosphere) into comprehensive Earth System Models (ESMs).Although originally the focus was solely on the efficient and flexible implementation of processes related to atmospheric chemistry into a circulation model of the atmosphere, the methodology turned out to be much more powerful (e.g.Pozzer et al., 2011).In particular, MESSy has been connected to the ECHAM5 general circulation model of the atmosphere (Roeckner et al., 2006), expanding it into the global chemistry climate model EMAC (Jöckel et al., 2006) with a wide variety of applications (see also Sect.1).
The basic idea behind MESSy is its four layer structure, as visualised in Fig. 1, consisting of -the basemodel layer (BML), usually a legacy domain model, e.g. of the atmosphere; it can be regarded as power supply, -the basemodel interface layer (BMIL), which hosts the MESSy infrastructure and can be regarded as multiple socket outlet, -the submodel interface layer (SMIL), which provides the communication of the submodels with the MESSy infrastructure (and therefore also with the basemodel); it can be regarded as the plug, -the submodel core layer (SMCL), which hosts the actual process or diagnostic formulation; it can be regarded as the machinery.
Each MESSy submodel is split into an interface and a core layer part.The MESSy infrastructure is coded as so-called generic submodels, which interfaces reside in the BMIL, with COSMO members.However, all users are required to sign an agreement with a COSMO national meteorological service and to respect certain conditions and restrictions on code usage.For questions concerning the request and the agreement, please contact the chairman of the COSMO Steering Committee.In the case of a planned operational or commercial use of the COSMO-Model package, special regulations will apply" (cited from the COSMO User's Guide, December 2009 http://www.cosmo-model.org/content/model/documentation/core/cosmoUserGuide.pdf). 6http://www.unidata.ucar.edu/software/netcdf/ 7Appendix A contains a Glossary explaining some terms repeatedly used here.The terms from the Glossary are written in italics throughout the article.
A. Kerkweg and P. Jöckel: COSMO/MESSy mate scale, a number of extensions have been implemented by the CLM-community, for instance, netCDF 6 output and a restart 7 facility.Hereafter, the COSMO model is named COSMO-CLM when features are discussed which are only used by the climate community.

MESSy
The Modular Earth Submodel System (MESSy, Jöckel et al., 2005Jöckel et al., , 2010) is a multi-institutional project providing the infrastructure to expand state-of-the-art geoscientific domain models (e.g., of the atmosphere) into comprehensive Earth System Models (ESMs).Although originally the focus was solely on the efficient and flexible implementation of processes related to atmospheric chemistry into a circulation model of the atmosphere, the methodology turned out to be much more powerful (e.g., Pozzer et al., 2011).In particular, MESSy has been connected to the ECHAM5 general circulation model of the atmosphere (Roeckner et al., 2006), expanding it into the global chemistry climate model EMAC (Jöckel et al., 2006) with a wide variety of applications (see also Sect.1).
The basic idea behind MESSy is its four layer structure, as visualised in Fig. 1  whereas the interfaces for the regular submodels belong to the SMIL.The core, i.e. the basemodel independent parts of all submodels, is located in the SMCL.Examples of generic submodels are (see also Sect.3), among others, the memory management and input/output interface, including a restart facility (CHANNEL, Jöckel et al., 2010), the tracer infrastructure for the handling of constituents (TRACER, Jöckel et al., 2008), or the time and event management (TIMER, Jöckel et al., 2010).For a more detailed overview of the currently available generic and regular submodels, we refer to Jöckel et al. (2010) and the project web-page8 .The next section (Sect.3) describes the implementation details of the connection of the MESSy infrastructure to the COSMO model.Section 4 documents the required generalisation (i.e.including an adapter into the plug) of the SMIL of the regular submodels.The necessity for such a generalisation arose mainly from the fact that in the Eulerian (grid-point) formulation of many processes the vertical direction needs to be particularly distinguished (in other words plays a particular role) from the horizontal directions.As a consequence, the rank of the vertical index in multi-dimensional array variables representing geo-located data might depend on the basemodel (and indeed it is different between ECHAM5 and COSMO).This dependency has been bypassed.

Implementation of the MESSy infrastructure
An important principle for the implementation of the MESSy interface in the COSMO model is to keep the changes in the COSMO model code as small as possible, i.e. minimally invasive.Additionally, nothing in the functionality of the COSMO model should be changed for the usual "COSMO only" user.Hence, all changes have been introduced using the pre-processor directive MESSY: Thus, the changes are only active if the model is configured with --enable-MESSY (default).Otherwise, with --disable-MESSY, the original COSMO model code is compiled.
One of the basic ideas of MESSy is to keep as much code as possible applicable for all models employed.In this way, code doubling is avoided, making it highly consistent and less error-prone.This is not always easily achieved, as for 3-D applications MESSy currently still depends on a basemodel providing the model domain decomposition, the dynamics and other distinct parts of the model physics.The global climate model ECHAM5 (Roeckner et al., 2006) was the basemodel chosen as basis for the development of MESSy (Jöckel et al., 2006).Before COSMO/MESSy was built, ECHAM5 was the only 3-D basemodel hosting the MESSy interface.Therefore, MESSy uses mostly the naming convention of ECHAM5.For instance, the names of variables used in the MESSy SMIL are in many cases the same as in ECHAM5.As long as MESSy was connected to one 3-D basemodel only, there was no need to take different realisations of 3-D model grids and/or representations into account.During the course of connecting MESSy to the COSMO model, MESSy was further generalised to allow for different realisations of model grids and representations.Note that all the changes apply to the interface layers only, as the core layer of the submodels is basemodel independent anyway.
The three most prominent extensions are: 1.The order of the dimensions in multi-dimensional array variables representing geo-located information: In principle, the three dimensions in space can be arbitrarily defined in a three dimensional model.But, in meteorological models the vertical is a distinguished spatial dimension.Therefore the position of the vertical axis becomes important when using three dimensional fields.In ECHAM5 the order of the three dimensions in space is (h1,v,h2), with h1,h2 the two horizontal dimensions and v the vertical dimension, whereas the order in COSMO is (h1,h2,v).Before the connection of MESSy to the COSMO model, the MESSy interfaces of the submodels worked under the assumption, that the order of dimensions is (h1,v,h2).This restriction has been eliminated now, the order of dimensions is automatically switched for different basemodels.This has been implemented using pre-processor directives, as explained in Sect. 4.

A big outer loop over the second horizontal dimension:
The specific order of dimensions in ECHAM5, (h1,v,h2), has been introduced for code optimisation.Within a loop over the second horizontal dimension (outer loop) arrays of rank 2 (h1,v) are processed.currently available generic submodels is provided in the "COSMO/MESSy Implementation Documentation" in the Supplement.Note that MESSy version 2 as documented by Jöckel et al. (2010) includes already the basis for the coupling of MESSy to COSMO.Therefore, according to the MESSy concept, the further developments discussed here, either do not interfere with ECHAM5/MESSy, because they are not relevant for it, or they are immediately applicable also for ECHAM5/MESSy.

SWITCH/CONTROL: switches and main entry points for individual submodels
SWITCH/CONTROL is the MESSy submodel managing the switching and calling of the individual submodels.In the basemodels, distinct entry points have been defined from which the MESSy submodels are called.One entry point consists of one call to a subroutine of CONTROL.This subroutine executes calls to the generic submodels and, if switched on, to the respective regular submodels.Figure 2 lists all current entry points for MESSy.
As not all SMILs of all submodels have been adapted so far, not all entry points available in ECHAM5/MESSy are also used in COSMO/MESSy.Those marked by hatched boxes are presently not used in COSMO/MESSy, the yellow boxes, in contrast, denote the active entry points in COSMO/MESSy.The MESSy submodels currently called from the respective subroutine are listed next to the box.Blue text highlights the generic submodels, whereas the regular submodels are written in black.In those entry points located within the outer loop over the second horizontal dimension in EMAC (indicated by the green box), i.e. messy_vdiff, messy_convec, messy_physc and messy_local_end, the loop is mimicked with loop variable jrow, in agreement with the naming convention.At first place in these loops, some rank-2 POINTERs are associated with the respective 2-D sub-arrays of 3-D fields by calling the subroutine main_data_2D_set_jrow:

CHANNEL: the memory management, output and restart control
The generic MESSy submodel CHANNEL is described in detail by Jöckel et al. (2010).It manages the memory usage, the output, the restart files and the model internal access to data.Internally, in CHANNEL all data are stored as concatenated list of structures, each containing a description of the data (meta information) and a rank-4 POINTER for the data  To make the COSMO model data fields available to all MESSy submodels, they have to be defined as channel objects.As one prerequisite the TARGET, ALLOCATABLE attribute of the COSMO arrays are replaced by the POINTER attribute.
The COSMO model arrays are allocated by creating the respective channel objects (in a channel called COSMO_ORI) in the subroutine messy_COSMO_create_channel instead of being allocated in the COSMO subroutine alloc_meteofields (compare flow chart in Fig. 3).For instance, the 3-D field of the density of the reference atmosphere (rho0) is defined by RHO0 is the name of the channel object.The memory space required for the rank-3 POINTER rho0 is indicated by the representation ID GP_3D_MID, which determines that the variable is defined at grid mid-points.The CHARAC-TER formal parameter callstr denotes the calling subroutine.Detailed information is provided in the CHANNEL user manual (available in the supplement of Jöckel et al., 2010).The channel COSMO_ORI is used for data management only, i.e. its contents should normally not be subject to data output (but the objects can be written, if required for debugging purposes).
CHANNEL provides its own output control.Within the channel namelist channel.nmldefaults for each channel and channel object can be set to create output e.g. for instantaneous values, averages, standard deviations, minimum and maximum values etc. (see CHANNEL manual, supplement of Jöckel et al., 2010).To keep the functionality of the &GRIBOUT namelists of the COSMO namelist file INPUT_IO, the requested output is redirected into MESSy channels.The channels are named COSMOXXXy with XXX the number of the &GRIBOUT namelist, e.g. if three &GRIBOUT namelists are specified in the INPUT_IO namelist file, XXX in the respective channels is replaced by the numbers 001, 002 and 003, respectively; y indicates the respective variable group and is defined as in COSMO: c for constants, m for variables on the model grid and p and z for output interpolated on pressure or altitude levels, respectively.In that way, the MESSy CHANNEL output of COSMO model variables can still be managed by the original COSMO model namelists, if desired.
In addition, the original COSMO output files can be written, if the LOGICAL switch L_BM_ORIG_OUPUT in the CHANNEL namelist is set .TRUE.9 .
The generic submodel CHANNEL also manages the model restarts.For a model restart all variables, required for the unambiguous continuation of the simulation, are dumped into so-called restart-files.Unambiguous in this context means that the results are binary identical compared with results from a continuous simulation without interruption.The restart-files are read during the initialisation phase of a restarted simulation.The climate version of the COSMO-CLM model already included a restart facility.Nevertheless, as CHANNEL replaces the complete COSMO model memory management, the restarts are also managed by CHANNEL.Thus, the read procedure of the restart file in the COSMO subroutine organize_input (src_input.f90) is skipped in COSMO/MESSy.

TIMER: the "heart-beat" and event management
The generic submodel TIMER became part of the MESSy system within the scope of the connection of MESSy to the COSMO model.Until now, ECHAM5/MESSy used the time and event management of ECHAM5 and all MESSy submodels utilising events depended on the ECHAM5 event management routines.Therefore, the functionality of ECHAM5 events had to be made available within COSMO/MESSy as well.Consequently, the generic MESSy submodel TIMER is based on the time and event managing routines of ECHAM5 written by Ingo Kirchner (Max-Planck-Institut for Meteorology, Hamburg, now at FU Berlin).A user manual for the TIMER submodel can be found in the supplement of Jöckel et al. (2010).
As MESSy provides the time management and because each model can apply only one time management, the COSMO model time management is overwritten by the MESSy generic submodel TIMER.Consequently, all time variables as time step, start date, stop date and the restart frequency are defined in the namelist file of TIMER (timer.nml), the entries in the COSMO namelist &RUNCTL in the namelist file INPUT_ORG for dt, hstart, hstop and ydate_ini are hence ignored, as is the entry restart_hour of the namelist INPUT_IO.
The most important time variables are initialised via namelist in messy_timer_setup, which is called from the organize_setup subroutine of the COSMO model via the CONTROL subroutine messy_setup (see Fig. 3 Additionally, the counters for the re-initialisation of the maximum 10 m wind velocity and the minimum and maximum temperature are re-set here.After this re-initialisation, all COSMO time variables are consistent with the TIMER setup.
Because the original COSMO-CLM restart facility is replaced by the CHANNEL restart interface in COSMO/MESSy, as mentioned above, the subroutine messy_timer_COSMO_reinit_time is additionally called during a restart after reading the restart attributes from the restart-files.
The generic submodel TIMER is described in more detail by Jöckel et al. (2010) and the corresponding Supplement.

DATA: the data transfer interface
The generic submodel DATA is part of the memory management and data transfer interface.The basemodel data fields are made available and easily accessible by DATA.One of the guidelines for the MESSy interface implementation into the COSMO model is to use as many parts as possible available for ECHAM5/MESSy also for COSMO/MESSy, i.e. to minimise the code that needs to be maintained independently for ECHAM5/MESSy and COSMO/MESSy.Therefore, DATA is the model unifying the data structure and names of the basemodel variables.Different aspects had to be taken into account: -When "US(E)ing" a (basemodel) variable provided by DATA, the MESSy submodels access this variable always with the same variable name independent of the basemodel.For instance, the variable name of the sea-land fraction in ECHAM5/MESSy is slf, whereas the same variable in the COSMO model is fr_land.stant).Examples for such variables are the geopotential at full and interface levels (geopot and geopoti, respectively), the mass contained in a grid box (grmass), or the volume of a grid box (grvol).

MPI: a high-level interface to the MPI library
As MESSy must be applicable in a parallel decomposed environment, the generic MESSy submodel MPI builds a high level application interface (API) for the use of the MPI library.Currently, for each basemodel the access to the MPI library is mirrored in the MESSy submodel MPI.In case of ECHAM5/MESSy, the ECHAM5 routines are simply USEd into the generic submodel MPI.Consequently, for the COSMO model, two aspects had to be taken into account: -Each basemodel uses the MPI library in its own way, thus, COSMO/MESSy uses the MPI library for the MESSy submodels similar to the COSMO model.
-The names of the high level routines organising the data exchange via the MPI library in the generic submodel MPI must be the same for different basemodels.As the MESSy submodels call the MPI routines (from their SMIL) using the original ECHAM5 names, these names must also be used for COSMO/MESSy.
To achieve this goal, the generic MESSy submodel MPI contains subroutines named as in ECHAM5/MESSy and structured, where necessary, either as the ECHAM5/MESSy or the COSMO model routines, respectively.In particular it contains the following subroutines/interfaces: -gather_gp: this is a threefold overloaded interface for gathering 4-D-, 3-D-or 2-D fields in grid-point space.The subroutines for 3-D-or 4-D fields reduce the field by one rank and call the next "lower" subroutine, thus gather_gp basically performs an exchange of parallel decomposed horizontal 2-D fields.In the COSMO model this is accomplished by the subroutine gather_field.Consequently, this subroutine is called from the subroutine gather_gp to actually perform the data transfer.
-scatter_gp: this subroutine provides the reverse action to gather_gp.It distributes fields to a parallel decomposed grid.Equivalently overloaded, scatter_gp basically calls the COSMO model subroutine distribute_field.
-p_bcast: the subroutine p_bcast is twelvefold overloaded for the transfer of 0D -messy_mpi_initialize: this MESSy subroutine is called directly from the basemodel to initialise the LOGICAL variables p_parallel_io and p_parallel indicating whether the actual PE is responsible for input and output, and whether the model is working in a parallel environment, respectively.
-p_abort: this subroutine manages the proper termination of a simulation in a parallel environment by calling the COSMO model routine model_abort.
Some further variables used in the MESSy submodels are defined accordingly.The "COSMO/MESSy Implementation Documentation" in the Supplement contains a description of these variables.

TRANSFORM: the interface for grid transformations
The generic MESSy submodel TRANSFORM hosts utility subroutines to transform variable fields from one representation into another (e.g.grid point, spectral or Lagrangian).
Most of these representations are meaningless in a regional grid point model.This leaves only one important subroutine, namely locate_in_decomp.This subroutine determines the corresponding horizontal indices and PE number in a parallel decomposed grid for a point given in geographical coordinates, or, for COSMO/MESSy optionally in coordinates of the rotated grid.The subroutine locate_in_decomp is for instance required in the regular MESSy submodel TR-EXP (see Sect. 5.1.2 and Jöckel et al., 2010) to correctly identify the geo-location of point sources.

TRACER: the management of constituents
The generic submodel TRACER (Jöckel et al., 2008) provides the interface for the management of constituents (e.g. chemical species) in the model system.
The tracer definition consists of two parts: -the meta-information defining the properties of the tracer, e.g.
-the quantity (amount-fraction, number-density or concentration) and the unit of the component, -the medium in which the tracer resides (air, aerosol, cloud, ocean, etc.), -switches for the processes the tracer should undergo within the simulation (e.g.advection, dry deposition, scavenging, etc.) and -some tracer specific constants as the molar mass, the henry coefficient, etc.
-The tracer data field itself.It is always defined for one tracer set, i.e. for a group of tracers with the same representation i.e. the same geometry in space and time.
From the rank-6 tracer field -rank 1, 2 and 4 span the spatial dimensions, -rank 3 is the tracer index, i.e. the number of an individual tracer in the tracer set, -rank 5 is of length 1 and -rank 6 contains the data instances, e.g. the different time levels of a time integration scheme.
For a detailed overview of the submodel TRACER, we refer to Jöckel et al. (2008).Here, we provide details about the implementation of the TRACER interface into the COSMO model.
One of the striking differences between the COSMO model and ECHAM5 is the treatment of the prognostic variables (see Sect. 3.4).In the COSMO model prognostic variables are allocated with an extra rank for the different time levels required by the time integration scheme.The indices of the respective time levels are nnow, nnew and nold, which are rotated each model time step in order to avoid the copying of the data from the "nnew" to the "nnow" time level and so forth.Therefore the data for one specific time level is not always located at the same memory space.This is different in the MESSy.Due to the association of the tracers to the channel memory management, each of the instances (rank 6) of the tracer field is permanently associated to one specific time level of the time integration scheme.Table 1 shows the definition of the instances for COSMO/MESSy as currently chosen10 .The number of instances depends on the time integration scheme: for the Runge-Kutta scheme only 5 instances are needed, whereas the leap-frog scheme requires 6 instances.The different time levels of the tracer fields are copied at the beginning of the new time step, i.e. in the subroutine initialize_loop in lmorg.f90.
Consequently, the access to the respective time levels of the tracer field is not as simple as in the COSMO model itself.To mimic the functionality of the access to the different time levels of a prognostic field (by the indices nnow, nnew and nold), additionally the POINTER ARRAY xt_array is allocated to the number of time levels in the integration scheme.In addition to copying the tracer time instances in initialize_loop, the individual POINTERs of xt_array are associated to the respective tracer fields: The transport of the tracers via advection and horizontal and vertical diffusion is included into the COSMO model in equivalence to the transport of water vapour (qv).Due to the meta-information associated to each tracer, each transport process can independently be switched on or off for each individual tracer.Further information are provided within the comprehensive description of the generic MESSy submodel TRACER by Jöckel et al. (2008).

Implementation of the MESSy submodels
A MESSy submodel consists of two layers (Jöckel et al., 2005, see also Sect.2.2): -The submodel core layer (SMCL) hosts the code parts describing the physical, chemical or diagnostic process dealt with by the submodel.Because the SMCL contains only the scientific content of a submodel, it is completely independent of the basemodel (global, regional, column or box model) and therefore stays unchanged for the implementation of MESSy into the COSMO model.-The submodel interface layer (SMIL) manages the communication (data flow) between the basemodel interface layer (i.e. the MESSy infrastructure) and the submodel core layer.It allocates the memory for the submodel specific fields, and organises the access to the fields defined by other MESSy submodels or the basemodel, which are required in the respective submodel core.Thus, a SMIL file provides the connector to the basemodel.

Geosci
So far, for the implementation of the MESSy regular submodels the SMIL files of each submodel had to be changed for different basemodels.Within the scope of the connection of MESSy to the COSMO model, the SMIL files have been generalised in a way that for additional new basemodels no further changes of the SMIL files are required.Most of the differences between the basemodels (here the COSMO model and ECHAM5) is accounted for by the generic MESSy submodels, especially by DATA (Sect.3.4).Only one big difference remains for all regular submodels, which has to be taken care of in the SMIL files: the order of the spatial dimensions.Whereas ECHAM5/MESSy arranges the spatial dimensions of a grid-point field (h1,z,h2), with h1,h2 the horizontal dimensions and z the vertical dimension, the COSMO model applies the order (h1,h2,z).This difference in order needs to be taken into account in the SMIL, when data fields are provided to the core layer subroutines, where the vertical dimension is usually particularly distinguished.This is implemented by applying a rankflipping based on pre-processor directives.This approach proves useful to avoid unnecessary doubling of code, implying a lower error-proneness and a reduced risk of inconsistencies, as changes in the submodel interfaces are directly valid for all grid-point basemodels.
The rank-flipping is implemented via pre-processor directives, the so-called "rank identifiers".Basically, the ranks that need to be flipped (dependent on the basemodel) are replaced by a directive in the code.As a lot of different combinations of indices and colons is possible, a naming convention for the Fortran95 index variables and the rank identifiers is required.All rank identifiers are defined basemodel dependent in the include file messy_main_ppd.inc,which is included in all SMIL files.For instance, with rank identifiers the calculation of the thickness of the lowest model layer in pressure units (dp) reads: pressi_3d is the pressure at the layer interfaces and kproma the length of the first (blocked) horizontal dimension.The rank identifiers _RI_YV_ and _RI_YVp1_ are defined as where nlev is the number of vertical levels and jrow is the loop variable of the second horizontal dimension, indicated by the V and the Y in the rank identifier, respectively.The pre-processor replaces the rank identifiers with the provided definition.In the example above, the code after running the pre-processor reads A detailed description of the naming convention of the rank identifiers and a list of the currently used rank identifiers is included in the "COSMO/MESSy Implementation Documentation" in the Supplement.

Tracer transport
At the beginning of this section, the previously published submodels for convective tracer transport (CVTRANS), simplified prognostic tracers (PTRAC), point sources and simplified chemistry (TREXP) and tracers of opportunity (DRADON) are shortly introduced, before they are used for the evaluation of the functionality of the TRACER interface of COSMO/MESSy and the tracer transport characteristics.

PTRAC
The MESSy submodel PTRAC (Prognostic TRACers) provides a namelist interface for the definition of tracers.Usually tracers are defined in the SMIL module of a submodel, which requires a recompilation of the code, if a new tracer is added.To provide an easy test-bed for tracer studies without the need of recompilation, PTRAC enables the definition of tracers via namelist.The submodel is described in detail by Jöckel et al. (2008).

TREXP
The regular submodel TREXP provides the possibility to define point sources of tracers and to define tracers experiencing one degradation reaction, either by zero order decay or by a first order reaction with another educt.The tracers, the reaction coefficients and the point sources are all specified in the &CPL namelist of TREXP.The point sources are also applicable to tracers defined elsewhere in MESSy submodels.Jöckel et al. (2010) provide a detailed description of the MESSy submodel TREXP.

DRADON
The regular submodel DRADON (Diagnostic RADON) provides an ideal test bed for the diagnosis of tracer transport.Jöckel et al. (2010) explain the concept of DRADON in detail.Different parameterisations of the radon emission fluxes can be applied.For the test shown in Sect.5.2.3 we applied a constant emission of 1 atom (cm 2 s) −1 over ice and snow free surfaces.

CVTRANS
Convective tracer transport is managed by the MESSy submodel CVTRANS (Tost et al., 2010).The submodel requires as input the convective mass fluxes, which have to be delivered by the applied convection scheme.Currently, in COSMO/MESSy this data delivery is only implemented for the Tiedtke scheme (Tiedtke, 1989), therefore, CVTRANS is presently only applicable in combination with the Tiedtke scheme.

Tracer transport tests
To verify the functionality of the TRACER interface in COSMO/MESSy and to evaluate the tracer transport characteristics, different tests have been performed.The model was used in the on-line coupled mode, as described in Kerkweg and Jöckel (2012) to ensure consistent boundary conditions for the tracers.All tests use either artificial tracers, or 222 Rn.Moisture in the COSMO model is transported by advection, horizontal and vertical diffusion.As the transport of tracers should be consistent with the moisture transport, the COSMO transport routines are adopted for tracer transport in COSMO/MESSy.
Even though horizontal diffusion of tracers is implemented in COSMO/MESSy, all tests showed that the numerical diffusion is larger than the calculated horizontal diffusion.Thus, we recommend to neglect the explicit calculation of horizontal diffusion and show no results for horizontal tracer diffusion here.In addition to the transport processes of the basemodel, the MESSy submodel CVTRANS (Tost et al., 2010) for the convective tracer transport is included.
The COSMO/MESSy model domain for the evaluation is located over Central Europe using a horizontal grid of approximately 40 km (0.36 • ) mesh size and 40 vertical levels.
A second, smaller model domain is located roughly over Germany utilising a horizontal grid of approx.7 km (0.0625 • ) and the same 40 vertical levels as in the 40 km COSMO domain.The simulations will be called COSMO-40 and COSMO-7 hereafter.ECHAM5/MESSy in T106L31 resolution provides the initial and boundary data for the COSMO-40 simulation and COSMO-40 the initial and boundary conditions for the COSMO-7 simulation.The simulated period starts at 1 March 2010.As these are simple tracer transport tests, ECHAM5/MESSy was not nudged to the analysed meteorology.
We first show transport tests for tracers which are initialised but have neither sources nor sinks (Sect.5.2.1), continue with tracers which are initialised with zero and emitted from a point source (Sect.5.2.2) and conclude with a simulation for 222 Rn (Radon).For all these tests the Runge-Kutta time integration and a Semi-Lagrangian advection scheme are used.

Tests with artificial passive tracers
Tracer transport in COSMO/MESSy is operator split into four different processes: advection, convective transport, vertical and horizontal diffusion.For the reason given above, we do hereafter not consider horizontal diffusion.For atmospheric chemistry simulations, each of these transport   sive tracers (i.e., without sources or sinks in the regional domain) is expected to be closed, implying that the tracer mass within the domain plus inflow minus outflow in/out of the domain is conserved.The latter (inflow and outflow) are determined, at least implicitly, by the boundary conditions.As there is no in-and outflux budgeting routine in COSMO/MESSy yet, our analyses are somewhat limited, but as a first step artificial passive tracers are used to diagnose, to the extent possible, the mass conservation, positive definiteness and monotonicity of the transport processes as implemented in COSMO/MESSy.Table 2 lists the initialisation patterns for the artificial tracers.H is initialised homogeneously to a mixing ratio of 10 −9 mol mol −1 .Figure 4 depicts the initialisation of the V 1 tracer.V 2 is initialised reversely to V 1, i.e., its mixing ratio increases from bottom to top.Thus, V 1 and V 2 are defined in a way that V 1 + V 2 is also homogeneously distributed.Since the artificial tracers do not have any internal sources or sinks, their abundance is only changed by transport processes including inflow / outflow across the domain boundaries.The MESSy TRACER interface nicely provides the possibility to switch (on or off) the specific transport processes for each tracer individually.Utilising this, tracers with the three initial conditions (Tab.2) have been combined with the five different transport process switch combinations as listed in Table 3.Thus, in total 15 tracers have been simulated: for  one tracer undergoes all three processes.

Monotonicity
For the analysis of the monotonicity we use the H-tracers, which are initially homogeneously distributed.In an ideal model, none of the transport processes must cause the tracer mixing ratio to deviate from its initial value at any time and any place within the domain (i.e., minimum and maximum value in the entire domain are equal).For a regional model, this monotonicity test also requires the same constant mixing ratio being prescribed at all domain boundaries, which was the case in our setup.The homogeneous tracers Hnt, Hv and Hc stay indeed exactly homogeneous (not shown), indicating that vertical diffusion and convective transport are monotone.The advected tracers exhibit a small but negligible violation of the monotonicity with a deviation from the initial mixing ratio of less than ±0.03‰within 31 days.

Positive definiteness
The generic MESSy sub-submodel TRACER PDEF (Jöckel operators needs to be mass conserving, positive definite and monotone.Inside a regional model domain the mass of a specific tracer is not expected to be conserved.With a perfect transport scheme, however, the mass budget of passive tracers (i.e.without sources or sinks in the regional domain) is expected to be closed, implying that the tracer mass within the domain plus inflow minus outflow in/out of the domain is conserved.The latter (inflow and outflow) are determined, at least implicitly, by the boundary conditions.As there is no in-and outflux budgeting routine in COSMO/MESSy yet, our analyses are somewhat limited, but as a first step artificial passive tracers are used to diagnose, to the extent possible, the mass conservation, positive definiteness and monotonicity of the transport processes as implemented in COSMO/MESSy.Table 2 lists the initialisation patterns for the artificial tracers.H is initialised homogeneously to a mixing ratio of 10 −9 mol mol −1 .Figure 4 depicts the initialisation of the V 1 tracer.V 2 is initialised reversely to V 1, i.e. its mixing ratio increases from bottom to top.Thus, V 1 and V 2 are defined in a way that V 1 + V 2 is also homogeneously distributed.Since the artificial tracers do not have any internal sources or sinks, their abundance is only changed by transport processes including inflow/outflow across the domain boundaries.
The MESSy TRACER interface nicely provides the possibility to switch (on or off) the specific transport processes for each tracer individually.Utilising this, tracers with the three initial conditions (Table 2) have been combined with the five different transport process switch combinations as listed in Table 3.Thus, in total 15 tracers have been simulated: for each initialisation pattern one tracer encounters no transport (nt), three tracers experience exactly one transport process (advection (a), vertical diffusion (v), or convection (c)) and one tracer undergoes all three processes.

Monotonicity
For the analysis of the monotonicity we use the H -tracers, which are initially homogeneously distributed.In an ideal model, none of the transport processes must cause the tracer mixing ratio to deviate from its initial value at any time and any place within the domain (i.e.minimum and maximum value in the entire domain are equal).For a regional model, this monotonicity test also requires the same constant mixing ratio being prescribed at all domain boundaries, which was the case in our setup.The homogeneous tracers H nt, H v and H c stay indeed exactly homogeneous (not shown), indicating that vertical diffusion and convective transport are monotone.The advected tracers exhibit a small but negligible violation of the monotonicity with a deviation from the initial mixing ratio of less than ±0.03 ‰ within 31 days.
Figure 5 displays the integrated (over the model domain) negative tracer masses for the H , V 1 and V 2 tracers, which have been corrected by TRACER PDEF.The H and the V 1 tracers do not become negative, whereas for V 2 advection produces small negative tracer mixing ratios at the beginning of the simulation; after 10 days no more negative values are produced.The fact that V 2 is initialised with zero at the bottom and that the tracer is well mixed in the lower layers after 10 days indicates that the negative values origin in the surface layer, but only if very small mixing ratios prevail.Even then, the negative mass produced for V 2 is 7 orders of magnitude smaller compared to the integrated (over the domain) mass of the tracer, which is in the order of 10 6 kg.

Conservation of mass
Since the air mass in the regional model domain is not conserved, e.g. when a low or high pressure system is moving across the model domain, the tracer mass variation within the model domain must be related to the air mass variation.For tracer mixing ratios in units of mol mol −1 (dry air) the dry air mass is the appropriate measure.Figure 6 (left) shows the dry air and tracer masses (of H and V 1 plus V 2) integrated over the model domain, normalised to their corresponding time averages.Due to this normalisation, the mass variations of the tracers and of the dry air mass are directly comparable.The lines for the variation of the dry air mass and the homogeneous tracer are on top of each other, implying that the transport operators are strictly mass conserving for homogeneously distributed tracers.Moreover, the sum of masses of V 1 and V 2 follows exactly the same variation.Indeed, the deviation of V 1 + V 2 from the initial condition (which is equivalent to H ) is negligibly small (less than ±0.02 ‰) throughout the domain at any time during the simulation (not shown), except for the very beginning, where the negative masses occur (see above).As conclusion the transport operators are strictly linear.The invariance of this linearity relation requires that each change of the mass (per mass of dry air) of V 1 needs to be counterbalanced by a corresponding change of V 2 with opposite sign.This is indeed the case, as is visible in the middle and right panel of Fig. 6: the individual tracer masses (normalised to their time average and divided by the mass of dry air normalised to its time average) of V 1 and V 2 show a complementary variation.As important implication, this must also comprise the in-and outfluxes across the domain boundaries, where, in our setup, the linearity relation V 1+V 2 = H also holds throughout the simulation (not shown).
This complementary variation also appears for the tracers experiencing no transport at all (V 1nt, V 2nt) or only vertical diffusion or convective transport (V 1v, V 2v; V 1c, V 2c), as shown in Fig. 6 (middle and left panel).The amplitudes are smaller, however, since they are not or only vertically redistributed and the only contribution to the in-and outflux is the Newtonian relaxation to the boundary condition at the domain boundary.
Although the analysis so far does not reveal any indication for potential mass conservation violations, they cannot be entirely ruled out: a mass conservation violation caused by mass-wind inconsistencies of the advection operator depends on the vertical tracer (z-axis) gradients (Jöckel et al., 2001).Since in our case, V 1 + V 2 = H throughout the domain at any time, it follows dV 1/dz = −dV 2/dz, because dH /dz = 0.As a consequence, spurious mass changes of V 1 and V 2 due to mass-wind inconsistencies would be of the same absolute value, but of opposite sign.As a consequence they would cancel out and -in our analysis -therefore be indistinguishable from in-and outfluxes across the domain boundaries.In summary, convective transport and vertical diffusion are positive definite, monotone, mass conserving 14 Fig. 6.Left: Integrated dry air and tracer masses normalised to their time averages (see text) for dry air, the homogeneous tracers H and V 1 + V 2. The lines are on top of each other.Middle and right panel: Tracer masses normalised to their time average and divided by the mass of dry air normalised to its time average for the passive tracers V 1 (middle) and V 2 (right).All panels are for the COSMO-7 simulation.mass conservation violation.This has been applied here.Figure 5 displays the integrated (over the model domain) negative tracer masses for the H, V 1 and V 2 tracers, which have been corrected by TRACER PDEF.The H and the V 1 tracers do not become negative, whereas for V 2 advection produces small negative tracer mixing ratios at the beginning of the simulation; after 10 days no more negative values are produced anymore.The fact, that V 2 is initialised with zero at the bottom and that the tracer is well mixed in the lower layers after 10 days indicates that the negative values origin in the surface layer, but only if very small mixing ratios prevail.Even then, the negative mass produced for V 2 is 7 orders of magnitude smaller compared to the integrated (over the domain) mass of the tracer, which is in the order of 10 6  and linear.The advection operator exhibits a small, negligible violation of the positive definiteness, but only for very small mixing ratios in the lowest model layer.Consequently, the correction in the submodel TRACER PDEF (required for numerical stability) results in these cases in a likewise negligible mass conservation violation.Other issues, like potential mass-wind inconsistencies could not be revealed with our tests.Further studies including the budgeting of in-and outflow across the domain boundaries will provide more details in future.Mainly the linearity of the advection operators let us conclude that the tracer transport in COSMO/MESSy is ready for atmospheric chemistry applications.

Tracer tests utilising point sources
As a second test case, tracers have been initialised with zero and are emitted by a point source using the MESSy submodel TREXP (see Sect. 5.1.2).Two artificial tracers (without sinks), emitted at different point sources, have been chosen for the experiment: 1. Tracer PNT: the emission point source of the first tracer is located at 20.00 • W and 50.00 • N in a pressure altitude of 900 hPa.The tracer is named PNT as abbreviation for PoiNT source.This location is chosen to be outside, but close to, the COSMO-40 model domain.
Therefore, the tracer needs to be transported into both COSMO model domains.
2. Tracer VOL: the second tracer is emitted at a pressure altitude of 800 hPa over Island (at 19.60 • W, 63.63 • N).
It is named VOL, as its location coincides with the volcano Eyjafjallajökull.This point resides in the COSMO-40, but not in the COSMO-7 model domain.Thus, VOL is, in contrast to PNT, emitted within the larger COSMO model domain, but is also advected into the COSMO-7 domain.
Figures 7 and 8 display snapshots of the 12th, 15th and 18th simulation day at 12:00 UTC and at the 900 hPa or 800 hPa pressure level for the PNT and VOL tracer, respectively.Pictures in one column correspond to the same simulation day.The first row displays the tracer distribution as simulated by ECHAM5/MESSy, the second and the third row show the tracer distributions in the larger (i.e.40 km) and the smaller (i.e. 7 km) COSMO model domain, respectively, while the last row depicts a composite of all three simulations.
To investigate the transport of the tracers into the regional model domain, the PNT tracer is emitted outside of both COSMO model domains.Comparing the results of the global and the COSMO-40/MESSy model simulation, PNT is advected correctly into the regional model domain.Obviously, details of the tracer distribution are much better resolved in the COSMO/MESSy model than in the ECHAM5/MESSy simulation.These findings are also valid for the comparison of COSMO-40/MESSy to COSMO-7/MESSy. Figure 9 displays the time series of maximum PNT mixing ratios within a rectangle ranging from 5 • W to 25 • E and 36 • to 56 • N at three pressure levels: 850, 900 and 950 hPa.This region is part of all three model domains.A two-peak structure is exhibited at all heights in all model domains.The peaks coincide with the dynamical patterns, when streamers of the PNT tracer are advected into the rectangular region.The lines indicating the maxima as simulated by the COSMO/MESSy models are very similar.In contrast, the maxima of PNT in ECHAM5/MESSy coincide with the other two models at 900 and 950 hPa, while they are -especially for the two peaks -lower at the 850 hPa   level.This is most probably due to less upward transport in the global model.
The VOL tracer is emitted at a point located in the COSMO-40 model domain.As for PNT, the more detailed structures resolved on the finer grids are clearly visible in Fig. 8. Here, the discrepancies between the ECHAM5/MESSy and the COSMO-40/MESSy simulations are larger as for PNT.The VOL mixing ratio is nearly always larger in the COSMO-40/MESSy compared to the ECHAM5/MESSy simulation.This finding is affirmed by Fig. 10, which displays the maximum mixing ratios of the tracer VOL in a similar way as Fig. 9, but for the 750, 800 and 850 hPa pressure levels.The main reason for these differences is the parameterisation of the point source emission in the MESSy submodel TREXP.As the emitted tracer mass is prescribed, the same tracer mass is emitted into grid boxes   1. Tracer PNT: the emission point source of the first tracer is located at 20.00 • W and 50.00 • N in a pressure altitude of 900 hPa.The tracer is named PNT as abbreviation for PoiNT source.This location is chosen to be outside, but close to, the COSMO-40 model domain.
Therefore, the tracer needs to be transported into both COSMO model domains.
2. Tracer VOL: the second tracer is emitted at a pressure altitude of 800 hPa over Island (at 19.60 • W, 63.63 • N).
It is named VOL, as its location coincides with the volcano Eyjafjallajökull.This point resides in the COSMO-40, but not in the COSMO-7 model domain.Thus, VOL is, in contrast to PNT, emitted within the larger COSMO model domain, but is also advected into the COSMO-7 domain.
Figures 7 and 8 display snapshots of the 12th, 15th and 18th simulation day at 12:00 UTC and at the 900 hPa or 800 hPa pressure level for the PNT and VOL tracer, respectively.Pictures in one column correspond to the same simulation day.The first row displays the tracer distribution as simulated by ECHAM5/MESSy, the second and the third row show the tracer distributions in the larger (i.e., 40 km) and the smaller (i.e., 7 km) COSMO model domain, respectively, while the last row depicts a composite of all three simulations.
To investigate the transport of the tracers into the regional model domain, the PNT tracer is emitted outside of both COSMO model domains.Comparing the results of the global and the COSMO-40/MESSy model simulation, PNT is advected correctly into the regional model domain.Obviously, details of the tracer distribution are much better resolved in the COSMO/MESSy model than in the ECHAM5/MESSy simulation.These findings are also valid for the comparison of COSMO-40/MESSy to COSMO-7/MESSy. Figure 9 displays the time series of maximum PNT mixing ratios within a rectangle ranging from 5 • W to 25 • E and 36 • to 56 • N at three pressure levels: 850, 900 and 950 hPa.This region is part of all three model domains.A two-peak structure is exhibited at all heights in all model domains.The peaks coincide with the dynamical patterns, when streamers of the PNT tracer are advected into the rectangular region.The lines indicating the maxima as simulated by the COSMO/MESSy models are very similar.In contrast, the maxima of PNT in ECHAM5/MESSy coincide with the other two models at 900 and 950 hPa, while they are -especially for the two peaks -lower at the 850 hPa level.This is most probably due to less upward transport in the global model.
The VOL tracer is emitted at a point located in the COSMO-40 model domain.As for PNT, the more detailed structures resolved on the finer grids are clearly visible in Fig. 8.Here the discrepancies between the ECHAM5/MESSy and the COSMO-40/MESSy simulations Fig. 9. Maximum mixing ratio (nmol mol −1 ) of tracer PNT.The maximum is calculated in a rectangle ranging from 5 • W to 25 • E and 36 • to 56 • N. The mixing ratios are displayed for three different pressure levels: 850, 900 and 950 hPa.The injection height is 900 hPa.The ECHAM5/MESSy result is displayed in black.The red and the blue lines show the maximum mixing ratios of the COSMO-40 and COSMO-7 domains, respectively.

Radon
As a further test, the MESSy submodel DRADON (Sect.5.1.3)is used, prescribing a constant emission flux of 1 atom (cm 2 s) −1 at the land surface (not covered by snow or ice).The constant emission flux was chosen in order to yield comparable results for all three model domains.
Figure 11 displays snapshots of vertical Radon distributions taken at a cross section at 50 • N at 00:00 UTC.The upper, middle and lower row display the results of the ECHAM5/MESSy, the COSMO-40/MESSy and the COSMO-7/MESSy simulation, respectively.The vertical are larger as for PNT.The VOL mixing ratio is nearly always larger in the COSMO-40/MESSy compared to the ECHAM5/MESSy simulation.This finding is affirmed by Fig. 10, which displays the maximum mixing ratios of the tracer VOL in a similar way as Fig. 9 for the 750, 800 and 850 hPa pressure levels.The main reason for these differences is the parameterisation of the point source emission in the MESSy submodel TREXP.As the emitted tracer mass is prescribed, the same tracer mass is emitted into grid boxes of very different volumes, leading to a higher degree of dilution and thus smaller mixing ratios in ECHAM5/MESSy.Thus the VOL distribution in COSMO-40/MESSy is more compact and shows higher peak values.As the emission point is not located in the COSMO-7 domain, the tracer VOL is only advected into this model domain.Because the lateral boundary data are interpolated from the COSMO-40/MESSy simulation, the maximum mixing ratios of the tracer VOL in COSMO-7/MESSy reflect mostly the maxima of the COSMO-40/MESSy simulation.

Radon
As a further test, the MESSy submodel DRADON (Sect.5.1.3)is used, prescribing a constant emission flux of 1 atom (cm 2 s) −1 at the land surface (not covered by snow or ice).The constant emission flux was chosen in order to yield comparable results for all three model domains.
Figure 11 displays snapshots of vertical Radon distributions taken at a cross section of 50 • N at 00:00 UTC.The upper, middle and lower row display the results of the ECHAM5/MESSy, the COSMO-40/MESSy and the COSMO-7/MESSy simulation, respectively.The vertical structure of the 222 Rn distribution is mainly driven by convective transport.Due to the trigger mechanisms of convection parameterisations, and due to the different grid box sizes the distribution is expected to differ in some places.It should be noted that the ECHAM5/MESSy and the COSMO/MESSy model, use convection schemes based on the Tiedtke scheme (Tiedtke, 1989), but including different further developments.On 6 March, the convective transport of the COSMO/MESSy simulations is very similar, with COSMO-7/MESSy showing some more individual updrafts than COSMO-40/MESSy.While in the COSMO/MESSy models the highest updrafts are located between 10 • W to 0 • W and 20 • E to 32 • E, the ECHAM5/MESSy simulation shows a prominent updraft around 5 • E. This is most probably due to the different trigger mechanisms employed in the respective convection schemes.On 11 March a bubble of Radon formed in the free troposphere, detached from the surface and the boundary layer.This bubble exists in all simulations, but while the global model predicts an almost even distribution within this bubble, the COSMO/MESSy simulations show much more structure inside.On 19 March all simulations provide relatively similar well mixed Radon distributions (despite the differences in the scales), suggest-  The implementation of the MESSy interface and some regular submodels into the regional weather prediction and climate model COSMO is documented, aiming at a new limited-area chemistry model.The most important code changes, both in the COSMO model code and in the MESSy code, which were required for the implementation of the MESSy interface into the COSMO model are briefly listed in Appendix B, while the details are provided in the Supplement.The further developments and modifications are implemented for the full model system, i.e. they either do not interfere with EMAC (because they are not relevant for it), or they are immediately available also for EMAC.
The code changes are entirely transparent for the COSMO model and EMAC users.For instance, the changes in the structure of MESSy required for the generalisation of the interface are all included in MESSy development cycle 2 as published by Jöckel et al. (2010).The changes in the COSMO model code are introduced via pre-processor directives, such that the standard COSMO code is compiled, if the pre-processor directive is not activated.
Due to the implemented generalisations it is possible to include the MESSy infrastructure into other models as well.The feasibility to apply specific regular MESSy submodels depends, however, on the processes they describe and if they are applicable to the domain and scale of the basemodel.For instance, the global ocean model MPIOM, implemented as a MESSy submodel by Pozzer et al. (2011), is not straightforwardly applicable to the regional domain and thus not used in COSMO/MESSy.Similarly, the application of the CVTRANS submodel for convective transport (Tost et al., 2010) is reasonable for COSMO/MESSy studies on scales requiring convection parameterisations, whereas it is inappropriate for convection resolving simulations.The choice of appropriate submodels for a specific scientific purpose is, at the end, the responsibility of the user.
One prerequisite for the next step, namely the inclusion of the submodels relevant for atmospheric chemistry applications with COSMO/MESSy, is the correct functionality of the TRACER infrastructure, which handles the data and meta-information of chemical constituents.This functionality has been successfully demonstrated with three different tests for artificial passive (i.e.without sources and sinks) tracers, with point sources and with 222 Rn with a standard source distribution.The tests involved three regular MESSy submodels (PTRAC, TREXP and DRADON).The results furthermore indicate that no severe violations of mass conservation, monotonicity or positive definiteness of the tracer transport are apparent.From that, we conclude that the model is ready for atmospheric chemistry applications, but further (regular) tests, which also budget the in-and outflow across the domain boundaries are desirable.
The implementation and the successful testing of the tracer transport forms an important first step and a prerequisite towards the development of a limited-area atmospheric chemistry model.The application and further development of the Modular Earth Submodel System (MESSy) to achieve this goal has several advantages: the same or similar process formulations can be used in both, the global (ECHAM5/MESSy) and the limited-area (COSMO/MESSy) model.Furthermore, the required dynamical and chemical boundary conditions for the limited-area model can be seamlessly taken from the global model.In this way, a highly consistent, scale-bridging model chain is constructed.As a powerful research tool, it can be applied for a wide variety of scientific questions: it can serve as a zooming option for atmospheric chemistry related campaign and mission support, for regional air pollution studies, for the dynamical downscaling of the chemical state of the atmosphere on various time scales.It can support the development of new parameterisations for global chemistry climate models, etc.Our ongoing developments therefore focus on the even closer, i.e. on-line coupling of the limited-area model to the global model, the evaluation of the coupled model system (at first w.r.t. the simulated meteorology) and on the complete integration of all atmospheric chemistry related processes in the limited-area model.
Consequently, this is only the first in a series of three companion articles about the new modeling system ECHAM5/MESSy -COSMO/MESSy (abbreviated MECO(n)).While the basis for a regional chemistry model, i.e. the implementation of the MESSy infrastructure and tracer transport into the COSMO model, is presented here, the second article (Kerkweg and Jöckel, 2012) is about a different issue, but based on the COSMO/MESSy model of this first part.It describes a newly developed coupling methodology, for the on-line nesting of COSMO/MESSy instances into ECHAM5/MESSy.Last but not least, the third article (Hofmann et al., 2012) provides a first meteorological evaluation of the newly developed on-line coupled model chain described technically in the second article.This evaluation, focusing on distinct meteorological events on synoptic scale, and on the question if and how they can be reproduced by MECO(n), is a prerequisite for further applications with chemistry, like chemical weather (air pollution) forecasts, measurement campaign analyses etc. with a focus on the representation of distinct meteorological situations.meta-information.The "relation" can be, for instance, the simple fact that the channel objects are defined by the same submodel.
-channel object: it represents a data field including its meta-information and its underlying geometric structure (representation), e.g. the 3-dimensional vorticity in spectral representation, the ozone mixing ratio in Eulerian representation, the pressure altitude of trajectories in Lagrangian representation.
-dimensions: they represent the basic geometry of one dimension, e.g. the number of latitude points, the number of trajectories, etc.
-event: this is a data type provided by the generic submodel TIMER, which is used to schedule processes at specific (regular) time intervals, e.g. to trigger regular output or input during a simulation.The event control is part of the MESSy generic submodel TIMER.The supplement of Jöckel et al. (2010) comprises a manual for TIMER and details about the event definition.
-rerun event: it triggers the output of restart files.
-restart: a restart is performed, if the computing time allowed by a job scheduler is limited and too short for the complete simulation.In this case, the simulation is interrupted in between and restarted in a new job.To achieve binary identical results for simulations with and without interruption, restart files are written, of which the contents fully determine the state of a model simulation.These files are read in the initialisation phase during a model restart.
-tracer set: a group of tracers in the same representation (see Jöckel et al., 2008).-Prognostic and diagnostic fields are (de)allocated and initialised by the MESSy memory management instead of in alloc_meteofields.
-The calculation of the diagnostic fields for the output is additionally called in those time steps, when MESSy output is triggered.
-Some local variables of COSMO are replaced by POINTER defined in the MESSy submodel DATA.The POINTER are used in COSMO instead of the local variables in order to make the information available to the MESSy submodels.
-The clock of COSMO is set (and synchronised) by (with) the MESSy submodel TIMER.

B2.1 Infrastructure
Naturally, most extensions for the implementation of MESSy into COSMO are required for the infrastructure submodels, because these form the adapter to the respective basemodel.
Here, we name only the major changes: -CONTROL Only those submodels applicable to COSMO are called.
-CHANNEL Representations are defined corresponding to the COSMO grid.Restart checks depend on the basemodel.
-DATA This code is naturally specific for COSMO (i.e.basemodel dependent; see manual).
-MPI Shared high level interfaces are used, which internally call the corresponding COSMO routines.

B2.2 Regular submodels -
The most important change is the introduction of the rank-identifier for flipping the ranks of the variables between the ECHAM and COSMO order of dimensions.
-Consequently, the definition of representations needs to be adapted to the dimension order of the respective basemodel.

Fig. 1 .
Fig. 1.Sketch of the four layer structure of MESSy.

Fig. 1 .
Fig. 1.Sketch of the four layer structure of MESSy.

Fig. 2 .
Fig. 2. List of all entry points available in ECHAM5/MESSy.The entry points highlighted with yellow boxes are currently implemented in COSMO/MESSy, the hatched ones are not yet required, as the dependent submodels are not used in COSMO/MESSy so far.Next to the boxes the generic (blue) and regular (black) submodels currently called in COSMO/MESSy by the respective CONTROL subroutine are listed.The green box surrounds those entry points located in the outer loop over the second horizontal dimension in ECHAM5/MESSy.

Figure 3
Figure 3 depicts a flow chart of the COSMO/MESSy model.It illustrates where the individual MESSy entry points are located in the COSMO basemodel flow.The MESSy entry points are highlighted by yellow boxes, whereas the COSMO model routines are coloured in blue.COSMO model routines which are obsolete for COSMO/MESSy are crossed out and substituted by the corresponding MESSy subroutine in the dark orange box directly below.The light orange boxes indicate direct calls of generic submodel subroutines.The "COSMO/MESSy Implementation Documentation" in the Supplement lists and explains all changes, which became necessary in the COSMO model code for the implementation of the MESSy interface.

Fig. 3 .
Fig. 3. Simplified flow chart of COSMO/MESSy.Listed are the main routine calls of COSMO and the MESSy entry points.Indented boxes are called by the non-indented box above.Blue boxes indicate COSMO routines; yellow boxes highlight MESSy entry points (CONTROL); orange boxes point to direct calls of MESSy infrastructure submodels and dark orange boxes indicate direct calls of COSMO specific MESSy infrastructure submodel routines.

Fig. 3 .
Fig. 3. Simplified flow chart of COSMO/MESSy.Listed are the main routine calls of COSMO and the MESSy entry points.Indented boxes are called by the non-indented box above.Blue boxes indicate COSMO routines; yellow boxes highlight MESSy entry points (CONTROL); orange boxes point to direct calls of MESSy infrastructure submodels and dark orange boxes indicate direct calls of COSMO specific MESSy infrastructure submodel routines.

Fig. 4 .
Fig. 4. Initialisation pattern of the passive tracer V 1.The horizontal axis shows rotated coordinates.

Fig. 5 .
Fig. 5. Corrected negative tracer mass (see text) in kg for the passive tracers H (left), V 1 (middle) and V 2 (right) in the COSMO-7 region.For H and V 1 all lines are on top of each other.

Fig. 4 .
Fig. 4. Initialisation pattern of the passive tracer V 1.The horizontal axis shows rotated coordinates.

Fig. 5 .
Fig. 5. Corrected negative tracer mass (see text) in kg for the passive tracers H (left), V 1 (middle) and V 2 (right) in the COSMO-7 region.For H and V 1 all lines are on top of each other.

Fig. 6 .
Fig. 6.Left: integrated dry air and tracer normalised to their time averages (see text) for dry air, the homogeneous tracers H and V 1 + V 2. The lines are on top of each other.Middle and right panel: tracer masses normalised to their time average and divided by the mass of dry air normalised to its time average for the passive tracers V 1 (middle) and V 2 (right).All panels show results from the COSMO-7 simulation.

Fig. 7 .
Fig. 7. Horizontal distribution at 900 hPa of the artificial tracer PNT.The location of the emission point is indicated by the light blue plus sign.Results are shown for the 12th, 15th and 18th simulation day at 12:00 UTC (columns).First row: ECHAM5/MESSy, second row COSMO-40/MESSy, third row COSMO-7/MESSy and last row composite of all three model domains.

Fig. 7 .
Fig. 7. Horizontal distribution at 900 hPa of the artificial tracer PNT.The location of the emission point is indicated by the light blue plus sign.Results are shown for the 12th, 15th and 18th simulation day at 12:00 UTC (columns).First row: ECHAM5/MESSy, second row COSMO-40/MESSy, third row COSMO-7/MESSy and last row composite of all three model domains.

Fig. 8 .
Fig. 8. Horizontal distribution at 800 hPa of the artificial tracer VOL.The location of the emission point is indicated by the light blue plus sign.Results are shown for the 12th, 15th and 18th simulation day at 12:00 UTC (columns).First row: ECHAM5/MESSy, second row COSMO-40/MESSy, third row COSMO-7/MESSy and last row composite of all three model domains.

Fig. 8 .
Fig. 8. Horizontal distribution at 800 hPa of the artificial tracer VOL.The location of the emission point is indicated by the light blue plus sign.Results are shown for the 12th, 15th and 18th simulation day at 12:00 UTC (columns).First row: ECHAM5/MESSy, second row COSMO-40/MESSy, third row COSMO-7/MESSy and last row composite of all three model domains.

Fig. 9 .
Fig. 9. Maximum mixing ratio (nmol/mol) of tracer PNT.The maximum is calculated in a rectangle ranging from 5 • W to 25 • E and 36 • to 56 • N. The mixing ratios are displayed at three different pressure levels: 850, 900 and 950 hPa.The injection height is 900 hPa.The ECHAM5/MESSy result is displayed in black.The red and the blue lines show the maximum mixing ratios of the COSMO-40 and COSMO-7 domains, respectively.

Fig. 10 .
Fig. 10.Maximum mixing ratio (nmol/mol) of tracer VOL.The maximum is calculated in a rectangle ranging from 5 • W to 25 • E and 36 • to 56 • N. The mixing ratios are displayed at three different pressure levels: 750, 800 and 850 hPa.The injection height is 800 hPa.The ECHAM5/MESSy result is displayed in black.The red and the blue lines show the maximum mixing ratios of the COSMO-40 and COSMO-7 domains, respectively.

Fig. 10 .
Fig. 10.Maximum mixing ratio (nmol mol −1 ) of tracer VOL.The maximum is calculated in a rectangle ranging from 5 • W to 25 • E and 36 • to 56 • N. The mixing ratios are displayed for three different pressure levels: 750, 800 and 850 hPa.The injection height is 800 hPa.The ECHAM5/MESSy result is displayed in black.The red and the blue lines show the maximum mixing ratios of the COSMO-40 and COSMO-7 domains, respectively.

-
Prognostic and diagnostic fields are transformed form ALLOCATABLE fields to POINTER in order to be allocated and used by the MESSy memory management.This also includes the change of IF (ALLOCATED(field)) into IF (ASSOCIATED(field)) throughout the code.

Geosci. Model Dev., 5, 87-110, 2012 www.geosci-model-dev.net/5/87/2012/
Therefore fr_land is renamed in DATA to slf by .e. at the end of the previous time step), mostly indicated by 'm1' at the end of the variable name.Additionally, the ECHAM5/MESSy model utilises variables for the tendencies applied to the prognostic variables during the current time step (mostly indicated by 'te' at the end of the variable name).In contrast to the ECHAM5/MESSy model, the COSMO model contains 2-and 3-time level integration schemes.To simplify the treatment for these different schemes and to avoid copying one time level of a prognostic variable to another time level of the same prognostic variable, the variable fields in the COSMO model are of rank 4, consisting of the three dimensions in space and an additional rank dimensioned by the number of time levels utilised in the integration scheme.To access the different time levels, index variables of type INTEGER (nnew, nnow and nold) are defined.Their values are rotated every time step, thus the 'm1' or 'nnow' values are not located at the same memory space all the time during the integration.The tendency of a variable in the COSMO model is defined in a similar way as in DATA also provides fields required by the MESSy submodels, which are not directly supplied by the basemodel.In DATA additional channel objects are defined and calculated during the integration phase (if they vary with time) or in the initialisation phase (if they are con- -An additional challenge is the treatment of the time levels of the prognostic variables in both models.As ECHAM5/MESSy uses the leap-frog time integration scheme, it provides variables for the values after time integration (iECHAM5/MESSy.It is mostly indicated by 'tens' at the end of the variable name: e.g. the specific humidity in ECHAM5 is determined by the two rank-3 variables qm1 and qte.In the COSMO model the specific humidity field is the rank-4 variable qv.The field after the last completed time step is accessed by qv(:,:,:,nnow) for the 2-time level scheme or by qv(:,:,:,nold) for the leap-frog scheme.The corresponding specific humidity tendency field in COSMO is defined by qvtens.To minimise the changes in the MESSy submodel interface layer, POINTERs, named as the corresponding ECHAM5/MESSy variables, are associated to the current time slice and/or the tendency variable, e.g.:REAL(DP),POINTER,PUBLIC, & DIMENSION(:,:,:) :: qm1 => NULL() pendent.Therefore DATA copies the COSMO model dimension variables to variables named as the respective MESSy variables.The "COSMO/MESSy Implementation Documentation" in the Supplement contains a list of the dimension variable names and their meaning.Geosci.Model Dev., 5, 87-110, 2012 www.geosci-model-dev.net/5/87/2012/ - and 1D LOGI-CAL, CHARACTER (LEN=*) and different KINDs of INTEGER or REAL variables.It is a renamed and expanded version of the COSMO model interface distribute_values, which requires six parameters: the buffer to be sent/received, the rank of the task in the group communicator, the group communicator, the buffer length, the MPI datatype and an error status.The corresponding MESSy routine (p_bcast as in ECHAM5), in contrast, requires only two parameters (buffer and rank of the sender), others are optional.To create a common interface, the group communicator, the buffer length and the MPI datatype have been made optional parameters in the COSMO/MESSy section of the generic submodel MPI.If they are not present, they are deduced from other input:

Table 1 .
Definition of tracer field instances in COSMO/MESSy.The middle column lists the variable names of the respective fields in TRACER.The abbreviations RK and LF denote the Runge-Kutta and Leap-frog scheme, respectively.

Table 2 .
Definition of initialisation patterns for passive tracers used in this study.

Table 3 .
Table of passive tracers for the transport test with PTRAC.X is a placeholder for either H, V 1 or V 2 and refers to the initialisation patterns (Tab.2).

Table 3 .
Table of passive tracers for the transport test with PTRAC.X is a placeholder for either H , V 1 or V 2 and refers to the initialisation patterns (Tab.2).