\def\myfiledate{2011-04-14}

\documentclass[a4paper,twoside]{article}
\usepackage{ucs}\usepackage[utf8x]{inputenc}
\usepackage{natbib}
\usepackage{chem}
%\usepackage{afterpage}
\usepackage{url}
\usepackage{color}
%\usepackage{multicol}
\usepackage{rotating} % loads graphicx
%\usepackage{longtable}
\usepackage{graphicx}
%\usepackage{eclclass}
%\usepackage{verbatim}
\usepackage[pdftex,bookmarks=false,colorlinks]{hyperref}
\hypersetup{anchorcolor=black,citecolor=black,filecolor=black,linkcolor=black,%
  menucolor=black,pagecolor=black,urlcolor=black}

% activate to show "draft" watermark:
%\usepackage{pdfdraftcopy}\draftstring{DRAFT}\draftfontsize{150}\draftangle{60}

\oddsidemargin-5mm
\evensidemargin-5mm
\topmargin-15mm
\textheight250mm
\textwidth175mm
\raggedbottom
\parindent0mm
\parskip1.0ex plus0.5ex minus0.5ex
\renewcommand{\arraystretch}{1}
\renewcommand{\topfraction}{0.95}
\renewcommand{\dbltopfraction}{0.95}
\renewcommand{\bottomfraction}{0.95}
\renewcommand{\floatpagefraction}{0.9}
\renewcommand{\dblfloatpagefraction}{0.9}
\renewcommand{\textfraction}{0.1}
\setcounter{topnumber}{3}
\setcounter{bottomnumber}{3}

\newcommand{\todo}[1]{{\uppercase{\bf ((#1))}}}
\newcommand{\egcite}[1]{\citep[e.g.][]{#1}}

\def\nosep{\setlength\parsep{0mm}\setlength\topsep{0mm}\setlength\itemsep{0mm}}

\def\mypageheader{Sander et al.: CAABA/MECCA User Manual}
\markboth{\mypageheader}{\mypageheader}
\pagestyle{myheadings}

\input{mecca_info.tex} % \def\meccaversion{...}

% line break after \paragraph:
\makeatletter
\renewcommand\paragraph{\@startsection{paragraph}{4}{\z@}%
  {-2.0ex\@plus -1ex \@minus -.2ex}%
  {1.0ex \@plus .2ex}%
  {\normalfont\normalsize\bfseries}}
\makeatother

\begin{document}

\thispagestyle{empty}

\begin{center}
  \includegraphics[height=0.3\textheight]{caaba_mecca_logo_print}\\[10mm]
  {\Huge\bf CAABA/MECCA-{\meccaversion}}\\[5mm]
  {\Huge\bf User Manual}\\[10mm]
  {\huge\em \underline{C}hemistry \underline{A}s \underline{A}
    \underline{B}oxmodel \underline{A}pplication /}\\[3mm]
  {\huge\em \underline{M}odule \underline{E}fficiently
    \underline{C}alculating the\\[5mm]
    \underline{C}hemistry of the
    \underline{A}tmosphere}\\[10mm]
  {\huge\bf Rolf Sander et al.}\\[10mm]
  \Large
  Air Chemistry Department\\
  Max-Planck Institute of Chemistry\\
  PO Box 3060, 55020 Mainz, Germany\\
  \url{rolf.sander@mpic.de}\\[10mm]
  {\huge\url{http://www.mecca.messy-interface.org}}

  \vfill

  {\large This manual is part of the electronic supplement of our
    article\\
    ``The atmospheric chemistry box model CAABA/MECCA-3.0''\\
    in Geosci.\ Model Dev.\ (2011), available at:\\
    \url{http://www.geosci-model-dev.net}}

  Date: \myfiledate

\end{center}

\clearpage

\setlength{\columnsep}{8mm}
\twocolumn\sloppy

\tableofcontents

\begin{figure*}%[htb]
  \begin{center}
  \includegraphics[width=0.9\textwidth]{modular-mecca}
  \end{center}
  \caption{Diagram showing MECCA as part of the CAABA box model or of a
    global model.}
  \label{fig:modular-mecca}
\end{figure*}

\section{Introduction}

MECCA (\underline{M}odule \underline{E}fficiently
\underline{C}alculating the \underline{C}hemistry of the
\underline{A}tmosphere) is an atmospheric chemistry module that contains
a comprehensive chemical mechanism with tropospheric and stratospheric
chemistry of both the gas and the aqueous phase \citep{1666}. For the
numerical integration, MECCA uses the KPP software \citep{1665}.

To apply the MECCA chemistry to atmospheric conditions, MECCA must be
connected to a base model. As shown in Fig.~\ref{fig:modular-mecca}, the
base model can be a complex, 3-dimensional model \egcite{1851} but it
can also be a simple box model. The connection is established via the
MESSy interface (\url{http://www.messy-interface.org}) developed by
\citet{1664}.

This manual describes how to install and work with MECCA when it is
connected to the box model CAABA (\underline{C}hemistry \underline{A}s
\underline{A} \underline{B}oxmodel \underline{A}pplication). This
combination will be referred to as ``CAABA/MECCA''. The main features of
the CAABA box model are shown in Fig.~\ref{fig:caaba_sketch}. In
addition to MECCA chemistry, CAABA also contains several other
submodels, e.g.\ JVAL for calculating J-values, SAPPHO for simplified
and parameterized photolysis rates, and SEMIDEP for simplified emission
and deposition.

\section{Installation}
\label{sec:install}

This section can be skipped if CAABA/MECCA is already installed on your
computer.

\subsection{System Requirements}

\subsubsection{Linux/Unix}

CAABA/MECCA has been tested successfully on several UNIX-like operating
systems. The easiest installation is probably on a Linux PC since
several auxiliary programs are already included in a typical Linux
distribution.

\subsubsection{MAC OS X}

CAABA/MECCA does not work with the version of the \verb|sed| program
that is shipped with MAC OS X. Instead, it is necessary to install the
GNU version of \verb|sed|, called \verb|gsed|. This can be done using
the MacPorts (\url{http://www.macports.org/}). To ensure that
\verb|gsed| is always executed when \verb|sed| is called, a symbolic
link from \verb|sed| to \verb|gsed| can be created, e.g.:
\begin{verbatim}
sudo ln -s /opt/local/bin/gsed
  /opt/local/bin/sed
\end{verbatim}
Also, there may be problems executing the command ``\verb|echo -n|''
under OS X. If this is the case, the script \verb|xmecca| needs to be
adjusted.

\subsubsection{Windows}

A native installation under Windows is neither recommended nor
supported. However, it is possible to execute the model in a virtual
machine running Linux on a Windows PC. The VMware Player, which is
needed to run the virtual machine, can be obtained from
\url{http://www.VMware.com}. Information how to obtain the Ubuntu Linux
image file with the CAABA/MECCA model (about 8 GB) from our ftp site is
given on the MECCA web page at
\url{http://www.mecca.messy-interface.org}.

\subsection{Prerequisites}
\label{sec:prerequisites}

\begin{description}
\item[A Fortran95 compiler (mandatory):] Several compilers have been
  tested successfully: g95 (for Linux), Lahey/Fujitsu (for Linux), Intel
  (for Linux), Compaq (Alpha UNIX). Other compilers can be used as well,
  if they accept standard Fortran95 code. It should be noted that the
  g95 compiler for Linux is free and can be downloaded from
  \url{http://www.g95.org/}.
\item[The Kinetic PreProcessor KPP (mandatory):] This flexible numerical
  integration package by \citet{1665} transforms the chemistry mechanism
  into a set of ordinary differential equations (ODEs) in Fortran95
  syntax. MECCA needs the KPP version that is provided in the
  \verb|mecca/kpp/| directory. Follow the instructions in
  \verb|mecca/kpp/readme| to install KPP.
\item[Perl, tcsh, gawk, sed, and gmake (mandatory):] These UNIX tools
  are standard on Linux systems. Check that recent versions of them are
  installed. Especially gawk may lead to strange error messages.
  To test gawk, type:\\
  \verb|gawk 'BEGIN {print match("X","[^a-z]")}'|\\
  The result should be ``1''. However, you may get ``0'' as the result
  on your system. Supposedly, this is not a bug in gawk but a feature.
  You can solve the problem by setting the environment variable
  \verb|LC_ALL| to ``\verb|C|'':\\
  \verb|export LC_ALL=C| $\qquad$ (if you use bash)\\
  \verb|setenv LC_ALL C| $\qquad$ (if you use tcsh)\\
  When you try the gawk test again, it should work fine.
\item[La\TeX\ (optional):] If you have La\TeX\ installed on your
  computer, you can print a table (including rate coefficients and
  references) of the currently selected mechanism (see
  Sect.~\ref{sec:latextable} for details).
\item[netCDF library (strongly recommended):] The netCDF library is
  needed to create model output in netCDF format. It can be obtained
  from \url{http://www.unidata.ucar.edu/software/netcdf/}. Note that the
  \verb|*.mod| files in the netCDF library are compiler-specific. Thus,
  it is necessary to create a netCDF library for each Fortran95 compiler
  and maybe also for each Fortran95 compiler version.

  Software for manipulating or displaying netCDF data is listed at:
  \url{http://www.unidata.ucar.edu/software/netcdf/software.html}. If
  you don't have the netCDF library, you can still run the model but
  produce only ASCII output.
\item[netCDF tools (optional):] Several tools are used to analyse the
  netCDF output when the model is run in the multirun
  (Sect.~\ref{sec:multirun}) or Monte-Carlo (Sect.~\ref{sec:montecarlo})
  mode. Specifically, the NCO programs \verb|ncpdq|, \verb|ncrcat|, and
  \verb|ncks| from \url{http://nco.sourceforge.net}, the program
  \verb|ncdump|, and the program \verb|ncclamp| from
  \url{http://ncclamp.sourceforge.net} are needed.
\item[ferret (optional):] The ferret plotting program is needed to plot
  the contents of the netCDF output using the ferret scripts in the
  \verb|jnl/| directory (see Sect.~\ref{sec:ferret} for details). To
  ensure that ferret finds all necessary files, you have to add
  ``\verb|./tools|'' to the \verb|FER_GO| environment variable. The
  \verb|ferret_paths*template| files show how to do this. For example,
  when using the tcsh, type:\\
  \verb|setenv FER_GO "$FER_GO ./tools"|
\item[graphviz (optional):] The graphviz program can be used to create
  graphical visualizations of the reaction mechanism like
  Fig.~\ref{fig:Br}.
\item[fpc Pascal compiler (optional):] To use the tagging diagnostics
  and isotope modeling features (see Sect.~\ref{sec:tagdbl} for
  details), version 2.2.0 or higher of the free pascal compiler
  \verb|fpc| from \url{http://www.freepascal.org/} is needed.
\end{description}

\begin{figure*}%[htb]
  \begin{center}
  \includegraphics[width=0.9\textwidth]{caaba_sketch}
  \end{center}
  \caption{The CAABA box model}
  \label{fig:caaba_sketch}
\end{figure*}

\begin{figure*}
  \begin{center}
    \includegraphics[width=\textwidth]{kpp_use_diagr}
  \end{center}
  \caption{Module structure of KPP-produced Fortran95 files. The arrows
    start at the module which is exporting the {\tt PUBLIC} variables
    and subroutines (which are shown in blue). They point to the module
    importing them via the Fortran95 USE instruction. A dotted line
    represents an optional connection.}
  \label{fig:kpp_use_diagr}
\end{figure*}

\begin{figure*}%[htb]
  \begin{center}
  \includegraphics[width=0.9\textwidth]{xcaaba_flowcontrol}
  \end{center}
  \caption{Illustration of the tasks performed by {\tt xcaaba}. {\tt
      xcaaba} and all scripts called by {\tt xcaaba} are shown on a blue
    background. User-generated (static) input files are shown on a
    yellow background whereas automatically generated temporary files
    are shown on a white background.}
  \label{fig:xcaaba_flowcontrol}
\end{figure*}

\begin{figure*}%[htb]
  \begin{center}
  \includegraphics[width=\textwidth]{caaba_flowcontrol}
  \end{center}
  \caption{Flow control of a CAABA box model simulation}
  \label{fig:caaba_flowcontrol}
\end{figure*}

\begin{figure*}
  \begin{center}
    \fbox{\includegraphics[width=0.9\textwidth]{Br}}
  \end{center}
  \caption{Visualization of the MECCA gas-phase bromine chemistry
    generated with the graphviz software.}
  \label{fig:Br}
\end{figure*}

\subsection{Installation}

Once all prerequisites are fulfilled, you can install CAABA/MECCA by
simply unpacking the zip archive:\\[2mm]
{\tt unzip caaba\_\meccaversion.zip}\\[2mm]

Next, you have to check that all settings in \verb|Makefile| are
correct. If necessary, edit the file: Choose a Fortran95 compiler
(\verb|COMPILER|), enter its name (\verb|F90|) and the compiler options
(\verb|F90FLAGS|). If you add a new compiler, check if you need to
activate the C-preprocessor option. To activate netCDF output, you also
have to edit the \verb|Makefile|:
\begin{itemize}\nosep
\item Check that the variable \verb|OUTPUT| is set to \verb|NETCDF| (not
  to \verb|ASCII|).
\item Enter the correct netCDF library information in
  \verb|NETCDF_INCLUDE| and \verb|NETCDF_LIB|.
\end{itemize}

\subsection{Troubleshooting}

Should there be any problems with the CAABA/MECCA installation, please
check the following:
\begin{itemize}\nosep
\item Confirm that all prerequisites (see above) are fulfilled!
\item Confirm that the perl path in the first line of
  \verb|sfmakedepend| is correct. It should be the same as the output of
  the command:\\
  \verb|which perl|
\item Confirm that the tcsh paths in the first lines of
  \verb|xcaaba| and \verb|xmecca| are correct.
\item Confirm that the model code was unzipped successfully from the zip
  archive. Check for potential problems during the unzipping process:
  \begin{itemize}
  \item Make sure that the directory structure has not changed.
    Unfortunately, some unzipping programs seem to put all files into
    one directory, ignoring the original directory structure.
  \item Make sure that links have not been converted to files. For
    example, the output of the command ``\verb|file caaba.nml|'' should
    tell you that \verb|caaba.nml| is a symbolic link to
    \verb|nml/simple/caaba.nml|.
  \end{itemize}
\end{itemize}

\section{Compiling and running the CAABA/MECCA box model with the shell
  script {\tt xcaaba}}
\label{sec:execute}

First, go to the base directory of the model code (note that all path
names given in this manual are relative to this base directory):\\[2mm]
{\tt cd caaba\_\meccaversion}\\[2mm]
Next, the tcsh script \verb|xcaaba| will guide you through the process
of running the box model, as illustrated in
Fig.~\ref{fig:xcaaba_flowcontrol}. To execute \verb|xcaaba|, type:
\begin{verbatim}
./xcaaba
\end{verbatim}
\verb|xcaaba| will ask several questions, and recommended answers are
given below. If you only press the Return key, you select the default.
\begin{verbatim}
Start xmecca?
\end{verbatim}
If you answer ``\verb|y|'', you can create a new chemical mechanism with
\verb|xmecca| as described in detail in Sect.~\ref{sec:xmecca}. However,
for the first tests with CAABA/MECCA it is recommended to answer
``\verb|n|'' and use the simple default mechanism.
\begin{verbatim}
Choose an option:
s = Start from scratch
c = Compile
r = Run existing executable
h = Help
q = Quit
\end{verbatim}
Choose ``\verb|c|'' to compile the Fortran95 code. After a successful
compilation, \verb|xcaaba| asks you to choose a namelist file:
\begin{verbatim}
Choose a namelist file from
the nml/ directory:
...
\end{verbatim}
Namelists control the behaviour of CAABA/MECCA during run-time, and
editing them allows to define the model setup at run-time (see
Sect.~\ref{sec:nmlfiles}). The default is to use the same namelist as
last time. For the first tests, the file \verb|caaba_simple.nml| can be
chosen. The active contents of the chosen namelist will be shown. Next,
\verb|xcaaba| asks if you want to run the model:
\begin{verbatim}
Run CAABA boxmodel with MECCA chemistry?
y = Yes (default)
n = No
q = Quit
\end{verbatim}
Answer ``\verb|y|'', and the CAABA/MECCA model simulation will start.
The flow control is illustrated in Fig.~\ref{fig:caaba_flowcontrol}. The
model day and the current solar zenith angle (sza) are printed on the
screen during the model simulation. The default is to integrate 8 days.
\begin{verbatim}
Save the output and model code
in output/ directory?
\end{verbatim}
Answer ``\verb|y|'', and \verb|xcaaba| will put the files into a
subdirectory with a name based on the date and time of the model
simulation, e.g.\ \verb|output/2011-04-14-12:43:35/|. It is recommended
to rename the subdirectory to a more descriptive name.

\section{Selecting a chemical mechanism with the shell script {\tt xmecca}}
\label{sec:xmecca}

MECCA contains a very comprehensive set of chemical reactions in both
the gas phase and the aqueous phase. For many applications, using the
complete mechanism will consume too much CPU time. Therefore, the shell
script \verb|xmecca| has been written which allows to create a
custom-made subset of the chemical mechanism interactively. Normally,
\verb|xmecca| is called via \verb|xcaaba|. However, you can also start
it manually:
\begin{verbatim}
cd mecca
./xmecca
\end{verbatim}
\verb|xmecca| will ask several questions, and recommended answers are
given below. If you only press the Return key, you select the default.
\begin{verbatim}
Select a batch file which defines
the chemistry mechanism that you
want to generate. 
\end{verbatim}
It is strongly recommended that you select a batch file here. Batch
files contain all the information that \verb|xmecca| needs to create a
chemical mechanism. Several batch files are available already, and it is
also possible to add your own batch files as explained at the end of
this section. If you do not select a valid batch file, you can continue
and answer all questions interactively as described here.
\begin{verbatim}
How many aerosol phases?
\end{verbatim}
For a gas-phase only mechanism, type ``0''. For a mechanism with
aqueous-phase chemistry in seasalt and in sulfate particles, type ``2''.
Other values are possible if they have been defined in subroutine
\verb|define_aerosol| in \verb|messy_mecca_box.f90|.
\begin{verbatim}
Modify gas.eqn with a replacement file?
\end{verbatim}
Answer ``0'' unless you have written your own replacement file. More
information about the replacement feature can be found in the file
\verb|rpl/gas.rpl-example|.
\begin{verbatim}
Choose a selection number or type a boolean
expression:
\end{verbatim}
Now you can choose a subset of chemical reactions. A few predefined
standard selections are available. For all other purposes, a batch file
should be created, as explained at the end of this section. Some of the
predefined selections are:
\begin{description}
\item [EVAL:] A mechanism that was used for the evaluation of the MECCA
  chemistry in the global model ECHAM5/MESSy \citep{1851}.
\item [Minimum tropospheric chemistry:] A very small tropospheric
  mechanism.
\item [Minimum MBL chemistry:] A small mechanism that contains
  aqueous-phase chemistry and should only be used if the number of
  aerosol phases is $>0$.
\end{description}
For details about the selection, see Sect.~\ref{sec:selectreactions}.
\begin{verbatim}
Add Monte-Carlo factor to all rate
coefficients? [y/n/q, default=n]
\end{verbatim}
Answer ``\verb|n|'' here unless you want to perform Monte-Carlo
calculations, as described in Sect.~\ref{sec:montecarlo}.
\begin{verbatim}
Add diagnostic tracers to gas.eqn?
[q/0/?, default=0]
\end{verbatim}
Answer ``0''. Diagnostic tracers are usually only used for 3-dimensional
model runs.
\begin{verbatim}
Calculate accumulated reaction rates
of all equations? [y/n/q, default=n]
\end{verbatim}
Answer ``\verb|y|'' if you want to have all accumulated reaction rates
in the model output. Otherwise, answer ``\verb|n|''.
\begin{verbatim}
Tagging, doubling, both, or
none? [t/d/b/n/q, default=n]
\end{verbatim}
Please use the default ``\verb|n|'' here. The recommended way to test
tagging diagnostics and isotope modeling is described in
Sect.~\ref{sec:tagdbl}.
\begin{verbatim}
Run KPP?
\end{verbatim}
Answer ``\verb|y|''.
\begin{verbatim}
Choose an integrator [q=quit,
default=rosenbrock_posdef]:
\end{verbatim}
The default integrator is strongly recommended (see
Sect.~\ref{sec:selectintegrator} for details). Next, KPP will create
several Fortran95 files.
\begin{verbatim}
Remove indirect indexing with decomp?
[y/n/q, default=n]
\end{verbatim}
If this question shows up, answer ``\verb|n|''.
\begin{verbatim}
Create LaTeX listing of selected mechanism?
[y/n/q, default=n]
\end{verbatim}
If you answer ``\verb|y|'' here, a table of the current reaction
mechanism will be produced. Only the selected reactions will be listed.
The table also contains the rate coefficients and their references, as
described in Sect.~\ref{sec:latextable}.
\begin{verbatim}
Create graphviz plots of selected
mechanism? [y/n/q, default=n]
\end{verbatim}
If you have the ``\verb|dot|'' program from the graphviz software
installed, you can create graphical visualizations of the reaction
mechanism. As an example, the graphviz-generated plot of gas-phase
bromine chemistry is shown in Fig.~\ref{fig:Br}. For more information,
look at the files \verb|xgraphvizall|, \verb|xgraphviz|, and
\verb|spc_extract.awk| in the \verb|mecca/graphviz/| directory.
\begin{verbatim}
Do you want to delete the temporary
xmecca files?
\end{verbatim}
It is okay to delete these temporary files unless you need them for
debugging purposes.

When \verb|xmecca| finishes successfully, the Fortran95 code of your
selected mechanism has been created. The KPP-produced Fortran95 files
(Tab.~\ref{tab:files}) are moved into the \verb|mecca/smcl/| directory
(with lower-case names). An exception is
\verb|messy_mecca_kpp_Model.f90|, which is produced by KPP but not
needed for MECCA. The modular structure of the KPP-produced Fortran95
files is shown in Fig.~\ref{fig:kpp_use_diagr}.

If you need to create a chemical mechanism very often, it is quite
tedious to answer all questions every time. To make this easier, you can
copy the template \verb|batch/example.bat| to a new name (e.g.\
\verb|batch/myfile.bat|) and then enter your answers into that batch
file. Now you can create a new chemical mechanism in batch mode with
\begin{verbatim}
./xmecca myfile
\end{verbatim}
It is also possible to add the name of the batch file to the
\verb|xcaaba| command:
\begin{verbatim}
./xcaaba myfile
\end{verbatim}

\subsection{Selecting a set of chemical reactions}
\label{sec:selectreactions}

All chemical reactions are marked. Each marker consists of several
labels, which contain information about the domain
(troposphere/stratosphere), the phase where the reaction occurs
(gas/aqueous), its relevant chemical elements, and more. See
Sect.~\ref{sec:markerslabels} for a complete list of labels. To define a
set of chemical reactions, you can either choose a pre-defined selection
by number or enter a boolean expression based on the labels. Boolean
expressions are typed in gawk syntax. The most important operators and
expressions are:\\
\begin{tabular}{l@{ = }l}
\verb|&&| & AND\\
\verb!||! & OR\\
\verb|!|  & NOT\\
\verb|()| & parentheses\\
\verb|1|  & TRUE\\
\verb|0|  & FALSE
\end{tabular}

For example, to select all gas-phase reactions (G) except for those
including halogens (Cl, Br, I), type:

\verb|G && !Cl && !Br && !I|.

It is important to understand the logic behind this selection mechanism.
The expression ``\verb|Cl && Br|'' selects only those reactions that
contain chlorine {\em and} bromine. Similarly, the expression
``\verb|G && Het|'' selects only those reactions that occur in the gas
phase {\em and} are heterogeneous. However, since no reaction has both
the ``\verb|G|'' {\em and} the ``\verb|Het|'' label, this results in an
empty mechanism. If you want a mechanism that contains both gas-phase
and heterogeneous reactions, you must select all reactions that contain
either the label ``\verb|G|'' {\em or} the label ``\verb|Het|'', i.e.\
you must use the expression ``\verb!G || Het!''.

\subsection{Selecting a numerical integrator}
\label{sec:selectintegrator}

Several numerical integrators are defined in the subdirectory
\verb|mecca/kpp/int/| and can be used with KPP. The default is the
positive definite Rosenbrock solver with automatic time-step control
(\verb|rosenbrock_posdef|). It is very robust and capable of integrating
very stiff sets of equations (e.g.\ chemical mechanisms including both
gas- and aqueous-phase chemistry). Although a Rosenbrock solver with
manual time-step control (\verb|ros2_manual|) is also available, it is
strongly recommended not to use it for stiff sets of equations. If you
choose it, you do so at your own risk!

\section{Plotting the model results with the ferret software}
\label{sec:ferret}

If you have chosen netCDF output, you can plot the model results with
the ferret program (\url{http://ferret.wrc.noaa.gov/Ferret/}). Change
into the \verb|jnl/| directory, then start the program by typing
``\verb|ferret|''. When ferret has started, you can plot the gas-phase
species of the latest model simulation with the ferret script
\verb|xxxg.jnl| by typing:
\begin{verbatim}
go xxxg.jnl
\end{verbatim}
Similarly, \verb|xxxa.jnl| can be used to plot aqueous-phase species:
\begin{verbatim}
go xxxa.jnl
\end{verbatim}
The file \verb|xxxa.jnl| accepts several parameters to modify the plots.
The first parameter should be ``\verb|0d|'' for plotting box model
results. The second parameter can be set to ``\verb|mpl|'' or
``\verb|mpm|'' in order to plot either aqueous-phase concentrations in
\underline{m}ole \underline{p}er \underline{l}iter or mixing ratios in
\underline{m}ole \underline{p}er \underline{m}ole, respectively. The
third parameter defines the aerosol bin. With two aerosol bins,
``\verb|A01|'' refers to sulfate particles, and ``\verb|A02|'' to
sea-salt particles. For example, type:
\begin{verbatim}
go xxxa.jnl 0d mpl A02
\end{verbatim}
Photolysis rate coefficients can be plotted with \verb|jval.jnl|:
\begin{verbatim}
go jval.jnl
\end{verbatim}
If the calculation of accumulated reaction rates had been switched on in
\verb|xmecca| (see Sect.~~\ref{sec:xmecca}), plots of the reaction rates
can be made. One possibility is to plot all reactions with:
\begin{verbatim}
go rxnrates.jnl
\end{verbatim}
Alternatively, it is possible to plot only the production and
destruction rates for a certain species, e.g.\ for OH:
\begin{verbatim}
go rxnrates_scaled.jnl OH
\end{verbatim}
To plot results from previous simulations which are saved in the
\verb|output/| directory, edit the file \verb|setmodelrun.jnl| and enter
the paths of the directories in the ``\verb|GO _define_sensi|'' command.
To compare model simulations, you can enter two or more
``\verb|GO _define_sensi|'' commands in \verb|setmodelrun.jnl|. To plot
the difference between model simulations, activate the line
``\verb|DEFINE SYMBOL diffplot TRUE|'' in \verb|setmodelrun.jnl|.

\section{Run CAABA/MECCA in special modes}
\label{sec:special}

In the base configuration described so far, CAABA/MECCA calculates the
temporal evolution of the chemistry inside an air parcel. This is ideal
for sensitivity studies analyzing the effect of individual reactions
inside a large chemical mechanism. For other applications, some special
modes exist as described below.

\subsection{Multiple model simulations and steady-state}
\label{sec:multirun}

The so-called ``multirun'' mode performs multiple model simulations,
each of them terminating when a steady-state has been reached. This is
useful to calculate the steady-state concentrations of short-lived
species when the concentrations of longer-lived species (e.g.\
non-methane hydrocarbons) are known from measurements. The default
termination condition is that the relative change of \chem{OH} and
\chem{HO_2} between two model time steps is less than
$10^{-6}$~\unit{s^{-1}}. If necessary, this can be changed in the
function \verb|steady_state_reached| in \verb|messy_mecca.f90|. To avoid
that the concentrations of long-lived species change from their initial
values, they can be fixed in the file \verb|mecca/messy_mecca_kpp.kpp|
by adding them to the ``\verb|#SETFIX|'' line. A tracer initialization
file (see Sect.~\ref{sec:tracfile}) and a photolysis initialization file
(see Sect.~\ref{sec:readjfile}) must be available in the
\verb|multirun/input/| directory. As an example, the files
\verb|example_*.nc| are available. To create such input netCDF files
from ASCII files, the script \verb|asc2ferret4nc.tcsh| can be used.
Finally, since the multirun mode needs ``\verb|ncks|'' from the netCDF
Operators (NCO) software, it must be ensured that this program is
available. After these preparations, the multirun mode can be entered by
running \verb|xcaaba| and answering ``Choose a namelist file'' with
``\verb|m|''. The user can either select one input file or make model
simulations for all input files in the \verb|multirun/input/| directory.
For each input netCDF file, the script \verb|multirun.tcsh| is called.
For each time step contained in the input file, \verb|multirun.tcsh|
first creates a suitable namelist file \verb|caaba.nml| and then
performs a CAABA/MECCA model simulation. The namelist file contains
values for temperature and pressure taken from the input netCDF file. In
addition, the steady-state option is switched on with
``\verb|l_steady_state_stop = T|'', and the submodel READJ is activated
and used with ``\verb|USE_READJ = T|'' and
``\verb|photrat_channel = 'readj'|''. After all model simulations have
finished, a summary of the output is placed in the output directory. The
name of the output directory will be based on the name of the input
netCDF file, e.g.\ when the file \verb|example_small.nc| is used, the
output will be in \verb|output/multirun/example_small/|.

\subsection{Monte-Carlo}
\label{sec:montecarlo}

\subsubsection{Performing Monte-Carlo simulations}

In the Monte-Carlo mode, several CAABA/MECCA simulations are performed,
with each individual simulation using slightly different rate
coefficients. To activate it, you first have to create a new chemistry
mechanism with \verb|xmecca| (see Sect.~\ref{sec:xmecca}) and answer the
question ``Add Monte-Carlo factor to all rate coefficients'' with
``\verb|y|'' (e.g.\ using the batch file \verb|mcfct.bat|). This will
start the gawk script \verb|mcfct.awk|, which adds Monte-Carlo factors
to the rate coefficients in the equation file. Next, the \verb|xcaaba|
script can be used to start the simulations. A suitable namelist for
this purpose is \verb|caaba_mcfct.nml|. The script
\verb|montecarlo.tcsh| in the directory \verb|montecarlo/| will now
perform the Monte-Carlo model simulations. The default is to make 100
model simulations. To choose another value (up to 9999), change the
definition of \verb|maxline| in \verb|montecarlo.tcsh|.

\subsubsection{Analyzing Monte-Carlo simulations}

After performing the model simulations, the resulting netCDF files are
merged (using the tools \verb|ncpdq|, \verb|ncclamp|, and \verb|ncrcat|)
and then stored in the output directory with a name based on the date
and time of the model simulations, e.g.\ \verb|$outputdir| =
\verb|output/montecarlo/2010-06-24-16:29:00|. The final concentrations
and rate coefficients of all simulations are summarized in
\verb|caaba_mecca_c_end.nc| and \verb|caaba_mecca_k_end.nc|. Results of
the individual simulations can be found in the directories
\verb|$outputdir/runs/*|.

\paragraph{Time series}

If the model is set up to run for a fixed length (e.g.\ using the
default of \verb|runtime_str| = 8 days), the time series of all
simulations can be plotted together with ferret by activating the lines
for Monte-Carlo in \verb|setmodelrun.jnl|. However, these plots become
illegible if more than about 5 simulations are made.

\paragraph{Steady-state calculations}

The most efficient way to analyze a large number of Monte-Carlo
simulations is to use the steady-state option and only compare the final
values of the different model simulations, not the individual time
series. The ferret script \verb|montecarlo.jnl| can be used to create
scatter plots of concentrations vs rate coefficients. It also plots
linear regression lines for all comparisons above a given threshold of
the coefficient of determination ($r^2$). Note, however, that $r^2$ is
only an indicator of the goodness of a {\em linear} correlation. It is
also possible that the dependence of a species on a rate coefficient is
{\em nonlinear}. For an exhaustive analysis of the model results, the
threshold in the file \verb|scatterplot_mc.jnl| can be set to zero, thus
creating scatter plots of {\em all} species against {\em all} rate
coefficients.

\subsubsection{Variation of rate coefficients}

In each individual Monte-Carlo simulation $j$, all rate coefficients
$k_i$ are varied by a Monte-Carlo factor:

\begin{equation}
  k^{\rm MC}_{i,j} = k_i \times f_i^{x_{i,j}}
\end{equation}

Here, $k^{\rm MC}_{i,j}$ is the rate coefficient of reaction $i$ used in
the Monte-Carlo simulation $j$. It is defined as the product of the
recommended value $k_i$ and the Monte-Carlo factor $f_i^{x_{i,j}}$. This
Monte-Carlo factor consists of two parts, the uncertainty factor $f_i$
and the exponent $x_{i,j}$:

\paragraph{The uncertainty factor $f_i$}

The uncertainty factor $f_i$ describes the uncertainty of the measured
(or estimated) rate coefficient $k_i$. Its value can usually be found in
publications of laboratory studies or summaries like the JPL evaluation
\citep{1945}.

The tables of the IUPAC evaluations \egcite{1759} list the decadic
logarithm $\lg(f_i)$ of the uncertainty factor, which they call
``$\Delta\log k$''.

Sometimes an absolute uncertainty is quoted instead of an uncertainty
factor, e.g.\ $k = 2 \pm 0.2$ or $k = 2 \pm 10\%$. In this case we
define $f_i$ such that the upper limit is reached when multiplied with
$k_i$, i.e.\ in the current example $f_i = (2+0.2)/2 = 1.1$.

The uncertainty factor is defined in the equation files (\verb|*.eqn|)
in a comment starting with the paragraph symbol. Three different syntax
types are possible:
\begin{itemize}
\item If there is just one \verb|§| sign, (e.g.\ ``\verb|{§1.1}|''), the
  value inside the curly braces is the uncertainty factor $f_i$.
\item With two \verb|§| signs, (e.g.\ ``\verb|{§§0.04}|''), the value
  inside the curly braces equals $\lg(f_i)$.
\item If there is only a \verb|§| sign (``\verb|{§}|'') but no number,
  the uncertainty factor is set to the default value of $f_i$~= 1.25.
\end{itemize}

\paragraph{The Monte-Carlo exponent $x_{i,j}$}

\begin{figure}
  \begin{center}
    \includegraphics[width=\columnwidth]{random}
  \end{center}
  \caption{Example histogram of normally-distributed random numbers.}
  \label{fig:random}
\end{figure}

There is one Monte-Carlo exponent $x_{i,j}$ (variable ``\verb|mcexp|''
in the code) for each rate coefficient $k_i$ and for each individual
Monte-Carlo simulation $j$. The values of $x_{i,j}$ are
normally-distributed random numbers centered around zero (see
Fig.~\ref{fig:random}), and produced with the Marsaglia polar method
(\url{http://en.wikipedia.org/wiki/Marsaglia_polar_method}). As input
for the Marsaglia polar method, uniformly distributed random numbers
between 0 and 1 calculated with either the standard Fortran95 function
{\tt RANDOM\_NUMBER} or the Mersenne Twister algorithm \citep{2404} are
used.

\subsubsection{Changing the uncertainty factors}

The uncertainty factors can be changed by modifying the equation files,
as shown in Sect.~\ref{sec:eqnfiles}. Note that predefined rate
coefficients (e.g.\ \verb|k_HO2_HO2|) already contain an uncertainty
factor and there must not be an additional factor in the reaction where
they are used.

In some cases, it may be useful to vary only one or a few rate
coefficients. To do this, it is first necessary to find the correct
indices of \verb|mcexp(...)| in \verb|mecca.eqn| (note that these
indices may change when creating a new mechanism with \verb|xmecca|). As
an example, to vary only the rate coefficients that use \verb|mcexp(40)|
and \verb|mcexp(50)|, the following lines can be added to subroutine
\verb|mecca_init| in \verb|messy_mecca_box.f90| after
\verb|CALL define_mcexp|:
\begin{verbatim}
DO i=1, MAX_MCEXP
  IF ((i/=40).OR.(i/=50)) mcexp(i) = 0.
ENDDO
\end{verbatim}

To verify that the rate coefficients are modified in the Monte-Carlo
simulations, it is possible to temporarily activate the subroutine
\verb|montecarlo_check| in \verb|template_messy_mecca_kpp.f90| and check
the output in \verb|caaba.log|. After these tests,
\verb|montecarlo_check| must be switched off again to allow normal model
simulations.

\subsection{Lagrangian trajectories}
\label{sec:lagrangian}

CAABA can be used as a Lagrangian trajectory box model \citep{2403}. The
usual combination of submodels for this purpose includes the CAABA
submodel TRAJECT for the processing of trajectory information, MECCA for
atmospheric chemistry, and JVAL for photolysis rate calculation. All
important settings for trajectory calculations can be made via the
namelist file plus a few external files.

\subsubsection{Namelist parameters}
\label{sec:trajnml}

A namelist template can be found in
\verb|nml/caaba_traject_example.nml|. After copying the namelist file to
a new name in the same directory and altering the settings, it will be
available when running CAABA via \verb|xcaaba|. There are standard and
trajectory-exclusive namelist parameters to be set:
\begin{itemize}\nosep
\item \noindent\textbf{Submodel switches (mandatory):} The trajectory
  mode of CAABA requires that \verb|USE_TRAJECT = T|,
  \verb|USE_MECCA = T|, and \verb|USE_JVAL = T| and/or
  \verb|USE_SAPPHO = T|. One of the two photolysis rate models is
  sufficient, but it is also possible to run them in parallel (see also
  \verb|photrat_channel|). If the use of external photolysis rates is
  desired, \verb|USE_JVAL = T| is mandatory since the scaling of
  photolysis rate coefficients is currently only implemented for
  \verb|JVAL|. For the application in
  \citet{2403}, the submodel \verb|SEMIDEP| was switched off.
\item \textbf{Scenarios (optional):} The use of scenarios is optional in
  trajectory mode. When using external input for chemical initialization
  and photolysis rates, they can be used to supplement default values
  for unspecified photolysis rate coefficients or initial mixing ratios.
\item \textbf{Photolysis rate channel:} Since the scaling of
  photolysis rate coefficients is currently only implemented for
  \verb|JVAL|, choose
  \verb|photrat_channel = 'jval'| when planning to prescribe photolysis
  rates.
\item \textbf{Trajectory input (mandatory):} The path to the netCDF file
  containing trajectory waypoints should be specified as
  \verb|input_physc = 'traject/example_traj.nc'|. For its structure, see
  section \ref{sec:trajfile}.
\item \textbf{Tracer initialization (optional):} Tracer mixing ratios
  can be initialized with an external netCDF file by specifying the path
  to it: \verb|init_spec = 'traject/example_init.nc'|. For its
  structure, please refer to section \ref{sec:tracfile}. Without an
  external file for tracer initialization, tracers mixing ratios are
  initialized by box model default or by the chosen scenario.
\item \textbf{Photolysis rates (optional):} To scale photolysis rates to
  prescribed $J$(\chem{NO_2}) values, specify the path to the respective
  file as \verb|input_jval = 'traject/example_jval.nc'|. For its
  structure, see \ref{sec:jvalfile}.
\item \textbf{Variable names (partly mandatory):} There are default
  trajectory variable names designated in CAABA. They can be selectively
  changed by providing alternative variable names. Here is a list of
  trajectory variables, their default name, and respective examples how
  to specify an alternative variable name:
  \begin{center}
    \begin{tabular}{lll}
      \hline
      quantity       & default       & alternative\\
      \hline
      longitude      & \verb|LON|    & \verb|vlon = 'LON_TR'|\\
      latitude       & \verb|LAT|    & \verb|vlat = 'LAT_TR'|\\
      pressure       & \verb|PRESS|  & \verb|vpress = 'P'|\\
      temperature    & \verb|TEMP|   & \verb|vtemp = 'TM1'|\\
      rel. humidity  & \verb|RELHUM| & \verb|vrelhum = 'rh'|\\
      spec. humidity &               & \verb|vspechum = 'sh'|\\
      \hline
    \end{tabular}
  \end{center}
  For only one of the two humidity variables, an alternate name may be
  active in the namelist. This is the variable that is subsequently used
  as humidity and that is internally converted to relative humidity if
  necessary. When specific humidity is provided, then both specific
  humidity and relative humidity are written to output
  \verb|caaba_physc.nc|, since CAABA uses relative humidity internally.
  When relative humidity is provided, only relative humidity will be
  written to output.
\item \textbf{Humidity definitions (optional):} To switch from the
  traditional relative humidity definition to the WMO definition,
  \verb|l_relhum_wmo = T| can be selected (see also
  Sect.~\ref{sec:trajout}).
\item \textbf{Use relative humidity? (optional):} By default, the
  relative humidity is used to calculate the concentration of water
  vapor $c$(\chem{H_2O}). If \verb|l_ignore_relhum = T|, then the
  relative humidity is ignored, and $c$(\chem{H_2O}) must be initialized
  in \verb|SUBROUTINE x0| (either directly or via external chemical
  initialization).
\item \textbf{Water vapor saturation pressure (optional):} By default,
  the parameterizations from \citet{2553} for both, vapor pressure over
  liquid and over ice are used. To use the saturation vapor pressure
  parameterizations from EMAC, choose the namelist entry
  \verb|l_psat_emac = T|,
\item \textbf{Integration time (optional):} \verb|timesteplen_str| sets
  the integration and output time step. See also
  Sect.~\ref{sec:caaba_nml}.
\item \textbf{Selecting trajectory parts (optional):} Two namelist
  parameters allow flexible selection of the model runtime along the
  trajectory. \verb|runlast| defines the start of the CAABA simulation
  counted backwards in time from the trajectory end, i.e.
  \verb|runlast = 4.5| means ``calculate the last 4.5 days of the
  trajectory''. The unit is days. The parameter \verb|runtime_str|
  defines the overall model simulation time. Thus, \verb|runlast| and
  \verb|runtime_str| combined select any desired part of the
  trajectory.
\end{itemize}

\subsubsection{Trajectory input file}
\label{sec:trajfile}

The trajectory information is provided to CAABA via an external netCDF
file specified in the namelist by \verb|input_physc|. A sample file is
available at \verb|traject/example_traj.nc|. The file should contain a
time origin in 'seconds/minutes/hours/days since yyyy-mm-dd hh:mm:ss',
where the seconds in the time string are optional, for example:
\verb|"MINUTES since 2000-01-19 08:00:00"|. The file must contain at
least two trajectory waypoints and the following time-dependent
variables:
\begin{center}
  \begin{tabular}{lll}
    \hline
    quantity       & default name &  unit\\
    \hline
    longitude      & \verb|LON|   &  degrees east\\
    latitude       & \verb|LAT|   &  degrees north\\
    pressure       & \verb|PRESS| &  Pascal\\
    temperature    & \verb|TEMP|  &  Kelvin\\
    (rel. humidity)  &            &  1\\
    (spec. humidity) &            &  kg/kg\\
    \hline
  \end{tabular}
\end{center}
Of the two humidity quantities, only one needs to be present. If
specific humidity is given, it is converted to relative humidity (RH)
within CAABA. Depending on the namelist parameter \verb|l_relhum_wmo|,
either the traditional definition of relative humidity
\begin{equation}
 \mbox{RH} = \frac{p_{H_2O}(T)}{p_{\rm sat}(T)},
\end{equation}
or the WMO definition \citep{1660}
\begin{equation}
\mbox{RH} = \frac{\omega_v}{\omega_{vs}} =
\frac{p_{\chem{H_2O}}(T)}{p - p_{\chem{H_2O}}(T)} \times
\frac{p - p_{\rm sat}(T)}{p_{\rm sat}(T)}
\end{equation}
is assumed, with p = atmospheric pressure, $p_{\chem{H_2O}}$ = water
vapor partial pressure, $p_{\rm sat}$ = water vapor saturation pressure,
$\omega_v$ = water vapor mass mixing ratio, and $\omega_{vs}$ =
saturation water vapor mass mixing ratio in dry air. We calculate the
latter two as:
\begin{eqnarray}
  \omega_v & = & \frac{q}{1-q}\\
  \omega_{vs} & = & \frac{M(\chem{H_2O})}{M({\rm air})} \times
  \frac{p_{\rm sat}(T)}{p - p_{\rm sat}(T)},
\end{eqnarray}
using the temperature-dependent saturation water vapor pressure $p_{\rm
  sat}(T)$, and the molar masses $M$ of water and dry air. There are two
parameterization schemes available for the water vapor saturation
pressure, which can be selected by namelist (see
Sect.~\ref{sec:trajnml}).

\subsubsection{Photolysis rate file}
\label{sec:jvalfile}

It is possible to scale photolysis rate coefficients via prescribed
$J$(\chem{NO_2}) values from a netCDF file. An example is available at
\verb|traject/example_jval.nc|. The name in the netCDF file must be
\verb|J_NO2|. The files specified in \verb|input_jval| (e.g.\
\verb|example_jval.nc|) and \verb|input_physc| (e.g.\
\verb|example_traj.nc|) must both refer to exactly the same trajectory
because the $J$ values are read into the model at the same times as
other trajectory information.

\subsubsection{Trajectory mode output}
\label{sec:trajout}

Output along the trajectory is written to \verb|caaba_messy.nc|. There
are some special variables written out in addition to the default
\verb|caaba_messy.nc| output. They are listed in the second part of
Tab.~\ref{tab:trajout}.
\begin{table*}%[htb]
  \begin{center}
    \caption{List of trajectory mode output variables.}
    \label{tab:trajout}
\begin{tabular}{lll}
  \hline
  variable         & unit             & notes\\
  \hline
  \verb|lon|       & \multicolumn{2}{l}{(dummy 3-D x-coordinate)}\\
  \verb|lat|       & \multicolumn{2}{l}{(dummy 3-D y-coordinate)}\\
  \verb|lev|       & \multicolumn{2}{l}{(dummy 3-D z-coordinate)}\\
  \verb|time|      & \textless time unit\textgreater\ since \textless time origin\textgreater & UTC\\
  \verb|lon_tr|    & degrees east     & longitude\\
  \verb|lat_tr|    & degrees north    & latitude\\
  \verb|press|     & Pa               & pressure\\
  \verb|temp|      & K                & temperature\\
  \verb|relhum|    & 1                & relative humidity (RH)\\
  \hline
  \verb|spechum|   & kg/kg            & specific humidity ($q$)\\
  \verb|sza|       & deg              & solar zenith angle\\
  \verb|J_NO2_x|   & 1/s              & \chem{NO_2} photolysis rate\\
  \verb|localtime| & (same as time)   & local time at trajectory position\\
  \verb|year_loc|  &                  & year of local time\\
  \verb|month_loc| &                  & month of local time\\
  \verb|day_loc|   &                  & day of local time\\
  \verb|hour_loc|  &                  & hour of local time\\
  \verb|min_loc|   &                  & minute of local time\\
  \verb|sec_loc|   &                  & second of local time\\
  \hline
\end{tabular}
\end{center}
\end{table*}

Specific humidity (q) is only written to output if it was provided
as input.

\begin{table*}%[htb]
  \begin{center}
    \caption{List of CAABA/MECCA Fortran95 files}
    \label{tab:files}
    \begin{tabular}{lp{0.55\textwidth}}
      \hline
      \multicolumn{2}{l}{CAABA box model related files}\\
      \hline
      \verb|caaba.f90|                         & main box model file\\
      \verb|caaba_io.f90|                      & input/output\\
      \verb|caaba_io_netcdf.inc|               & netCDF input/output\\
      \verb|caaba_io_ascii.inc|                & ASCII input/output\\
      \verb|caaba_mem.f90|                     & declaration of CAABA variables\\
      \verb|messy_main_control_cb.f90|         & flow control\\
      %\verb|messy_e4chem_box.f90|              & connection of E4CHEM to CAABA\\
      \verb|messy_jval_box.f90|                & connection of JVAL to CAABA\\
      \verb|messy_mecca_box.f90|               & connection of MECCA to CAABA\\
      \verb|messy_mecca_dbl_box.f90|           & for isotope modeling\\
      \verb|messy_mecca_tag_box.f90|           & for isotope modeling\\
      \verb|messy_readj_box.f90|               & connection of READJ to CAABA\\
      \verb|messy_sappho_box.f90|              & connection of SAPPHO to CAABA\\
      \verb|messy_semidep_box.f90|             & simplified emission and
                                                 deposition, including
                                                 connection of SEMIDEP to CAABA\\
      \verb|messy_traject_box.f90|             & trajectory calculations\\
      \hline
      \multicolumn{2}{l}{static core files (generic)}\\
      \hline
      \verb|messy_cmn_gasaq.f90|               & common definitions for
                                                 gas/aqueous mass transfer\\
      \verb|messy_cmn_photol_mem.f90|          & common definitions for photolysis\\
      \verb|messy_main_blather.f90|            & print utilities\\
      \verb|messy_main_constants_mem.f90|      & physical constants\\
      \verb|messy_main_rnd.f90|                & random number generation\\
      \verb|messy_main_rnd_lux.f90|            & (file exists but is not
                                                 not used with CAABA)\\
      \verb|messy_main_rnd_mtw.f90|            & Mersenne-Twister random numbers\\
      \verb|messy_main_timer.f90|              & timer\\
      \verb|messy_main_tools.f90|              & auxiliary functions\\
      \verb|messy_main_tools_kp4_compress.f90| & (file exists but is not
                                                 not used with CAABA)\\
      \hline
      \multicolumn{2}{l}{static core files (submodel-specific)}\\
      \hline
      %\verb|messy_e4chem.f90|                  & \\
      \verb|messy_jval.f90|                    & calculation of J-values\\
      \verb|messy_jval_jvpp.inc|               & include file for JVAL\\
      \verb|messy_readj.f90|                   & read J-values\\
      \verb|messy_sappho.f90|                  & simplified and parameterized
                                                 photolysis rate coefficients\\
      \hline
      \multicolumn{2}{l}{static MECCA core files}\\
      \hline
      \verb|messy_mecca.f90|                   & MECCA core\\
      \verb|messy_mecca_aero.f90|              & aerosol chemistry\\
      \verb|messy_mecca_khet.f90|              & (file exists but is not
                                                 used with CAABA)\\
      \verb|messy_mecca_dbl_parameters.inc|    & for isotope modeling\\
      \verb|messy_mecca_tag_parameters.inc|    & for isotope modeling\\
      \hline
      \multicolumn{2}{l}{KPP- and {\tt xmecca}-produced files}\\
      \hline
      \verb|messy_mecca_kpp.f90|               & a wrapper for the KPP files\\
      \verb|messy_mecca_kpp_function.f90|      & ODE function\\
      \verb|messy_mecca_kpp_global.f90|        & global data headers\\
      \verb|messy_mecca_kpp_initialize.f90|    & initialization\\
      \verb|messy_mecca_kpp_integrator.f90|    & numerical integration\\
      \verb|messy_mecca_kpp_jacobian.f90|      & ODE Jacobian\\
      \verb|messy_mecca_kpp_jacobiansp.f90|    & Jacobian sparsity\\
      \verb|messy_mecca_kpp_linearalgebra.f90| & sparse linear algebra\\
      \verb|messy_mecca_kpp_monitor.f90|       & equation info\\
      \verb|messy_mecca_kpp_parameters.f90|    & model parameters\\
      \verb|messy_mecca_kpp_precision.f90|     & arithmetic precision\\
      \verb|messy_mecca_kpp_rates.f90|         & user-defined rate laws\\
      \verb|messy_mecca_kpp_util.f90|          & utility input-output\\
      \hline
    \end{tabular}
  \end{center}
\end{table*}

\subsection{Tagging diagnostics and isotope modeling}
\label{sec:tagdbl}

Before using the tagging feature, please make sure that version 2.2.0 or
higher of \verb|fpc| (Free Pascal compiler,
\url{http://www.freepascal.org/}) is installed and configured on your
system. To create an isotopically tagged chemical mechanism with
\verb|xmecca|, the batch file \verb|iso_example.bat| can be used. Note
that \verb|tagdbl| is set to ``\verb|d|'' (doubling) here. The
alternative option ``\verb|t|'' (tagging) is not intended to be used in
box-model simulations and thus currently disabled. Before running
CAABA/MECCA, the execution of the doubling code has to be enabled by
setting \verb|l_dbl=T| in the \verb|&CTRL| namelist in \verb|mecca.nml|.
In this configuration, the model simulation will produce additional
files \verb|caaba_mecca_dbl_*.nc| containing output from the doubling,
e.g., concentrations of isotopologues and isotope ratios. The batch file
\verb|iso_example.bat| contains a simple example of the
\chem{^{12}C}/\chem{^{13}C} isotopes of methane tagging. Further tagging
configurations can be found by running \verb|xmecca| interactively
(without using a batch file) and answering ``\verb|d|'' to the question
about doubling/tagging. To obtain further details, please contact Sergey
Gromov \url{<sergey.gromov@mpic.de>}.

\section{Modifying CAABA/MECCA}
\label{sec:modifying}

The CAABA/MECCA model simulation can be modified by changing the input
files (\verb|*.nc|), the namelist files (\verb|*.nml|), the species
files (\verb|*.spc|), the equation files (\verb|*.eqn|) and the
Fortran95 files (\verb|*.f90|). Frequently applied changes are:
\begin{itemize}
\item Change the model time (start, stop), location (lon, lat), and/or
  meteorology ($p$, $T$, RH) $\TO$ Sect.~\ref{sec:caaba_nml}
\item Add, delete, and change individual chemical reactions $\TO$
  Sect.~\ref{sec:expandchem}
\item Create or change the set of photolysis rate coefficients $\TO$
  Sect.~\ref{sec:scenarios}
\item Create or change the set of emission fluxes $\TO$
  Sect.~\ref{sec:scenarios}
\item Create or change the set of deposition velocities $\TO$
  Sect.~\ref{sec:scenarios}
\item Create or change the set of initial mixing ratios $\TO$
  Sect.~\ref{sec:scenarios}
\item Change the model initialization using input files $\TO$
  Sect.~\ref{sec:inputfiles}
\end{itemize}
The following sections describe typical (mostly minor) changes to the
model that are probably needed by many users. If you would like to make
extensive changes/additions to the model, please follow also the general
programming guidelines described in Sect.~\ref{sec:guidelines}.

\subsection{Namelist files}
\label{sec:nmlfiles}

Fortran95 namelist files allow modifications of the model simulation
without having to recompile the source code.

\subsubsection{The CAABA namelist file {\tt caaba.nml}}
\label{sec:caaba_nml}

The file \verb|caaba.nml| contains the namelist \verb|&CAABA|. Here
individual parts of the CAABA model (the so-called ``MESSy submodels'')
can be switched on or off. It is important that the following switches
are set to ``\verb|T|'' (=true):
\begin{verbatim}
USE_MECCA   = T
USE_SAPPHO  = T
USE_SEMIDEP = T
\end{verbatim}
To use the photolysis rate coefficients from SAPPHO in MECCA, set:
\begin{verbatim}
photrat_channel = 'sappho'
\end{verbatim}
Alternatively, you can switch on the JVAL submodel with
\verb|USE_JVAL = T| and then select \verb|photrat_channel = 'jval'|. It
is fine to switch on both the JVAL and the SAPPHO submodel, which can be
useful for a comparison. However, only the values selected by
\verb|photrat_channel| are used for the MECCA chemistry.

You can define the model start, runtime, and time step, e.g.:
\begin{verbatim}
model_start_day = 90.
runtime_str     = '10 days'
timesteplen_str = '15 minutes'
\end{verbatim}
If you don't set these, the default is a model start on day-of-year
(``Julian Day'') 80, a model simulation duration of 8 days, and an
output time step of 20 minutes. The \verb|timesteplen_str| value can be
given as any integer or floating point and in the units of
\verb|seconds|, \verb|minutes|, or \verb|hours|.

As an alternative, it is possible to stop the model simulation when a
steady state has been reached. This is normally used in the ``multirun''
mode (Sect.~\ref{sec:multirun}):
\begin{verbatim}
l_steady_state_stop = T
\end{verbatim}
The default location of the model simulation (latitude, longitude) is
45$^\circ$ N and 0$^\circ$ E. It can be changed with e.g.:
\begin{verbatim}
degree_lat = 82.  ! Alert
degree_lon = 297. ! Canada
\end{verbatim}
Changing the latitude only affects photolysis calculations (via the
zenith angle calculations). The longitude is only used for Lagrangian
trajectory calculations.

The values of temperature (\verb|temp|), pressure (\verb|press|),
relative humidity (\verb|relhum|), and the height of the boundary layer
(\verb|zmbl|) can be defined, e.g.:
\begin{verbatim}
temp   = 293.    ! [K]
press  = 101325. ! [Pa]
relhum = 0.81    ! [1]
zmbl   = 1000.   ! [m]
\end{verbatim}
The values shown here are the default values as defined in
\verb|caaba_mem.f90|.

In the submodel SEMIDEP, emissions are distributed evenly into the
well-mixed boundary layer of height \verb|zmbl|. Note that CAABA is a
box model, and changing \verb|zmbl| has no other effects than this.

It is possible to initialize chemical species from a netCDF file (see
Sect.~\ref{sec:tracfile}). To activate this feature, define a valid
input file name, e.g.:
\begin{verbatim}
init_spec = 'traject/example_init.nc'
\end{verbatim}
If the submodel READJ is switched on, a netCDF file containing J-values
must be defined (see Sect.~\ref{sec:readjfile}). In addition, an index
can be defined, if the netCDF file contains data for more than one time
step, e.g.:
\begin{verbatim}
init_j = 'multirun/input/example_small.nc'
init_j_index = 25
\end{verbatim}

To facilitate running CAABA under different boundary conditions,
so-called ``scenarios'' can be defined. Currently, there are scenarios
for photolysis (\verb|photo_scenario|), initialization
(\verb|init_scenario|), emission (\verb|emission_scenario|), and dry
deposition (\verb|drydep_scenario|):
\begin{verbatim}
photo_scenario    = 'MBL'
init_scenario     = 'MBL'
emission_scenario = 'MBL'
drydep_scenario   = 'MBL'
\end{verbatim}
Possible values (``\verb|MBL|'' is used here as an example) for the
scenarios can be found in the variable \verb|list_of_scenarios| in
\verb|caaba.f90|. Note that values have not yet been assigned for all
scenarios. See Sect.~\ref{sec:scenarios} for details about scenarios.

Finally, it is possible to skip the chemistry calculations with KPP
completely. This is only useful for debugging:
\begin{verbatim}
l_skipkpp = T
\end{verbatim}

\subsubsection{The MECCA namelist file {\tt mecca.nml}}
\label{sec:mecca_nml}

The file \verb|mecca.nml| contains the namelists \verb|&CTRL_KPP| and
\verb|&CTRL| (the namelist \verb|&CPL| is not used in connection with
CAABA). \verb|&CTRL_KPP| is used for fine-tuning the numerical
integration. The default selection \verb|icntrl(3) = 2| should normally
be suitable.

\subsubsection{The JVAL namelist file {\tt jval.nml}}
\label{sec:jval_nml}

The file {\tt jval.nml} contains several namelists. When JVAL is run as
a submodel for CAABA, only the control namelist (\verb|&CTRL|) is used.
Normally, there is no need to change the default settings.

\subsection{The species files {\tt gas.spc} and {\tt aqueous.spc}}
\label{sec:spcfiles}

The files \verb|*.spc| declare chemical species for KPP. All species
that may occur in an equation must be declared here. Additional dummy
species may also be declared here.

Gas-phase species are declared in \verb|gas.spc|. Examples for gas-phase
species are \verb|O2|, \verb|O1D|, and \verb|NO2|. The names of lumped
species start with ``\verb|L|''. The names of species with two or more
carbon atoms are taken from the master chemical mechanism (MCM).

MECCA also includes aqueous species which are declared in
\verb|aqueous.spc|. The names of cations end with ``\verb|p|'' for plus.
The names of single-charge anions end with ``\verb|m|'' for minus.
Doubly-charged anions end with ``\verb|mm|''. Examples for aqueous
species are ``\verb|H2O2|'' (\chem{H_2O_2}), ``\verb|Hp|'' (\chem{H^+}),
``\verb|NO3m|'' (\chem{NO_3^-}), and ``\verb|SO4mm|''
(\chem{SO_4^{2-}}).

All aqueous-phase species have the suffix ``\verb|_a##|'', which is a
placeholder for the aerosol phase. \verb|xmecca| replaces it by a
two-digit number, e.g., ``\verb|_a01|'' for the accumulation or
``\verb|_a02|'' for the coarse mode. This allows separate chemistry
calculations for aerosol particles of different size and composition.

All species are defined here with \verb|#DEFVAR|, i.e.\ KPP considers
them as prognostic variables. To treat a species as a constant (e.g.\
\chem{CO_2}), it can be added to the \verb|#SETFIX| command in the file
\verb|messy_mecca_kpp.kpp|. The declarations in \verb|#DEFVAR| can and
should remain unchanged when altering \verb|#SETFIX|.

\subsection{The equation files {\tt gas.eqn} and {\tt aqueous.eqn}}
\label{sec:eqnfiles}

The equation files \verb|*.eqn| define the chemical reaction mechanism
for KPP. After making any changes to the equation files, it is always
necessary to execute KPP via \verb|xmecca| again. Each reaction occupies
one line in this file. An example is:
\begin{verbatim}
<G1000> O2 + O1D = O3P + O2 : {%StTrG}
   3.3E-11*EXP(55./temp); {&1945}{§1.1}
\end{verbatim}
Note that, although split here in this manual, this is only 1 line in
the equation file. The line starts with the reaction number, which is
enclosed in angle brackets ``\verb|<...>|'' (see
Sect.~\ref{sec:rxnnumbers}). The second part (up to the colon) defines
the reaction, and the third part (between the colon and the semicolon)
defines the rate coefficient. The lines may also contain comments.
Comments in equation files are either enclosed in curly braces, or the
comment line starts with ``\verb|//|''. When using \verb|xmecca|, some
comments have a special meaning. Comments starting with the percent
symbol ``\verb|{%...}|''
  are markers (see Sect.~\ref{sec:markerslabels}). Comments starting
  with the ampersand ``\verb|{&...}|'', the ``at''-symbol
  ``\verb|{@...}|'', or the dollar ``\verb|{$...}|'' are used to store
  information for the listing of reactions, as explained in
  Sect.~\ref{sec:latextable}. Comments starting with the paragraph
  symbol ``\verb|{§...}|'' are defining uncertainties of the rate
  coefficients for Monte-Carlo calculations (see
  Sect.~\ref{sec:montecarlo}).

If the definition of a rate coefficient is very complex, it can be
stored in a Fortran95 variable and the variable is put into the
\verb|gas.eqn| file. For example, the rate of the self reaction of
\chem{HO_2} is quite complex since it depends on humidity. It is
predefined and the reaction line can be simplified to:
\begin{verbatim}
HO2 + HO2 = H2O2 : k_HO2_HO2;
\end{verbatim}
The declaration and definition of \verb|k_HO2_HO2| are also in the
\verb|gas.eqn| file. They can be found in the so-called KPP ``inline
types'' \verb|F90_GLOBAL| and \verb|F90_RCONST|, e.g.:
\begin{verbatim}
#INLINE F90_GLOBAL
  REAL :: k_HO2_HO2
#ENDINLINE

#INLINE F90_RCONST
  k_HO2_HO2 = (1.5E-12*EXP(19./temp)+ &
              1.7E-33*EXP(1000./temp)*cair)* &
              (1.+1.4E-21* &
              EXP(2200./temp)*C(KPP_H2O))
#ENDINLINE
\end{verbatim}
Another method to add reaction rates with complex dependencies are
Fortran95 functions. This is done for example for the oxidation of
\chem{S(IV)} by \chem{H_2O_2} (\verb|k_SIV_H2O2|). A function call is
given as the rate expression in the \verb|*.eqn| file. These functions
are defined with the ``inline type'' \verb|F90_RATES|:
\begin{verbatim}
#INLINE F90_RATES
  ELEMENTAL REAL(dp) FUNCTION k_SIV_H2O2 &
    (k_298,tdep,cHp,temp)
    ...
  END FUNCTION k_SIV_H2O2
#ENDINLINE
\end{verbatim}

\subsubsection{Reaction numbers}
\label{sec:rxnnumbers}

Each reaction in an equation file has a unique reaction ``number''
(number is not quite correct, since letters are included as well), which
is enclosed in angle brackets, e.g.: ``\verb|<G1000>|''. The reaction
number starts with one or more upper case letters denoting the type of
reaction. The following types exist:

\begin{tabular}{lp{0.8\columnwidth}}
  A   & aqueous-phase reactions\\
  H   & Henry's law (dissolution and evaporation)\\
  EQ  & equilibria in the aqueous phase (forward and backward reactions of
  acid/base and other equilibria)\\
  G   & gas-phase reactions\\
  J   & J-values of photolysis reactions\\
  HET & heterogeneous reactions (e.g.\ on polar stratospheric clouds)
\end{tabular}

The type is followed by a sequence of 3 or 4 digits. The first digit is
the number of the main element of the reaction. The following numbers
are used:

\begin{tabular}{lll}
  1) & \chem{O}  & Oxygen   \\
  2) & \chem{H}  & Hydrogen \\
  3) & \chem{N}  & Nitrogen \\
  4) & \chem{C}  & Carbon   \\
  5) & \chem{F}  & Fluorine \\
  6) & \chem{Cl} & Chlorine \\
  7) & \chem{Br} & Bromine  \\
  8) & \chem{I}  & Iodine   \\
  9) & \chem{S}  & Sulfur   \\
 10) & \chem{Hg} & Mercury  \\
\end{tabular}

Out of those elements that occur in a reaction, the one with the highest
number is called the main element. Accordingly, the second digit is
determined by the element with the second highest number (or set to zero
if there is no second element in the reaction). There is one exception
in this numbering scheme: For the carbon group, the second digit is the
number of C atoms in the largest organic molecule.

The following digits have no special meaning. If a reaction branches
into several pathways, a suffix ``\verb|a|'', ``\verb|b|'',
``\verb|c|'', \dots\ is added.

\subsubsection{Markers and labels}
\label{sec:markerslabels}

Each reaction must contain a marker. A marker contains several
case-sensitive labels. The syntax is ``\verb|{%...}|'' where the dots
represent the labels. Labels are used to select specific reactions, as
described above (Sect.~\ref{sec:selectreactions}). The labels are
placed in the marker without separators. The following labels are
available and should appear in this order:
\begin{enumerate}
\item the domain, i.e.\ altitudes at which the reaction occurs
  (mandatory, include at least one)\\
  \begin{tabular}{l@{ = }p{0.7\columnwidth}}
  \verb|St|   & Reactions relevant in the stratosphere\\
  \verb|Tr|   & Reactions relevant in the troposphere
  \end{tabular}
\item phase (mandatory, include exactly one)\\
  \begin{tabular}{l@{ = }p{0.7\columnwidth}}
    \verb|Aa##| & Aqueous, aerosol (\verb|##| is a placeholder for the
    2-digit aerosol phase number)\\
    \verb|Ara|  & Aqueous, rain (not used for the box model)\\
    \verb|G|    & Gas phase reactions\\
    \verb|Het|  & heterogeneous reactions (e.g.\ on polar stratospheric
                  clouds)\\
  \end{tabular}
\item elements (include all elements that occur in the reaction, except
  for H and O)\\
  \begin{tabular}{l@{ = }p{0.7\columnwidth}}
  \verb|N|    & Nitrogen\\
  \verb|C|    & Carbon with $>$ 1 C atom (only used for C,N,O species but
                not for halogenated or sulfur-containing organics)\\
  \verb|F|    & Fluorine\\
  \verb|Cl|   & Chlorine\\
  \verb|Br|   & Bromine\\
  \verb|I|    & Iodine\\
  \verb|S|    & Sulfur\\
  \verb|Hg|   & Mercury
  \end{tabular}
\item other\\
  \begin{tabular}{l@{ = }p{0.7\columnwidth}}
  \verb|J|    & Photolysis reactions\\
  \verb|Mbl|  & Minimum reaction mechanism for MBL chemistry\\
  \verb|Sc|   & Scavenging chemistry mechanism\\
  \verb|Scm|  & Scavenging chemistry mechanism, minimum selection
  \end{tabular}
\end{enumerate}

See Sect.~\ref{sec:newlabel} for a description how to add new labels to
\verb|xmecca|.

\subsubsection{Creating a table of the chemical mechanism}
\label{sec:latextable}

To ensure that the documentation of the chemical mechanism is always up
to date, the necessary information is contained inside the species and
equation files. If you have the programs pdfLa\TeX\ and Bib\TeX\
installed on your system, you can generate a table of the chemical
mechanism in pdf format.

The gawk scripts \verb|spc2tex.awk| and \verb|eqn2tex.awk| convert
information from the selected reactions into a La\TeX\ table. Bib\TeX\
citations are included in comments starting with an ampersand
``\verb|{&...}|''. If there is a second ampersand ``\verb|{&&...}|'',
additional information about reactions can be found in
\verb|meccanism.tex| as a footnote to the tables. Comments starting with
the ``at'' symbol ``\verb|{@...}|'' or the dollar ``\verb|{$...}|'' can
be used to put La\TeX\ commands directly into the \verb|*.eqn| files.
\verb|eqn2tex.awk| produces several \verb|*.tex| files which are
included into \verb|meccanism.tex|.

\begin{figure*}%[htb]
  \begin{center}
    \includegraphics[width=\textwidth]{caaba_use_diagr}
  \end{center}
  \caption{Module structure of MECCA when it is connected to the CAABA
    box model. The box-model related files are in the colored layers
    marked with BML, BMIL, and SMIL. The submodel core layer (SMCL) of
    MECCA is independent of the box model (see \citet{1664} for details
    about the MESSy layers). The arrows start at the module which is
    exporting the variables and subroutines. They point to the module
    importing them via the Fortran95 USE instruction. Here, the box {\tt
      messy\_mecca\_kpp*} represents all KPP-generated files. The
    KPP-internal structure is shown in Fig.~\ref{fig:kpp_use_diagr}.}
  \label{fig:caaba_use_diagr}
\end{figure*}

\subsection{Fortran95 files}
\label{sec:f95files}

The CAABA/MECCA simulations can be modified by changing the Fortran95
files (see Tab.~\ref{tab:files} for a list of files). The modular
structure of the Fortran95 files is shown in
Fig.~\ref{fig:caaba_use_diagr}. Most of the files need only be changed
by model developers. Those that are also interesting for model users,
are briefly explained below.

\subsubsection{{\tt caaba.f90}}

This file contains the main Fortran95 program (``\verb|PROGRAM caaba|'').

\subsubsection{{\tt caaba\_mem.f90}}

This file contains variable declarations which are needed by several
CAABA files.

\subsubsection{{\tt messy\_main\_control\_cb.f90}}

Flow control. Editing this file is only necessary when a new submodel is
added.

\subsubsection{{\tt messy\_jval\_box.f90}}

This file contains the connection of JVAL to CAABA.

\subsubsection{{\tt messy\_jval.f90} and {\tt messy\_jval\_jvpp.inc}}

These are the core files of JVAL, which calculate the J-values.

\subsubsection{{\tt messy\_mecca\_box.f90}}

This file contains the connection of MECCA to CAABA. The chemical
composition of seawater is defined in \verb|SUBROUTINE mecca_init|.

Aerosol properties (radius, liquid water content (LWC), and their
chemical composition) are defined in \verb|SUBROUTINE define_aerosol|.
In addition, the factor \verb|cvfac| is defined here, which converts the
aqueous-phase unit \unit{mol/L} (refering to the volume of the liquid)
to the gas-phase unit \unit{molecules/cm^3} (referring to the gas-phase
volume).

Initial mixing ratios of chemical species are defined in
\verb|SUBROUTINE x0|. Depending on which scenario was chosen in the
CAABA namelist file (see Sect.~\ref{sec:caaba_nml}), one of the
initialization subroutines \verb|x0_*| will be used.

\subsubsection{{\tt messy\_sappho\_box.f90}}

This file contains the connection of SAPPHO to CAABA.

\subsubsection{{\tt messy\_sappho.f90}}

This is the core file of SAPPHO, which defines simplified parameterized
photolysis rate coefficients.

\subsubsection{{\tt messy\_semidep\_box.f90}}

This file contains the connection of SEMIDEP to CAABA. Simplified
emission fluxes and deposition velocities are defined here.

\subsubsection{{\tt messy\_cmn\_gasaq.f90}}

This file contains gas- and aqueous-phase data of the chemicals: molar
masses (\verb|M|), Henry's law coefficients (\verb|Henry_T0|,
\verb|Henry_Tdep|), and accommodation coefficients (\verb|alpha_T0|,
\verb|alpha_Tdep|).

\subsubsection{{\tt messy\_mecca\_aero.f90}}

Several variables needed to calculate multiphase rate coefficients are
defined in \verb|messy_mecca_aero.f90|. The accommodation coefficients
from \verb|messy_cmn_gasaq.f90| are used for the calculation of the mass
transfer coefficients (\verb|kmt|). Together with the Henry's law
coefficients (also from \verb|messy_cmn_gasaq.f90|), they are needed to
calculate equilibria between the gas and the aqueous phase.
Heterogeneous reactions are described with the forward (\verb|k_exf|)
and backward (\verb|k_exb|) rate coefficients. The variable \verb|xaer|
is set to 1 or 0 to switch aqueous-phase chemistry on or off,
respectively.

\subsection{How to expand the chemical mechanism}
\label{sec:expandchem}

Small changes to the chemical mechanism can be made using ``replacement
files''. More information about this feature can be found in the file
\verb|mecca/rpl/gas.rpl-example|. This section contains brief
descriptions for experienced model developers explaining where to make
changes to the code for certain model expansions. In the descriptions,
``\verb|xyz|'' is used as an example for the name of the addition.

\subsubsection{Adding a new gas-phase species}

\begin{itemize}\nosep
\item \verb|gas.spc|:\\
  Add the new species, its elemental composition,
  the name in La\TeX\ syntax, and a comment, e.g.:\\
  \verb|NC4H10 = 4C + 10H ; {@C_4H_<10>} {n-butane}|\\
  Note that curly brackets needed by La\TeX\ must be entered as angle
  brackets.
\item \verb|jnl/xxxg.jnl|:\\
  Add one line per new species. Check if the new species is part of an
  existing familiy, e.g.\ add new reactive bromine species to
  \verb|Brx|.
\item \verb|jnl/tools/_kppvarg.jnl|:\\
  Add one line per new species.
\end{itemize}

\subsubsection{Adding a new aqueous-phase species}

\begin{itemize}\nosep
\item \verb|aqueous.spc|:\\
  Add the new species, the name in La\TeX\ syntax, and a comment,
  e.g.:\\
  \verb|SO4mm_a## = IGNORE; {@SO_4^<2->\aq} {sulfate}|\\
  The suffix \verb|_a##| is mandatory. The elemental composition is
  currently ignored. Note that curly brackets needed by La\TeX\ must be
  entered as angle brackets.
\item \verb|jnl/xxxa.jnl|:\\
  Add one line per new species.
\item \verb|jnl/_families_a.jnl|:\\
  Check if the new species is part of an existing familiy, e.g.\ add new
  bromine species to \verb|Brtot|.
\item \verb|jnl/tools/_kppvara.jnl|:\\
  Add one line per new species.
\end{itemize}

\subsubsection{Adding a new gas-phase reaction}
\label{sec:addgprxn}

First, choose an appropriate reaction number. To avoid that several
developers assign the same number to different new reactions, it is
strongly recommended that a preliminary reaction number is used
initially. This can be done by adding the developer's initials as a
suffix, e.g.\ John Doe would use \verb|G0001JD|, \verb|G0002JD|,
\verb|G0003JD|, and so on. When the new code is merged with other
development branches, the final reaction numbers will be assigned.

\begin{itemize}\nosep
\item \verb|gas.eqn|:
  \begin{itemize}
  \item Add one line per new reaction.
  \item Add Monte-Carlo uncertainty factor. If unknown, only add
    ``\verb|{§}|''.
  \end{itemize}
\item \verb|latex/meccanism.tex|:\\
  If necessary, add a footnote about the new reaction here.
\end{itemize}

\subsubsection{Adding a new gas-phase photolysis reaction}

First, add the reaction, as explained in Sect.~\ref{sec:addgprxn}.

In addition, check if the necessary photolysis rate coefficient is
already provided by SAPPHO, READJ, and/or JVAL. If not, add it:

\begin{itemize}\nosep
\item \verb|messy_cmn_photol_mem.f90|:
  \begin{itemize}\nosep
  \item Add a new index of photolysis \verb|ip_XYZ| at the end of the
    list.
  \item Increase \verb|IP_MAX|.
  \item Add the name to \verb|jname|.
  \end{itemize}
\item \verb|messy_sappho_box.f90|:\\
  Add \verb|XYZ| to \verb|CALL open_output_file| and
  \verb|CALL write_output_file|.
\item \verb|messy_sappho.f90|:\\
  Add a simple definition for \verb|jx(ip_XYZ)|.
\item \verb|messy_jval_box.f90|:\\
  Add one line.
\item \verb|messy_jval_jvpp.inc|:\\
  Calculate the definition with \verb|jvpp| or add it manually here.
\end{itemize}

\subsubsection{Adding a new aqueous-phase reaction}

First, choose an appropriate reaction number, as explained in
Sect.~\ref{sec:addgprxn}.

\begin{itemize}\nosep
\item \verb|aqueous.eqn|:
  \begin{itemize}
  \item Add one line per new reaction.
  \item Add Monte-Carlo uncertainty factor. If unknown, only add
    ``\verb|{§}|''.
  \end{itemize}
\item \verb|latex/meccanism.tex|:\\
  If necessary, add a footnote about the new reaction here.
\end{itemize}

\subsubsection{Adding a new Henry's law equilibrium}

First, choose an appropriate reaction number, as explained in
Sect.~\ref{sec:addgprxn}.

\begin{itemize}\nosep
\item \verb|aqueous.eqn|:
  \begin{itemize}
  \item Add two lines per new equilibrium, one for the forward and one
    for the backward reaction.
  \item Add Monte-Carlo uncertainty factors. If unknown, only add
    ``\verb|{§}|''.
  \end{itemize}
\item \verb|messy_cmn_gasaq.f90|:
  \begin{itemize}
  \item Add molar mass:\\
    \verb|CALL add_species('XYZ', ...)|
  \item Add the Henry's law coefficient:\\
    \verb|CALL add_henry('XYZ', ...)|
  \item Add the accommodation coefficient:\\
    \verb|CALL add_alpha('XYZ', ...)|
  \end{itemize}
  Using these data, the subroutines
  \verb|mecca_aero_trans_coeff|, \verb|mecca_aero_henry|,
  and \verb|mecca_aero_calc_k_ex| in
  \verb|messy_mecca_aero.f90| will automatically calculate
  the rate coefficients \verb|k_exf| and \verb|k_exb| for
  \verb|aqueous.eqn|.
\item \verb|latex/meccanism.tex|:\\
  If necessary, add a footnote about the new equilibrium here.
\end{itemize}

\subsubsection{Adding a new acid-base equilibrium}

First, choose an appropriate reaction number, as explained in
Sect.~\ref{sec:addgprxn}.

\begin{itemize}\nosep
\item \verb|aqueous.eqn|:
  \begin{itemize}
  \item Add two lines per new equilibrium, one for the forward and one
    for the backward reaction.
  \item Add Monte-Carlo uncertainty factors. If unknown, only add
    ``\verb|{§}|''.
  \end{itemize}
\item \verb|latex/meccanism.tex|:\\
  If necessary, add a footnote about the new equilibrium here.
\end{itemize}

\subsubsection{Adding a new label}
\label{sec:newlabel}

First, choose a name for the new label. The name must start with an
upper case letter and can be followed by one or more lower case letters
or numbers. Element symbols must not be used because they are reserved
for reactions of that element. For example, since S is sulfur, the
symbol S could not be used for the stratosphere. To avoid that several
developers introduce new labels with the same name for different
purposes, it is strongly recommended that a preliminary label is used
initially. This can be done by adding the developer's initials as a
prefix, e.g.\ John Doe would use \verb|Jd1|, \verb|Jd2|, \verb|Jd3|, and
so on. When the new code is merged with other development branches, a
final label name can be assigned.

\begin{itemize}\nosep
\item \verb|xmecca|:\\
  In the generation of \verb|awkfile1|, add another \verb|locate|
  function, and print the new label to the logfile.
\end{itemize}

\subsubsection{Adding a new emission}

\begin{itemize}\nosep
\item \verb|messy_semidep_box.f90|:\\
  Add one line to \verb|emission_default| (or one of the other
  \verb|emission_*| subroutines).
\end{itemize}

\subsubsection{Adding a new deposition}

\begin{itemize}\nosep
\item \verb|messy_semidep_box.f90|:\\
  Add one line to \verb|drydep_default| (or one of the other
  \verb|drydep_*| subroutines).
\end{itemize}

\subsubsection{Enlarging KPP}

If the selected chemistry mechanism has too many reactions, it may
become necessary to increase some limits of KPP. The main changes are
listed here:
\begin{itemize}\nosep
\item Increase \verb|MAX_EQN| and \verb|MAX_SPECIES| in \verb|gdata.h|.
\item With a large mechanism, some lines of the generated Fortran95 code
  become very long. The variable \verb|MAX_NO_OF_LINES| in
  \verb|code_f90.c| defines the maximum number of continuation lines
  that KPP will create. Unfortunately, if \verb|MAX_NO_OF_LINES| is too
  small, KPP may incorrectly split commands into two parts. Therefore, a
  large value of \verb|MAX_NO_OF_LINES| is needed for large mechanisms.
  Currently, \verb|MAX_NO_OF_LINES| is set to 100. For very large
  mechanisms, it may be necessary to increase it even further. Note,
  however, that the Fortran95 standard only allows a maximum of 39
  continuation lines and that you need a compiler that allows to exceed
  this limit.
\item To increase the maximum size of inlined code (\verb|F90_GLOBAL|
  etc.), change \verb|MAX_INLINE| in \verb|scan.h|.
\end{itemize}

To allow longer Fortran95 expressions for the rate coefficients in the
\verb|*.eqn| file, change consistently:
\begin{itemize}\nosep
\item \verb|crtToken|, \verb|nextToken|, \verb|crtFile|, and
  \verb|crt_rate| in \verb|scan.l|
\item \verb|MAX_K| in \verb|gdata.h| (note that \verb|MAX_EQNLEN| here
  only determines how long the printout of the equation in the monitor
  file will be)
\item \verb|union| in \verb|scan.y|
\end{itemize}

After each change in the source code, run \verb|gmake| again.

\subsection{How to modify and add new scenarios}
\label{sec:scenarios}

The photolysis (\verb|photo_scenario|), initialization
(\verb|init_scenario|), emission (\verb|emission_scenario|), and dry
deposition (\verb|drydep_scenario|) can be modified for the different
scenarios by entering appropriate values into the following subroutines:
\begin{center}
  \begin{tabular}{lll}
    \hline
    scenario & subroutine & file\\
    \hline
    \verb|photo|    & \verb|jvalues|   & \verb|messy_sappho.f90|\\
                    & \verb|jval_init| & \verb|messy_jval_box.f90|\\
    \verb|init|     & \verb|x0|        & \verb|messy_mecca_box.f90|\\
    \verb|emission| & \verb|emission|  & \verb|messy_semidep_box.f90|\\
    \verb|drydep|   & \verb|drydep|    & \verb|messy_semidep_box.f90|\\
    \hline
  \end{tabular}
\end{center}
To define a new scenario, it is also necessary to add its name to the
\verb|list_of_scenarios| in \verb|caaba.f90| and increase the value of
\verb|MAX_SCENARIOS| accordingly.

\subsection{How to modify the aerosol modes}

Aerosol properties are defined in \verb|SUBROUTINE define_aerosol| in
the file \verb|messy_mecca_box.f90|. Depending on the total number of
aerosol phases (\verb|APN|), different setups are defined. For example,
if \verb|APN = 2|, then an accumulation mode with sulfate (index 1) and
a coarse mode with sea salt are defined (index 2). For each aerosol
phase, the following properties must be defined inside the appropriate
\verb|CASE| section:

\begin{tabular}{lll}
  \hline
  variable        & definition                       & unit\\
  \hline
  \verb|xaer|     & off/on (\verb|0.| or \verb|1.|)  & \\
  \verb|radius|   & radius                           & \unit{m}\\
  \verb|lwc|      & liquid water content             & \unit{m^3/m^3}\\
  \verb|csalt|    & total salt concentration         & \unit{mcl/cm^3(air)}\\
  \verb|c0_NH4p|  & $c$(\chem{NH_4^+})               & \unit{mcl/cm^3(air)}\\
  \verb|c0_HCO3m| & $c$(\chem{HCO_3^-})              & \unit{mcl/cm^3(air)}\\
  \verb|c0_NO3m|  & $c$(\chem{NO_3^-})               & \unit{mcl/cm^3(air)}\\
  \verb|c0_Clm|   & $c$(\chem{Cl^-})                 & \unit{mcl/cm^3(air)}\\
  \verb|c0_Brm|   & $c$(\chem{Br^-})                 & \unit{mcl/cm^3(air)}\\
  \verb|c0_Im|    & $c$(\chem{I^-})                  & \unit{mcl/cm^3(air)}\\
  \verb|c0_IO3m|  & $c$(\chem{IO_3^-})               & \unit{mcl/cm^3(air)}\\
  \verb|c0_SO4mm| & $c$(\chem{SO_4^{2-}})            & \unit{mcl/cm^3(air)}\\
  \verb|c0_HSO4m| & $c$(\chem{HSO_4^-})              & \unit{mcl/cm^3(air)}\\
  \verb|exchng|   & $1/\tau$ = exchange              & \unit{1/s}\\
  \hline
\end{tabular}

These values can be changed to modify the properties (including the
initial chemical composition) of fresh aerosols. To add a new aerosol
mode, another \verb|CASE| section can be added. If a definition for the
chosen value of \verb|APN| exists already, the original code needs to be
deactivated. Otherwise, a block defining a new total number of aerosol
phases can simply be added.

Note that \verb|APN = 3| is a special case where the third mode reprents
sea water at the ocean surface (not aerosol particles). This is under
construction and not further described here.

\subsection{Input files}
\label{sec:inputfiles}

Data in netCDF files can be used to initialize chemical tracers and
photolysis rates as described below. In addition, time-dependent
photolysis rates of \chem{NO_2} from a netCDF input file can be used to
scale all other photolysis rates from the JVAL submodel. This has
already been described in Sect.~\ref{sec:jvalfile}.

\subsubsection{Tracer initialization file}
\label{sec:tracfile}

MECCA provides several initialization scenarios for a variety of species
and simulation aims (see \verb|init_scenario| in
Sect.~\ref{sec:scenarios}). However, for several consecutive simulations
with changing initialization, the possibility to define an
initialization file using \verb|init_spec| is convenient. The tracer
initialization file is a netCDF file with one point in time, at which
selected species' mixing ratios are defined in \unit{mol/mol}. The point
in time itself is not important and not checked when reading the file.
All tracers that are not specified in the input file are initialized
according to default or a chosen \verb|init_scenario|. An example for
the initialization file can be found at \verb|traject/example_init.nc|.

\subsubsection{Photolysis initialization file}
\label{sec:readjfile}

The submodel READJ reads J-values from a netCDF file. These values are
used throughout the whole model run. Thus, READJ is normally only used
when the model should run into steady-state (see
Sect.~\ref{sec:multirun}). It is possible to use one comprehensive
netCDF file containing values for both the tracer initialization
(Sect.~\ref{sec:tracfile}) and the photolysis initialization.

\subsection{How to add a new MESSy submodel}

Below you can find a brief description how to add a new MESSy submodel.
For more detailed information about programming in the MESSy system,
please refer to \citet{2400}.

\begin{itemize}
\item Choose a name (up to 7 lowercase alphanumerical characters,
  starting with a letter). Here, ``\verb|xyz|'' is used as an example.
\item \verb|caaba_mem.f90|:\\
  \verb|LOGICAL :: USE_XYZ = .FALSE.|
\item \verb|messy_xyz.f90|:\\
  Put all generic subroutines here, i.e.\ all subroutines that are used
  for the CAABA boxmodel as well as for a global model.
\item \verb|messy_xyz_box.f90|:\\
  Put CAABA-specific code here. Generic code is not included directly
  here. Instead, the generic subroutines in \verb|messy_xyz.f90| are
  called from here. This file contains up to four subroutines:
  \begin{itemize}
  \item If the submodel needs an initialization, put subroutine
    \verb|xyz_init| here.
  \item If the submodel performs calculations during the time loop, put
    subroutine \verb|xyz_physc| here.
  \item If the submodel prints results to the screen or netCDF files,
    put subroutine \verb|xyz_result| here.
  \item If the submodel needs to close any open files at the end of the
    model simulation, put subroutine \verb|xyz_finish| here.
  \end{itemize}
\item \verb|messy_main_control_cb.f90|:
  \begin{itemize}
  \item Add ``\verb|USE_XYZ|'' to ``\verb|USE caaba_mem|''
  \item If subroutine \verb|xyz_init| exists, add:\\
    \verb|USE messy_xyz_box, ONLY: xyz_init|\\
    \verb|IF (USE_XYZ) CALL xyz_init|\\
    to subroutine \verb|messy_init|.
  \item If subroutine \verb|xyz_physc| exists, add:\\
    \verb|USE messy_xyz_box, ONLY: xyz_physc|\\
    \verb|IF (USE_XYZ) CALL xyz_physc|\\
    to subroutine \verb|messy_physc|.
  \item If subroutine \verb|xyz_result| exists, add:\\
    \verb|USE messy_xyz_box, ONLY: xyz_result|\\
    \verb|IF (USE_XYZ) CALL xyz_result|\\
    to subroutine \verb|messy_result|.
  \item If subroutine \verb|xyz_finish| exists, add:\\
    \verb|USE messy_xyz_box, ONLY: xyz_finish|\\
    \verb|IF (USE_XYZ) CALL xyz_finish|\\
    to subroutine \verb|messy_finish|.
  \end{itemize}
\item \verb|caaba.f90|:\\
  Edit subroutine \verb|caaba_read_nml|:
  \begin{itemize}
  \item Add ``\verb|USE_XYZ|'' to ``\verb|USE caaba_mem|''.
  \item Add ``\verb|USE_XYZ|'' to namelist \verb|/CAABA/|.
  \item Print value of \verb|USE_XYZ| (see ``selected MESSy submodels'')
  \item If applicable, perform consistency checks for interaction of new
    submodel with other submodels.
  \end{itemize}
\item \verb|xyz.nml|:\\
  Add sensible default values for controlling the new submodel
  (optional).
\item \verb|manual/caaba_mecca_manual.tex|:\\
  Mention new submodel in this user manual (Sect.~\ref{sec:f95files},
  Tab.~\ref{tab:files}, and Fig.~\ref{fig:caaba_flowcontrol}).
\end{itemize}

\subsection{General programming guidelines}
\label{sec:guidelines}

To achieve some consistency of the model code, the following guidelines
should be adhered to when writing new code for the CAABA/MECCA modeling
system:

\begin{itemize}
\item All Fortran code must be compatible with the Fortran95 standard.
  Compiler-specific extensions should not be used.
\item Shell scripts should be written in \verb|tcsh| syntax unless there
  is a specific reason for using \verb|bash|, \verb|ksh|, or any other
  shell.
\item The file names of shell scripts should start with the letter
  \verb|x| for all scripts that are executed directly by the model user.
\item The file names of shell scripts should start with the underscore
  character (``\verb|_|'') if the script is called via another script
  but not executed directly be the user.
\item The names of temporary files should start with the prefix
  \verb|tmp_|.
\end{itemize}

\bibliographystyle{egu} % bst file
\bibliography{literat}  % bib files

\end{document}
