Index: /trunk/anuga_user_manual/source/anuga_user_manual.tex
===================================================================
--- /trunk/anuga_user_manual/source/anuga_user_manual.tex	(revision 9726)
+++ /trunk/anuga_user_manual/source/anuga_user_manual.tex	(revision 9727)
@@ -129,5 +129,6 @@
 cells, the topography and bathymetry, frictional resistance, initial
 values for water level (called \emph{stage}\index{stage} within \anuga),
-boundary conditions and forces such as rainfall, stream flows, windstress or pressure gradients if applicable.
+boundary conditions and operators such as rainfall, 
+stream flows, windstress or pressure gradients if applicable.
 
 \anuga tracks the evolution of water depth and horizontal momentum
@@ -160,12 +161,12 @@
 inundation software system, describe what it can do and give step-by-step
 instructions for setting up and running hydrodynamic simulations.
-The stable release of \anuga and this manual are available on sourceforge at
-\url{http://sourceforge.net/projects/anuga}. A snapshot of work in progress is
-available through the \anuga software repository at
-\url{https://anuga.anu.edu.au/svn/anuga/trunk/anuga_core/source/anuga}
-where the more adventurous reader might like to go.
+The stable releases of \anuga and this manual are available on github at
+\url{http://github.com/GeoscienceAustralia/anuga_core/releases}. 
+
+A snapshot of the current work in progress is at  
+\url{http://github.com/GeoscienceAustralia/anuga_core}. 
 
 This manual describes \anuga version \version. To check for later versions of this manual
-go to \url{https://anuga.anu.edu.au}.
+go to \url{http://github.com/GeoscienceAustralia/anuga_core/doc/anuga_user_manual.pdf}.
 
 \section{Scope}
@@ -178,5 +179,5 @@
 
 The latest installation instructions may be found at:
-\url{http://anuga.anu.edu.au/raw-attachment/wiki/WikiStart/anuga_installation_guide-1.2.1.pdf}.
+\url{http://github.com/GeoscienceAustralia/anuga_core/INSTALL.rst}.
 
 \section{Audience}
@@ -184,6 +185,6 @@
 Readers are assumed to be familiar with the Python Programming language and
 its object oriented approach.
-Python tutorials include
-\url{http://docs.python.org/tut} and \url{http://www.sthurlow.com/python}.
+A good Python tutorial is
+\url{https://docs.python.org/2/tutorial}.
 
 Readers also need to have a general understanding of scientific modelling,
@@ -195,5 +196,5 @@
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
-\chapter{Background}
+\section{Background}
 
 Modelling the effects on the built environment of natural hazards such
@@ -217,13 +218,13 @@
 method to accommodate discontinuities in the solution\footnote{
 While \anuga works with discontinuities in the conserved quantities stage,
-xmomentum and ymomentum, it does not allow discontinuities in the bed elevation.}.
+xmomentum and ymomentum, it does not allow discontinuities in the bed elevation.
+The next version 1.4 works with discontinuous bed elevation, and the current version 
+can use discontinuous bed if the flow algorithm is set to \texttt{DE0} or \texttt{DE1}}.
 
 To set up a particular scenario the user specifies the geometry
 (bathymetry and topography), the initial water level (stage),
-boundary conditions such as tide, and any forcing terms that may
-drive the system such as rainfall, abstraction of water, wind stress or atmospheric pressure
-gradients. Gravity and frictional resistance from the different
-terrains in the model are represented by predefined forcing terms.
-See section \ref{sec:forcing terms} for details on forcing terms available in \anuga.
+boundary conditions such as tide, and any operators  that may
+drive the system such as rainfall, abstraction of water,  erosion, culverts
+See section \ref{sec:operators} for details of operators available in \anuga.
 
 The built-in mesh generator, called \code{graphical\_mesh\_generator},
@@ -238,11 +239,19 @@
 and can be readily adapted to changing requirements throughout its
 lifetime.  Computationally intensive components are written for
-efficiency in C routines working directly with Python numeric
-structures.  The animation tool developed for \anuga is based on
+efficiency in C routines working directly with Python numpy
+structures.  
+
+The visualisation tool developed for \anuga is based on
 OpenSceneGraph, an Open Source Software (OSS) component allowing high
 level interaction with sophisticated graphics primitives.
 See \cite{nielsen2005} for more background on \anuga.
 
-\chapter{Restrictions and limitations on \anuga}
+
+ 
+
+ 
+
+
+\section{Restrictions and limitations on \anuga}
 \label{ch:limitations}
 
@@ -270,4 +279,34 @@
 \end{itemize}
 
+\pagebreak
+\section{Citing \anuga}
+
+When citing \anuga cite this manual:
+
+Bibtex entry:
+\begin{verbatim}
+@manual{anugamanual,
+  title={ANUGA User Manual},
+  author={Roberts, S and Nielsen, O. and Gray, D. and Sexton, J.},
+  organization={Geoscience Australia},
+  year={2015}
+}
+\end{verbatim}
+
+or this article
+
+Bibtex entry:
+\begin{verbatim}
+@article{nielsen2005hydrodynamic,
+  title={Hydrodynamic modelling of coastal inundation},
+  author={Nielsen, O and Roberts, S and Gray, D and McPherson, A and Hitchman, A},
+  journal={MODSIM 2005 International Congress on Modelling and Simulation},
+  pages={518--523},
+  year={2005}
+}
+\end{verbatim}
+
+
+
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
@@ -283,6 +322,4 @@
 \section{A Simple Example}
 \label{sec:simpleexample}
-
-\subsection{Overview}
 
 What follows is a discussion of the structure and operation of a
@@ -320,5 +357,5 @@
 \end{verbatim}
 
-\subsection{Outline of the Program}
+\section{Outline of the Program}
 
 In outline, \file{runup.py} performs the following steps:
@@ -336,5 +373,5 @@
 \end{enumerate}
 
-\subsection{The Code}
+\section{The Code}
 
 For reference we include below the complete code listing for
@@ -347,5 +384,5 @@
 
 
-\subsection{Establishing the Domain}\index{domain, establishing}
+\section{Establishing the Domain}\index{domain, establishing}
 
 The very first thing to do is import the various modules, of which the 
@@ -376,5 +413,7 @@
    \item a list \code{vertices} specifying the three vertices of each triangle, and
    \item a dictionary \code{boundary} that stores the edges on
-         the boundary and associates with each a symbolic tag. The edges are represented as pairs (i, j) where i refers to the triangle id and j to the edge id of that triangle. 
+         the boundary and associates with each a symbolic tag. 
+         The edges are represented as pairs (i, j) where i refers 
+         to the triangle id and j to the edge id of that triangle. 
          Edge ids are enumerated from 0 to 2 based on the id of the vertex opposite. 
 \end{itemize}
@@ -411,5 +450,5 @@
 However, this is not necessary, as the above is the default behaviour.
 
-\subsection{Initial Conditions}
+\section{Initial Conditions}
 
 The next task is to specify a number of quantities that we wish to
@@ -427,5 +466,5 @@
 forms in which we can invoke \method{set\_quantity}.
 
-\subsubsection{Elevation}
+\subsection{Elevation}
 
 The elevation, or height of the bed, is set using a function
@@ -454,5 +493,5 @@
 compatible. For example, using square root will not work.
 
-\subsubsection{Friction}
+\subsection{Friction}
 
 The assignment of the friction quantity (a forcing term)
@@ -467,5 +506,5 @@
 at every mesh point.
 
-\subsubsection{Stage}
+\subsection{Stage}
 
 The stage (the height of the water surface) is related to the
@@ -509,5 +548,5 @@
 details.
 
-\subsection{Boundary Conditions}\index{boundary conditions}
+\section{Boundary Conditions}\index{boundary conditions}
 
 The boundary conditions are specified as follows:
@@ -629,5 +668,5 @@
 \end{verbatim}
 
-\subsection{Evolution}\index{evolution}
+\section{Evolution}\index{evolution}
 
 The final statement: 
@@ -644,5 +683,5 @@
 outputs.  Behind the scenes more time steps are generally taken.
 
-\subsection{Output}
+\section{Output}
 
 The output is a NetCDF file with the extension \code{.sww}. It
@@ -701,8 +740,8 @@
 
 
-\section{A slightly more complex example}
+%====================================================
+\chapter{A slightly more complex example}
 \label{sec:channelexample}
 
-\subsection{Overview}
 
 The next example is about water-flow in a channel with varying boundary conditions and
@@ -711,5 +750,5 @@
 The example will be built up through three progressively more complex scripts.
 
-\subsection{Overview}
+\section{Overview}
 
 As in the case of \file{runup.py}, the actions carried
@@ -727,5 +766,5 @@
 \end{enumerate}
 
-\subsection{The Code}
+\section{The Code}
 
 Here is the code for the first version of the channel flow \file{channel1.py}:
@@ -736,5 +775,5 @@
 given above, discussing each major step of the code in turn.
 
-\subsection{Establishing the Mesh}\index{mesh, establishing}
+\section{Establishing the Mesh}\index{mesh, establishing}
 
 In this example we use a similar simple structured triangular mesh 
@@ -796,5 +835,5 @@
 Here is the code for the second version of the channel flow \file{channel2.py}:
 
-\verbatiminputunderscore{demos/channel2.py}
+\verbatiminputunderscore{../../anuga_core/examples/channel2.py}
 
 This example differs from the first version in that a constant outflow boundary condition has
@@ -830,5 +869,5 @@
 used to model the breakdown of a sluice door when water exceeds a certain level.
 
-\subsection{Output}
+\section{Output}
 
 The text output from this example looks like this:
@@ -845,5 +884,5 @@
 \end{verbatim}
 
-\subsection{Flow through more complex topographies}
+\section{Flow through more complex topographies}
 
 Here is the code for the third version of the channel flow \file{channel3.py}:
@@ -1351,141 +1390,4 @@
 \end{figure}
 
-%=================================
-\chapter{Parallel Simulation}
-
-The examples from the previous chapters were run just using one processor. As you will have noticed, the simulations can take a long time to run, especially if you are using a fine mesh. Indeed as you halve the spacing of the mesh, the number of triangles goes up by a factor of 4 and the timestep halves, so generally halving the spacing of the mesh increases the run time by approximately a factor of 8. 
-
-One way to alleviate this is to run your simulation in parallel and use more processors to spread the computational cost. 
-
-
-\anuga has the option of running in parallel using \mpi. A description of the method used to parallelize \anuga is given in section~\ref{theory:parallel}.
-
-Such jobs are run using the command
-
-\begin{verbatim}
-mpirun -np n python runscript.py
-\end{verbatim}
-where \code{n} is the total number of processors being used for the parallel run. 
-
-Essentially we can expect speedups comparable to the number of cores available. This is measured via scalability. Our current experiments have demonstrated a scalability of 70\% or above when the number of processors are chosen so that the local partitioned meshes contain around 2000 triangles. 
-
-For instance, suppose we have a problem with a mesh of 500,000 triangles, and use 250 processors. Then the mesh will be partitioned into sub meshes of approximately 2000 triangles. 
-We would expect 70\% scalability, and so expect a speedup of $250 \times 0.7 = 175$. 
-
-\section{The Code}
-
-Here is the code for \file{runParallelCairns.py} which is found in the \file{demos/cairns} directory:
-
-\verbatiminputunderscore{../../anuga_core/demos/cairns/runParallelCairns.py}
-
-
-\section{Structure of the Code}
-
-The code is very similar to the sequential code. The same procedures are used to setup the domain, setup the inital conditions, boundary conditions and evolve. 
-
-
-We first import a few procedures  need for the parallel code.
-\begin{verbatim}
-from anuga import distribute, myid, numprocs, finalize, barrier
-\end{verbatim}
-
-\code{myid} returns the id of the current processor running the code. 
-
-\code{numprocs} returns the total number of processors involved in this parallel jobs (the \code{n} in the original \code{mpirun} command. 
-
-\code{finalize} is called at the end of a script to close down the parallel job.
-
-\code{distribute} is used to partition and setup the parallel domains.
-
-\code{barrier} is a command for processors to wait as all the other processor to catchup to this point.
-
-
-
-
-The creation of the \code{domain} is only done on processor 0. Hence we have the structure:
-
-\begin{verbatim}
-#-------------------------------------
-# Do the domain creation on processor 0
-#-------------------------------------
-if myid == 0: 
-    ....
-    domain = ...
-	
-else:
-    domain = None
-\end{verbatim}
-
-We only need to create the original domain on one processor, otherwise we will have multiple copies of the full domain (which will easily eat up our memory). 
-
-Once we have our code{domain} setup we partition it and send the partitions to each of the other processors, via the command
-
-\begin{verbatim}
-#-------------------------------------
-# Now produce parallel domain
-#-------------------------------------
-domain = distribute(domain)
-\end{verbatim}
-
-This takes the \code{domain} on processor 0 and distributes that domain to each of the processors, (it overwrites the full domain on processor 0).  From this point of the code, there is a different domain on each processor, with each domain comunicating with the other domains to ensure required transfer of inforation to allow flow over the combined domains. 
-
-It is important to apply the boundary conditions after the \code{distribute}
-
-
-
-As the partitoined domain evolve, they will store their data to individual sww files, named as \code{domain_name_Pn_m.sww}, where \code{n} is the total number of processors being used and \code{m} is the specific processor id. 
-
-We have a procedure to merge these individual sww files via the command 
-
-\begin{verbatim}
-domain.sww_merge()
-\end{verbatim}
-
-And we close down the parallel job by issuing the command 
-
-\begin{verbatim}
-finalize()
-\end{verbatim}
-
-
-%=================================
-\chapter{Checkpointing}
-
-When running large and long running simulations, it is useful to provide a mechanism to allow the simulation to restarted if it is interrupted mid simulation. Maybe there is a power cut, or when using batch queues there may be a time restriction on the length of any one running job. Then the use of checkpointing becomes useful. 
-
-The idea is that at regular intervals the system makes a copy of the current state of the computation. Then if there is an interruption the computation can be started from the last checkpoint time and not original start time. 
-
-\section{The Code}
-
-Here is the code for \file{runCheckpoint.py}:
-
-\verbatiminputunderscore{../../anuga_core/demos/checkpointing/runCheckpoint.py}
-
-\section{Structure of the Code}
-
-As usual the code needs to import the required modules, and setup parameters and in this case procedures for setting up the stage and elevation. 
-
-Then we use a \code{try: except:} statement, where we try to open any appropriate checkpoint files, and if not successful, create and distribute a domain as usual. In the creation of the domain, we setup checkpointing just before we start the evolve loop. 
-
-
-\begin{verbatim}
-try:	
-    from anuga import load_checkpoint_file
-	
-    domain = load_checkpoint_file(domain_name = domain_name, checkpoint_dir = checkpoint_dir)
-
-except:
-    # create the domain as usual
-\end{verbatim}
-
-The \code{load\_checkpoint\_file} needs to know where to look for the checkpoint files (by default the current directory) and the domain name (default ``domain'').
-
-When the domain is originally created, we need to setup checkpointing via the command:
-\begin{verbatim}
-    if useCheckpointing:
-        domain.set_checkpointing(checkpoint_time = 5)
-\end{verbatim}
-
-Here we are setting the wall time between saving of checkpoint files to 5 sec (this is just for testing). Normally we would set the checkpointing to something large, say 15 minutes for a run of multiple hours. 
 
 
@@ -1502,5 +1404,5 @@
 will run the tests in verbose mode (produces output) and in parallel using 6 processors (if parallel version of anuga has been setup). 
 
-A subset of these validation tests can be run using \code{run_auto_validation_tests.py} from the \module{anuga_core} directory. 
+A subset of these validation tests can be run using \code{run_validations.py} from the \module{anuga_core} directory. 
 
 
@@ -1585,111 +1487,6 @@
 
 
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-\chapter{\anuga Public Interface}
-\label{ch:interface}
-
-This chapter gives an overview of the features of \anuga available
-to the user at the public interface. These are grouped under the
-following headings, which correspond to the outline of the examples
-described in Chapter \ref{ch:getstarted}:
-\begin{itemize}
-    \item Establishing the Mesh: Section \ref{sec:establishing the mesh}
-    \item Initialising the Domain: Section \ref{sec:initialising the domain}
-%    \item Specifying the Quantities: Section \ref{sec:quantities}
-    \item Initial Conditions: Section \ref{sec:initial conditions}
-    \item Boundary Conditions: Section \ref{sec:boundary conditions}
-    \item Operators: Section \ref{sec:operators}
-    \item Evolution: Section \ref{sec:evolution}
-\end{itemize}
-
-\section{Documentation}\index{Documentation}
-
-The listings here are intended merely to give the reader an idea of what
-each feature is and how it can be used -- they do
-not give full specifications; for these the reader
-may consult the programmer's guide. The code for every function or class contains
-a documentation string, or 'docstring', that specifies the precise
-syntax for its use. This appears immediately after the line
-introducing the code, between two sets of triple quotes.
-
-Python has a handy tool that lets you easily navigate this documentation,
-called \code{pydoc}. In Linux, it runs as a server, which serves the
-documentation up to your web browser:
-
-1. Open a terminal at your \code{anuga} folder
-2. Start the python documentation server with \code{pydoc -p 6767}
-3. Open a browser and type in \code{http://localhost:6767/}
-
-Now you have a real-time programmers' guide for \anuga, and an easy way to find 
-the functions you are interested in. Pydoctor and other Python doc generators
-look nicer and have graphs, etc, but pydoc works straight out of the box. 
-
-\section{Public vs Private Interface}\index{public vs private interface}
-
-To simplify the process of writing scripts, \anuga has a public API which packages
-up the commonly-used functionality of \anuga in the one place. To use it, simply
-import the anuga module like so:
-
-\begin{verbatim}
-import anuga
-\end{verbatim}
-
-You can now use the public API like so. Note the \code{anuga.} prefix:
-
-\begin{verbatim}
-anuga.sww2dem('in.sww', 'out.asc')
-\end{verbatim}
-
-If you wish to delve "under the hood" and modify the way \anuga runs at a more
-advanced level, you need to specify the full location of the module like so:
-
-\begin{verbatim}
-from anuga.fit_interpolate.interpolate import Interpolation_interface
-\end{verbatim}
-
-All modules are in the folder \file{anuga} or one of its subfolders, and the
-location of each module is described relative to \file{anuga}. Rather
-than using pathnames, whose syntax depends on the operating system,
-we use the format adopted for importing the function or class for
-use in Python code. For example, suppose we wish to specify that the
-function \function{create\_mesh\_from\_regions} is in a module called
-\module{mesh\_interface} in a subfolder of \module{anuga} called
-\code{pmesh}. In Linux or Unix syntax, the pathname of the file
-containing the function, relative to \file{anuga}, would be:
-
-\begin{verbatim}
-pmesh/mesh_interface.py
-\end{verbatim}
-
-\label{sec:mesh interface}
-while in Windows syntax it would be:
-
-\begin{verbatim}
-pmesh\mesh_interface.py
-\end{verbatim}
-
-Rather than using either of these forms, in this chapter we specify
-the location simply as \code{anuga.pmesh.mesh_interface}, in keeping with
-the usage in the Python statement for importing the function,
-namely:
-
-\begin{verbatim}
-from anuga.pmesh.mesh_interface import create_mesh_from_regions
-\end{verbatim}
-
-
-The following parameters are common to many functions and classes
-and are omitted from the descriptions given below:
-
-%\begin{tabular}{ll}
-\begin{tabular}{p{2.0cm} p{14.0cm}}
-  \emph{use\_cache} & Specifies whether caching is to be used for improved performance. See Section \ref{sec:caching} for details\\
-                    & on the underlying caching functionality\\                      
-  \emph{verbose}    & If \code{True}, provides detailed terminal output to the user\\
-\end{tabular}
-
-%=======================================================================
-\section{Mesh Generation}\index{Mesh!generation}
+%=================================
+\chapter{Mesh Generation}\index{Mesh!generation}
 \label{sec:establishing the mesh}
 Before discussing the part of the interface relating to mesh
@@ -1707,6 +1504,4 @@
   \label{fig:simplemesh}
 \end{figure}
-
-\clearpage
 
 The variables \code{points}, \code{triangles} and \code{boundary}
@@ -1774,6 +1569,6 @@
 
 
-%===========================================
-\subsection{High Level Mesh Generation Functions}
+%------------------------------------------------------------------
+\section{High Level Mesh Generation Functions}
 % \refmodindex[pmesh.meshinterface]{pmesh.mesh\_interface}
 \label{sec:meshgeneration}
@@ -1912,5 +1707,6 @@
 \end{funcdesc}
 
-\subsection{Advanced mesh generation}
+%-----------------------------------------------------
+\section{Advanced mesh generation}
 
 For more control over the creation of the mesh outline, use the
@@ -1938,5 +1734,5 @@
 \end{classdesc}
 
-\subsubsection{Key Methods of Class Mesh}
+\subsection{Key Methods of Class Mesh}
 
 \begin{methoddesc}{\emph{<mesh>}.add_hole}{x, y, geo_reference=None}
@@ -2103,6 +1899,7 @@
 \end{methoddesc}
 
-%==============================================
-\section{Initialising the Domain}\index{Initialising the Domain}
+
+%=================================
+\chapter{The Domain}\index{Initialising the Domain}
 \label{sec:initialising the domain}
 
@@ -2126,6 +1923,7 @@
 Module: \refmodule{abstract_2d_finite_volumes.domain}
 
-This class is used to create an instance of a structure used to
-store and manipulate data associated with a mesh. The mesh is
+This class isused to
+store and manipulate data associated with a mesh and to setup the particulars for 
+the evolution of the specific problem. The mesh is
 specified either by assigning the name of a mesh file to
 \code{source} or by specifying the points, triangle and boundary of the
@@ -2133,5 +1931,6 @@
 \end{classdesc}
 
-\subsection{Key Methods of Domain}
+%----------------------------------------------
+\section{Key Methods of Domain}
 
 \begin{methoddesc}{\emph{<domain>}.set_name}{name}
@@ -2365,7 +2164,322 @@
 \end{methoddesc}
 
-
-\section{Initial Conditions}\index{Initial Conditions}
+%-----------------------------------------
+\section{Evolution}\index{evolution}
+\label{sec:evolution}
+
+\begin{methoddesc}{\emph{<domain>}.evolve}{yieldstep=None,
+                                           finaltime=None,
+                                           duration=None,
+                                           skip_initial_step=False}
+Module: \module{abstract_2d_finite_volumes.domain}
+
+This method is invoked once all the
+preliminaries have been completed, and causes the model to progress
+through successive steps in its evolution, storing results and
+outputting statistics whenever a user-specified period
+\code{yieldstep} is completed.  Generally during this period the
+model will evolve through several steps internally
+as the method forces the water speed to be calculated
+on successive new cells. 
+
+\code{yieldstep} is the interval in seconds between yields where results are
+stored, statistics written and the domain is inspected or possibly modified.
+If omitted an internal predefined \code{yieldstep} is used.  Internally, smaller
+timesteps may be taken.
+
+\code{duration} is the duration of the simulation in seconds.
+
+\code{finaltime} is the time in seconds where simulation should end. This is currently
+relative time, so it's the same as \code{duration}.  If both \code{duration} and
+\code{finaltime} are given an exception is thrown.
+
+\code{skip_initial_step} is a boolean flag that decides whether the first
+yield step is skipped or not. This is useful for example to avoid
+duplicate steps when multiple evolve processes are dove tailed.
+
+The code specified by the user in the block following the evolve statement is
+only executed once every \code{yieldstep} even though 
+\anuga typically will take many more internal steps behind the scenes.
+
+You can include \method{evolve} in a statement of the type:
+
+\begin{verbatim}
+for t in domain.evolve(yieldstep, finaltime):
+    <Do something with domain and t>
+\end{verbatim}
+\end{methoddesc}
+
+%-----------------------------------------
+\section{Diagnostics}
+\label{sec:diagnostics}
+
+\begin{methoddesc}{\emph{<domain>}.statistics}{}
+Module: \module{abstract\_2d\_finite\_volumes.domain}
+
+Outputs statistics about the mesh within the \code{Domain}.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<domain>}.timestepping_statistics}{track_speeds=False, triangle_id=None}
+Module: \module{abstract_2d_finite_volumes.domain}
+
+Returns a string of the following type for each timestep:\\
+\code{Time = 0.9000, delta t in [0.00598964, 0.01177388], steps=12 (12)}
+
+Here the numbers in \code{steps=12 (12)} indicate the number of steps taken and
+the number of first-order steps, respectively.
+
+The optional keyword argument \code{track_speeds} will
+generate a histogram of speeds generated by each triangle if set to \code{True}. The
+speeds relate to the size of the timesteps used by \anuga and
+this diagnostics may help pinpoint problem areas where excessive speeds
+are generated.
+
+The optional keyword argument \code{triangle_id} can be used to specify a particular
+triangle rather than the one with the largest speed.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<domain>}.boundary_statistics}{quantities=None,
+                                                      tags=None}
+Module: \module{abstract_2d_finite_volumes.domain}
+
+Generates output about boundary forcing at each timestep.
+
+\code{quantities} names the quantities to be reported -- may be \code{None},
+a string or a list of strings.
+
+\code{tags} names the tags to be reported -- may be either None, a string or a list of strings.
+
+When \code{quantities = 'stage'} and \code{tags = ['top', 'bottom']}
+will return a string like:
+
+\begin{verbatim}
+Boundary values at time 0.5000:
+    top:
+        stage in [ -0.25821218,  -0.02499998]
+    bottom:
+        stage in [ -0.27098821,  -0.02499974]
+\end{verbatim}
+\end{methoddesc}
+
+\begin{methoddesc}{Q = \emph{<domain>}.get_quantity}{name,
+                                               location='vertices',
+                                               indices=None}
+Module: \module{abstract_2d_finite_volumes.domain}
+
+This function returns a Quantity object Q.
+Access to its values should be done through \code{Q.get_values()} documented on Page \pageref{pg:get values}.
+
+\code{name} is the name of the quantity to retrieve.
+
+\code{location} is
+
+\code{indices} is
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<domain>}.set_quantities_to_be_monitored}{quantity,
+                                             polygon=None,
+                                             time_interval=None}
+Module: \module{abstract\_2d\_finite\_volumes.domain}
+
+Selects quantities and derived quantities for which extrema attained at internal timesteps
+will be collected.
+
+\code{quantity} specifies the quantity or quantities to be monitored and must be either:
+\begin{itemize}
+  \item the name of a quantity or derived quantity such as 'stage-elevation',
+  \item a list of quantity names, or
+  \item \code{None}.
+\end{itemize}
+
+\code{polygon} can be used to monitor only triangles inside the polygon. Otherwise
+all triangles will be included.
+
+\code{time_interval} will restrict monitoring to time steps in the interval. Otherwise
+all timesteps will be included.
+
+Information can be tracked in the evolve loop by printing \code{quantity_statistics} and
+collected data will be stored in the SWW file.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<domain>}.quantity_statistics}{precision='\%.4f'}
+Module: \module{abstract_2d_finite_volumes.domain}
+
+Reports on extrema attained by selected quantities.
+
+Returns a string of the following type for each timestep:
+
+\begin{verbatim}
+Monitored quantities at time 1.0000:
+  stage-elevation:
+    values since time = 0.00 in [0.00000000, 0.30000000]
+    minimum attained at time = 0.00000000, location = (0.16666667, 0.33333333)
+    maximum attained at time = 0.00000000, location = (0.83333333, 0.16666667)
+  ymomentum:
+    values since time = 0.00 in [0.00000000, 0.06241221]
+    minimum attained at time = 0.00000000, location = (0.33333333, 0.16666667)
+    maximum attained at time = 0.22472667, location = (0.83333333, 0.66666667)
+  xmomentum:
+    values since time = 0.00 in [-0.06062178, 0.47886313]
+    minimum attained at time = 0.00000000, location = (0.16666667, 0.33333333)
+    maximum attained at time = 0.35103646, location = (0.83333333, 0.16666667)
+\end{verbatim}
+
+The quantities (and derived quantities) listed here must be selected at model
+initialisation time using the method \code{domain.set_quantities_to_be_monitored()}.
+
+The optional keyword argument \code{precision='\%.4f'} will
+determine the precision used for floating point values in the output.
+This diagnostics helps track extrema attained by the selected quantities
+at every internal timestep.
+
+These values are also stored in the SWW file for post-processing.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<quantity>}.get_values}{interpolation_points=None,
+                   location='vertices',
+                   indices=None,
+                   use_cache=False,
+                   verbose=False}
+\label{pg:get values}
+Module: \module{abstract_2d_finite_volumes.quantity}
+
+Extract values for quantity as a numeric array.
+
+\code{interpolation_points} is a list of (x, y) coordinates where the value is
+sought (using interpolation). If \code{interpolation_points} is given, values
+for \code{location} and \code{indices} are ignored.
+Assume either an absolute UTM coordinates or geospatial data object.
+    
+\code{location} specifies where values are to be stored.
+Permissible options are 'vertices', 'edges', 'centroids' or 'unique vertices'.
+
+The returned values will have the leading dimension equal to length of the \code{indices} list or
+N (all values) if \code{indices} is \code{None}.
+
+If \code{location} is 'centroids' the dimension of returned
+values will be a list or a numeric array of length N, N being
+the number of elements.
+      
+If \code{location} is 'vertices' or 'edges' the dimension of
+returned values will be of dimension \code{Nx3}.
+
+If \code{location} is 'unique vertices' the average value at
+each vertex will be returned and the dimension of returned values
+will be a 1d array of length "number of vertices" 
+      
+\code{indices} is the set of element ids that the operation applies to.
+
+The values will be stored in elements following their internal ordering.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<quantity>}.set_values}{numeric=None,
+                               quantity=None,
+                               function=None,
+                               geospatial_data=None,
+                               filename=None,
+                               attribute_name=None,
+                               alpha=None,
+                               location='vertices',
+                               polygon=None,
+                               indices=None,
+                               smooth=False,
+                               verbose=False,
+                               use_cache=False}
+Module: \module{abstract_2d_finite_volumes.quantity}
+
+Assign values to a quantity object.
+
+This method works the same way as \code{set_quantity()} except that it doesn't take
+a quantity name as the first argument since it is applied directly to the quantity.
+Use \code{set_values} is used to assign
+values to a new quantity that has been created but which is
+not part of the domain's predefined quantities.
+
+\code{location} is ??
+
+\code{indices} is ??
+
+The method \code{set_values()} is always called by \code{set_quantity()}
+behind the scenes.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<quantity>}.get_integral}{}
+Module: \module{abstract_2d_finite_volumes.quantity}
+
+Return the computed integral over the entire domain for the quantity.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<quantity>}.get_maximum_value}{indices=None}
+Module: \module{abstract_2d_finite_volumes.quantity}
+
+Return the maximum value of a quantity (on centroids).
+
+\code{indices} is the optional set of element \code{id}s that
+the operation applies to.
+
+We do not seek the maximum at vertices as each vertex can
+have multiple values -- one for each triangle sharing it.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<quantity>}.get_maximum_location}{indices=None}
+Module: \module{abstract_2d_finite_volumes.quantity}
+
+Return the location of the maximum value of a quantity (on centroids).
+
+\code{indices} is the optional set of element \code{id}s that
+the operation applies to.
+
+We do not seek the maximum at vertices as each vertex can
+have multiple values -- one for each triangle sharing it.
+
+If there are multiple cells with the same maximum value, the
+first cell encountered in the triangle array is returned.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<domain>}.get_wet_elements}{indices=None}
+Module: \module{shallow_water.shallow_water_domain}
+
+Returns the indices for elements where h $>$ minimum_allowed_height
+
+\code{indices} is the optional set of element \code{id}s that
+the operation applies to.
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<domain>}.get_maximum_inundation_elevation}{indices=None}
+Module: \module{shallow_water.shallow_water_domain}
+
+Return highest elevation where h $>$ 0.
+
+\code{indices} is the optional set of element \code{id}s that
+the operation applies to.
+
+Example to find maximum runup elevation:
+\begin{verbatim}
+z = domain.get_maximum_inundation_elevation()
+\end{verbatim}
+\end{methoddesc}
+
+\begin{methoddesc}{\emph{<domain>}.get_maximum_inundation_location}{indices=None}
+Module: \module{shallow_water.shallow_water_domain}
+
+Return location (x,y) of highest elevation where h $>$ 0.
+
+\code{indices} is the optional set of element \code{id}s that
+the operation applies to.
+
+Example to find maximum runup location:
+\begin{verbatim}
+x, y = domain.get_maximum_inundation_location()
+\end{verbatim}
+\end{methoddesc}
+
+
+
+
+
+%=================================
+\chapter{Initial Conditions}\index{Initial Conditions}
 \label{sec:initial conditions}
+
 In standard usage of partial differential equations, initial conditions
 refers to the values associated to the system variables (the conserved
@@ -2504,5 +2618,5 @@
 \end{methoddesc}
 
-\begin{funcdesc}{slump_tsunami}{length, depth, slope, width=None, thickness=None,
+\begin{funcdesc}{\emph{<callable_object>} =slump_tsunami}{length, depth, slope, width=None, thickness=None,
                                 radius=None, dphi=0.48, x0=0.0, y0=0.0, alpha=0.0,
                                 gravity=9.8, gamma=1.85,
@@ -2631,5 +2745,8 @@
 
 
-\section{Boundary Conditions}\index{boundary conditions}
+
+%=================================
+
+\chapter{Boundary Conditions}\index{boundary conditions}
 \label{sec:boundary conditions}
 
@@ -2659,5 +2776,6 @@
 \end{methoddesc}
 
-\subsection{Predefined boundary conditions}
+%----------------------------------------------------------
+\section{Predefined boundary conditions}
 
 \begin{classdesc}{Reflective_boundary}{domain=None}
@@ -2806,5 +2924,5 @@
 \end{classdesc}
 
-\subsection{User-defined boundary conditions}
+\section{User-defined boundary conditions}
 \label{sec:roll your own}
 
@@ -2821,4 +2939,542 @@
 Please refer to the source code for the existing boundary conditions
 for examples of how to implement boundary conditions.
+
+
+%=================================
+\chapter{Operators and Structures}\index{operators}\index{structures}
+\label{sec:operators}
+
+A way to effect the evolution is via fractional step \code{operators}.  The idea is that at each inner timestep we can use 
+an \code{operator} to change the centroid values of quantities such as the stage and the xmomentum and the ymomentum. T
+
+If you are using one of the older flow argorithms then the elevation can also be changed but more care needs to be applied as the \code{elevation} quantity needs to remain continuous. 
+The moral, play with \code{elevation} at your peril if using the old algorithms. 
+
+The newer algorithms (check that you are using discontinuous elevation with the domain member function \code{<domain>.get_using_discontinuous_elevation()}
+
+Currently predefined operators:
+
+
+\begin{classdesc}{Set_elevation_operator}{domain,
+                 elevation=None,
+                 indices=None,
+                 polygon=None,
+                 center=None,
+                 radius=None,
+                 description = None,
+                 label = None,
+                 logging = False,
+                 verbose = False}
+Module: \module{anuga.operators.set_elevation_operator}
+
+Set the elevation in a region (careful to maintain continuity of elevation)
+
+
+\code{domain} is the domain being evolved.
+
+\code{elevation} a function of $t$ or $x,y$ or $x,y,t$ written to take numpy arrays $x$ and $y$. $t$ is a scalar.
+
+\code{indices} is a list of triangle indices where the function \code{elevation} will be applied.
+
+\code{polygon} is a polygon. The function \code{elevation} will be applied to triangles with in the polygon.
+
+\code{center} is the center of a circle where the function \code{elevation} will be applied. (\code{radius} needs to be set).
+
+\code{radius} is the radius of a circle where the function \code{elevation} will be applied. (\code{center} needs to be set).
+
+\code{description} is a description of this operator.
+
+\code{label} is the label used when reporting on the operator.
+
+\code{logging} is the boolean flag to set logging to a file. 
+
+\code{verbose} is the boolean flag to setup extended output from the operator. 
+
+Example:
+
+\begin{verbatim}
+op0 = anuga.Set_elevation_operator(domain, elevation = lambda t : 0.01*t)
+\end{verbatim}
+would setup an operator which raises the elevation over the whole domain 1 cm per second.
+
+
+While 
+\begin{verbatim}
+P = [[x0, y0], [x1, y0], [x1, y1], [x0, y1]] 
+
+op0 = anuga.Set_elevation_operator(domain, \
+                       elevation = lambda t : 0.01*t, polygon=P)
+\end{verbatim}
+would setup an operator which raises the elevation over the polygon \code{P} 1 cm per second.
+
+\end{classdesc}
+
+%--------------------------------------------
+\section{Culvert operators}
+
+A culvert is a structure that allows water to flow under a road, rail road, trail, or similar obstruction that would otherwise build up behind the embankment. But culverts are also used to drain basins, dams or lakes. \anuga now has the ability to model a culvert by transferring flows from one point in the 2D domain to another. It should be noted that flow from the culvert inlet to the outlet occurs instantaneously with no lagging of the flows.
+
+Only box culverts and pipe culverts can be modelled at present using the Boyd Method (Generalised head-discharge equations for culverts, Proceedings of the fourth national local government engineering conference pp. 161-165, Perth, August, 1987).
+
+Usage
+
+In \anuga, a culvert can be modelled as two points (an inlet point and an outlet point) or along two lines (and inlet line and an outlet line) to account for skew.
+
+The user needs to nominate the location of the culvert inlet and outlet, any culvert losses (i.e. inlet, outlet, bends etc) and the culvert dimensions (in metres). 
+
+The user can also model the effects of momentum jetting at the outlet of the culvert and also account for velocity head of the flow at the inlet for more realistic results.
+
+Culvert hydraulics information can be logged by switching logging on or off.
+
+To use the Boyd Method in \anuga, the user can just copy the code below into their script. So if your model run has five pipe culverts, copy the pipe culvert routine five times into your run script and fill in the details for each of your five culverts.
+\begin{verbatim}
+#------------------------------------- -----------------
+# Box Culvert modelled along two lines 
+#------------------------------------- -----------------
+losses = {'inlet':0.5, 'outlet':1.0, 'bend':0.0, 'grate':0.0, 'pier': 0.0, 'other': 0.0} 
+el0 = numpy.array([[297750.2,6181379.1], [297753.5,6181380.9]]) 
+el1 = numpy.array([[297731.8,6181416.4], [297735.1,6181418.3]])  
+culvert = anuga.Boyd_box_operator(domain, 
+    losses=losses, 
+    width=1.5, 
+    height=1.5, 
+    end_points=[ep0, ep1], 
+    use_momentum_jet=True, #False to switch it off
+    use_velocity_head=True, #False to switch it off
+    manning=0.013, #culvert manning's n
+    logging=True, #False to switch it off
+    label='culvert_log_file', 
+    verbose=False)  
+\end{verbatim}
+
+\begin{verbatim}
+#------------------------------------- ------------------
+# Pipe Culvert modelled with a single point 
+#------------------------------------- ------------------
+losses = {'inlet':0.5, 'outlet':1.0, 'bend':0.0, 'grate':0.0, 'pier': 0.0, 'other': 0.0} 
+ep0 = numpy.array([296653.0,6180014.9]) 
+ep1 = numpy.array([296642.5,6180036.3])      
+culvert = anuga.Boyd_pipe_operator(domain, 
+    losses=losses, 
+    diameter=1.5, 
+    end_points=[ep0, ep1], 
+    use_momentum_jet=True, #False to switch it off
+    use_velocity_head=True, #False to switch it off
+    manning=0.013, #culvert manning's n
+    logging=True, #False to switch it off
+    label='culvert_log_file', 
+    verbose=False) 
+\end{verbatim}
+
+\section{User developed operators}
+
+Suppose we want to add water to our domain at a steady rate. we can create a class based on the \code{Operator} class which has a \code{__call__} function to do the required update to the centroid values of the stage variable.
+
+Here is some code to do this:
+\begin{verbatim}
+from anuga import Operator
+
+class My_operator(Operator):
+    """
+   An operator which adds water to the domain at a given rate (m/sec)
+    """
+
+    def __init__(self,
+                 domain,
+                 rate=0.0,
+                 description = None,
+                 label = None,
+                 logging = False,
+                 verbose = False):
+	  """Initialize the user defined operator
+	  """
+
+        Operator.__init__(self, domain, description, label, logging, verbose)
+       
+        self.rate = rate
+
+
+    def __call__(self):
+        """
+        Apply rate to all triangles 
+        """
+
+        timestep = self.domain.get_timestep()
+
+        self.stage_c[:] = self.stage_c[:]  + self.rate*timestep
+\end{verbatim}
+
+And then in your script you would asosciated \code{My_operator} with the \code{domain} with a call 
+\begin{verbatim}
+My_operator(domain,rate=1.0)
+\end{verbatim}
+
+
+%=================================
+\chapter{Querying SWW model output files}
+
+After a model has been run, it is often useful to extract various information from the SWW
+output file (see Section \ref{sec:sww format}). This is typically more convenient than using the
+diagnostics described in Section \ref{sec:diagnostics} which rely on the model running -- something
+that can be very time consuming. The SWW files are easy and quick to read and offer information
+about the model results such as runup heights, time histories of selected quantities,
+flow through cross sections and much more.
+
+\begin{funcdesc}{elevation = get_maximum_inundation_elevation}{filename,
+                                     polygon=None,
+                                     time_interval=None,
+                                     verbose=False}
+Module: \module{shallow_water.data_manager}
+
+Return the highest elevation where depth is positive ($h > 0$).
+
+\code{filename} is the path to a NetCDF SWW file containing \anuga model output.
+
+\code{polygon} restricts the query to the points that lie within the polygon.
+
+\code {time_interval} restricts the query to within the time interval.
+
+Example to find maximum runup elevation:
+
+\begin{verbatim}
+max_runup = get_maximum_inundation_elevation(filename)
+\end{verbatim}
+
+If no inundation is found (within the \code{polygon} and \code{time_interval}, if specified)
+the return value is \code{None}. This indicates "No Runup" or "Everything is dry".
+\end{funcdesc}
+
+\begin{funcdesc}{(x, y) = get_maximum_inundation_location}{filename,
+                                    polygon=None,
+                                    time_interval=None,
+                                    verbose=False}
+Module: \module{shallow_water.data_manager}
+
+Return location (x,y) of the highest elevation where depth is positive ($h > 0$).
+
+\code{filename} is the path to a NetCDF SWW file containing \anuga model output.
+
+\code{polygon} restricts the query to the points that lie within the polygon.
+
+\code {time_interval} restricts the query to within the time interval.
+
+Example to find maximum runup location:
+
+\begin{verbatim}
+max_runup_location = get_maximum_inundation_location(filename)
+\end{verbatim}
+
+If no inundation is found (within the \code{polygon} and \code{time_interval}, if specified)
+the return value is \code{None}. This indicates "No Runup" or "Everything is dry".
+\end{funcdesc}
+
+\begin{funcdesc}{sww2timeseries}{swwfiles,
+                                 gauge_filename,
+                                 production_dirs,
+                                 report=None,
+                                 reportname=None,
+                                 plot_quantity=None,
+                                 generate_fig=False,
+                                 surface=None,
+                                 time_min=None,
+                                 time_max=None,
+                                 output_centroids=False,
+                                 time_thinning=1,
+                                 time_unit=None,
+                                 title_on=None,
+                                 use_cache=False,
+                                 verbose=False}
+Module: \module{abstract_2d_finite_volumes.util}
+
+Read a set of SWW files and plot the time series for the prescribed quantities
+at defined gauge locations and prescribed time range.
+
+\code{swwfiles} is a dictionary of SWW files with keys of \code{label_id}.
+
+\code{gauge_filename} is the path to a file containing gauge data.
+
+THIS NEEDS MORE WORK.  WORK ON FUNCTION __DOC__ STRING, IF NOTHING ELSE!
+\end{funcdesc}
+
+\begin{funcdesc}{(time, Q) = get_flow_through_cross_section}{filename, polyline, verbose=False}
+Module: \module{shallow_water.data_manager}
+
+Obtain flow ($m^3/s$) perpendicular to specified cross section.
+
+\code{filename} is the path to the SWW file.
+
+\code{output_centroids} set to true to sample at the centre of the triangle containing the point.
+This may be useful for debugging. (new in 1.2)
+
+\code{polyline} is the representation of the desired cross section -- it may contain
+multiple sections allowing for complex shapes. Assumes absolute UTM coordinates.
+
+Returns a tuple \code{time, Q} where:
+
+\code{time} is all the stored times in the SWW file.
+
+\code{Q} is a hydrograph of total flow across the given segments for all stored times.
+
+The normal flow is computed for each triangle intersected by the \code{polyline} and
+added up.  If multiple segments at different angles are specified the normal flows
+may partially cancel each other.
+
+Example to find flow through cross section:
+
+\begin{verbatim}
+cross_section = [[x, 0], [x, width]]
+time, Q = get_flow_through_cross_section(filename, cross_section)
+\end{verbatim}
+\end{funcdesc}
+
+
+%=================================
+\chapter{Parallel Simulation}
+
+The examples from the previous chapters were run just using one processor. As you will have noticed, the simulations can take a long time to run, especially if you are using a fine mesh. Indeed as you halve the spacing of the mesh, the number of triangles goes up by a factor of 4 and the timestep halves, so generally halving the spacing of the mesh increases the run time by approximately a factor of 8. 
+
+One way to alleviate this is to run your simulation in parallel and use more processors to spread the computational cost. 
+
+
+\anuga has the option of running in parallel using \mpi. A description of the method used to parallelize \anuga is given in section~\ref{theory:parallel}.
+
+Such jobs are run using the command
+
+\begin{verbatim}
+mpirun -np n python runscript.py
+\end{verbatim}
+where \code{n} is the total number of processors being used for the parallel run. 
+
+Essentially we can expect speedups comparable to the number of cores available. This is measured via scalability. Our current experiments have demonstrated a scalability of 70\% or above when the number of processors are chosen so that the local partitioned meshes contain around 2000 triangles. 
+
+For instance, suppose we have a problem with a mesh of 500,000 triangles, and use 250 processors. Then the mesh will be partitioned into sub meshes of approximately 2000 triangles. 
+We would expect 70\% scalability, and so expect a speedup of $250 \times 0.7 = 175$. 
+
+\section{The Code}
+
+Here is the code for \file{runParallelCairns.py} which is found in the \file{examples/cairns} directory:
+
+\verbatiminputunderscore{../../anuga_core/examples/cairns/runParallelCairns.py}
+
+
+\section{Structure of the Code}
+
+The code is very similar to the sequential code. The same procedures are used to setup the domain, setup the inital conditions, boundary conditions and evolve. 
+
+
+We first import a few procedures  need for the parallel code.
+\begin{verbatim}
+from anuga import distribute, myid, numprocs, finalize, barrier
+\end{verbatim}
+
+\code{myid} returns the id of the current processor running the code. 
+
+\code{numprocs} returns the total number of processors involved in this parallel jobs (the \code{n} in the original \code{mpirun} command. 
+
+\code{finalize} is called at the end of a script to close down the parallel job.
+
+\code{distribute} is used to partition and setup the parallel domains.
+
+\code{barrier} is a command for processors to wait as all the other processor to catchup to this point.
+
+
+
+
+The creation of the \code{domain} is only done on processor 0. Hence we have the structure:
+
+\begin{verbatim}
+#-------------------------------------
+# Do the domain creation on processor 0
+#-------------------------------------
+if myid == 0: 
+    ....
+    domain = ...
+	
+else:
+    domain = None
+\end{verbatim}
+
+We only need to create the original domain on one processor, otherwise we will have multiple copies of the full domain (which will easily eat up our memory). 
+
+Once we have our code{domain} setup we partition it and send the partitions to each of the other processors, via the command
+
+\begin{verbatim}
+#-------------------------------------
+# Now produce parallel domain
+#-------------------------------------
+domain = distribute(domain)
+\end{verbatim}
+
+This takes the \code{domain} on processor 0 and distributes that domain to each of the processors, (it overwrites the full domain on processor 0).  From this point of the code, there is a different domain on each processor, with each domain comunicating with the other domains to ensure required transfer of inforation to allow flow over the combined domains. 
+
+It is important to apply the boundary conditions after the \code{distribute}
+
+
+
+As the partitoined domain evolve, they will store their data to individual sww files, named as \code{domain_name_Pn_m.sww}, where \code{n} is the total number of processors being used and \code{m} is the specific processor id. 
+
+We have a procedure to merge these individual sww files via the command 
+
+\begin{verbatim}
+domain.sww_merge()
+\end{verbatim}
+
+And we close down the parallel job by issuing the command 
+
+\begin{verbatim}
+finalize()
+\end{verbatim}
+
+
+%=================================
+\chapter{Checkpointing}
+
+When running large and long running simulations, it is useful to provide a mechanism to allow the simulation to restarted if it is interrupted mid simulation. Maybe there is a power cut, or when using batch queues there may be a time restriction on the length of any one running job. Then the use of checkpointing becomes useful. 
+
+The idea is that at regular intervals the system makes a copy of the current state of the computation. Then if there is an interruption the computation can be started from the last checkpoint time and not original start time. 
+
+\section{The Code}
+
+Here is the code for \file{runCheckpoint.py}:
+
+\verbatiminputunderscore{../../anuga_core/examples/checkpointing/runCheckpoint.py}
+
+\section{Structure of the Code}
+
+As usual the code needs to import the required modules, and setup parameters and in this case procedures for setting up the stage and elevation. 
+
+Then we use a \code{try: except:} statement, where we try to open any appropriate checkpoint files, and if not successful, create and distribute a domain as usual. In the creation of the domain, we setup checkpointing just before we start the evolve loop. 
+
+
+\begin{verbatim}
+try:	
+    from anuga import load_checkpoint_file
+	
+    domain = load_checkpoint_file(domain_name = domain_name, 
+                                  checkpoint_dir = checkpoint_dir)
+except:
+    # create the domain as usual
+\end{verbatim}
+
+The \code{load\_checkpoint\_file} needs to know where to look for the checkpoint files (by default the current directory) and the domain name (default ``domain'').
+
+When the domain is originally created, we need to setup checkpointing via the command:
+\begin{verbatim}
+    if useCheckpointing:
+        domain.set_checkpointing(checkpoint_time = 5)
+\end{verbatim}
+
+Here we are setting the wall time between saving of checkpoint files to 5 sec (this is just for testing). Normally we would set the checkpointing to something large, say 15 minutes for a run of multiple hours. 
+
+
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+\chapter{\anuga Public Interface}
+\label{ch:interface}
+
+This chapter gives an overview of the features of \anuga available
+to the user at the public interface. These are grouped under the
+following headings, which correspond to the outline of the examples
+described in Chapter \ref{ch:getstarted}:
+\begin{itemize}
+    \item Establishing the Mesh: Section \ref{sec:establishing the mesh}
+    \item Initialising the Domain: Section \ref{sec:initialising the domain}
+%    \item Specifying the Quantities: Section \ref{sec:quantities}
+    \item Initial Conditions: Section \ref{sec:initial conditions}
+    \item Boundary Conditions: Section \ref{sec:boundary conditions}
+    \item Operators: Section \ref{sec:operators}
+    \item Evolution: Section \ref{sec:evolution}
+\end{itemize}
+
+\section{Documentation}\index{Documentation}
+
+The listings here are intended merely to give the reader an idea of what
+each feature is and how it can be used -- they do
+not give full specifications; for these the reader
+may consult the programmer's guide. The code for every function or class contains
+a documentation string, or 'docstring', that specifies the precise
+syntax for its use. This appears immediately after the line
+introducing the code, between two sets of triple quotes.
+
+Python has a handy tool that lets you easily navigate this documentation,
+called \code{pydoc}. In Linux, it runs as a server, which serves the
+documentation up to your web browser:
+
+1. Open a terminal at your \code{anuga} folder
+2. Start the python documentation server with \code{pydoc -p 6767}
+3. Open a browser and type in \code{http://localhost:6767/}
+
+Now you have a real-time programmers' guide for \anuga, and an easy way to find 
+the functions you are interested in. Pydoctor and other Python doc generators
+look nicer and have graphs, etc, but pydoc works straight out of the box. 
+
+\section{Public vs Private Interface}\index{public vs private interface}
+
+To simplify the process of writing scripts, \anuga has a public API which packages
+up the commonly-used functionality of \anuga in the one place. To use it, simply
+import the anuga module like so:
+
+\begin{verbatim}
+import anuga
+\end{verbatim}
+
+You can now use the public API like so. Note the \code{anuga.} prefix:
+
+\begin{verbatim}
+anuga.sww2dem('in.sww', 'out.asc')
+\end{verbatim}
+
+If you wish to delve "under the hood" and modify the way \anuga runs at a more
+advanced level, you need to specify the full location of the module like so:
+
+\begin{verbatim}
+from anuga.fit_interpolate.interpolate import Interpolation_interface
+\end{verbatim}
+
+All modules are in the folder \file{anuga} or one of its subfolders, and the
+location of each module is described relative to \file{anuga}. Rather
+than using pathnames, whose syntax depends on the operating system,
+we use the format adopted for importing the function or class for
+use in Python code. For example, suppose we wish to specify that the
+function \function{create\_mesh\_from\_regions} is in a module called
+\module{mesh\_interface} in a subfolder of \module{anuga} called
+\code{pmesh}. In Linux or Unix syntax, the pathname of the file
+containing the function, relative to \file{anuga}, would be:
+
+\begin{verbatim}
+pmesh/mesh_interface.py
+\end{verbatim}
+
+\label{sec:mesh interface}
+while in Windows syntax it would be:
+
+\begin{verbatim}
+pmesh\mesh_interface.py
+\end{verbatim}
+
+Rather than using either of these forms, in this chapter we specify
+the location simply as \code{anuga.pmesh.mesh_interface}, in keeping with
+the usage in the Python statement for importing the function,
+namely:
+
+\begin{verbatim}
+from anuga.pmesh.mesh_interface import create_mesh_from_regions
+\end{verbatim}
+
+
+The following parameters are common to many functions and classes
+and are omitted from the descriptions given below:
+
+%\begin{tabular}{ll}
+\begin{tabular}{p{2.0cm} p{14.0cm}}
+  \emph{use\_cache} & Specifies whether caching is to be used for improved performance. See Section \ref{sec:caching} for details\\
+                    & on the underlying caching functionality\\                      
+  \emph{verbose}    & If \code{True}, provides detailed terminal output to the user\\
+\end{tabular}
+
+
+
 
 
@@ -3071,599 +3727,10 @@
 
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-\section{Operators}\index{operators}
-\label{sec:operators}
-
-A way to effect the evolution is via fractional step \code{operators}.  The idea is that at each inner timestep we can use 
-an \code{operator} to change the centroid values of quantities such as the stage and the xmomentum and the ymomentum. T
-
-If you are using one of theolder flow argorithms then the elevation can also be changed but more care needs to be applied as the \code{elevation} quantity needs to remain continuous. 
-The moral, play with \code{elevation} at your peril if using the old algorithms. 
-
-The newer algorithms (check that you are using discontinuous elevation with the domain member function \code{<domain>.get_using_discontinuous_elevation()}
-
-Currently predefined operators:
-
-
-\begin{classdesc}{Set_elevation_operator}{domain,
-                 elevation=None,
-                 indices=None,
-                 polygon=None,
-                 center=None,
-                 radius=None,
-                 description = None,
-                 label = None,
-                 logging = False,
-                 verbose = False}
-Module: \module{anuga.operators.set_elevation_operator}
-
-Set the elevation in a region (careful to maintain continuity of elevation)
-
-
-\code{domain} is the domain being evolved.
-
-\code{elevation} a function of $t$ or $x,y$ or $x,y,t$ written to take numpy arrays $x$ and $y$. $t$ is a scalar.
-
-\code{indices} is a list of triangle indices where the function \code{elevation} will be applied.
-
-\code{polygon} is a polygon. The function \code{elevation} will be applied to triangles with in the polygon.
-
-\code{center} is the center of a circle where the function \code{elevation} will be applied. (\code{radius} needs to be set).
-
-\code{radius} is the radius of a circle where the function \code{elevation} will be applied. (\code{center} needs to be set).
-
-\code{description} is a description of this operator.
-
-\code{label} is the label used when reporting on the operator.
-
-\code{logging} is the boolean flag to set logging to a file. 
-
-\code{verbose} is the boolean flag to setup extended output from the operator. 
-
-Example:
-
-\begin{verbatim}
-op0 = anuga.Set_elevation_operator(domain, elevation = lambda t : 0.01*t)
-\end{verbatim}
-would setup an operator which raises the elevation over the whole domain 1 cm per second.
-
-
-While 
-\begin{verbatim}
-P = [[x0, y0], [x1, y0], [x1, y1], [x0, y1]] 
-
-op0 = anuga.Set_elevation_operator(domain, \
-                       elevation = lambda t : 0.01*t, polygon=P)
-\end{verbatim}
-would setup an operator which raises the elevation over the polygon \code{P} 1 cm per second.
-
-\end{classdesc}
-
-%--------------------------------------------
-\subsection{Culvert operators}
-
-A culvert is a structure that allows water to flow under a road, rail road, trail, or similar obstruction that would otherwise build up behind the embankment. But culverts are also used to drain basins, dams or lakes. \anuga now has the ability to model a culvert by transferring flows from one point in the 2D domain to another. It should be noted that flow from the culvert inlet to the outlet occurs instantaneously with no lagging of the flows.
-
-Only box culverts and pipe culverts can be modelled at present using the Boyd Method (Generalised head-discharge equations for culverts, Proceedings of the fourth national local government engineering conference pp. 161-165, Perth, August, 1987).
-
-Usage
-
-In \anuga, a culvert can be modelled as two points (an inlet point and an outlet point) or along two lines (and inlet line and an outlet line) to account for skew.
-
-The user needs to nominate the location of the culvert inlet and outlet, any culvert losses (i.e. inlet, outlet, bends etc) and the culvert dimensions (in metres). 
-
-The user can also model the effects of momentum jetting at the outlet of the culvert and also account for velocity head of the flow at the inlet for more realistic results.
-
-Culvert hydraulics information can be logged by switching logging on or off.
-
-To use the Boyd Method in \anuga, the user can just copy the code below into their script. So if your model run has five pipe culverts, copy the pipe culvert routine five times into your run script and fill in the details for each of your five culverts.
-\begin{verbatim}
-#------------------------------------- -----------------
-# Box Culvert modelled along two lines 
-#------------------------------------- -----------------
-losses = {'inlet':0.5, 'outlet':1.0, 'bend':0.0, 'grate':0.0, 'pier': 0.0, 'other': 0.0} 
-el0 = numpy.array([[297750.2,6181379.1], [297753.5,6181380.9]]) 
-el1 = numpy.array([[297731.8,6181416.4], [297735.1,6181418.3]])  
-culvert = anuga.Boyd_box_operator(domain, 
-    losses=losses, 
-    width=1.5, 
-    height=1.5, 
-    end_points=[ep0, ep1], 
-    use_momentum_jet=True, #False to switch it off
-    use_velocity_head=True, #False to switch it off
-    manning=0.013, #culvert manning's n
-    logging=True, #False to switch it off
-    label='culvert_log_file', 
-    verbose=False)  
-\end{verbatim}
-
-\begin{verbatim}
-#------------------------------------- ------------------
-# Pipe Culvert modelled with a single point 
-#------------------------------------- ------------------
-losses = {'inlet':0.5, 'outlet':1.0, 'bend':0.0, 'grate':0.0, 'pier': 0.0, 'other': 0.0} 
-ep0 = numpy.array([296653.0,6180014.9]) 
-ep1 = numpy.array([296642.5,6180036.3])      
-culvert = anuga.Boyd_pipe_operator(domain, 
-    losses=losses, 
-    diameter=1.5, 
-    end_points=[ep0, ep1], 
-    use_momentum_jet=True, #False to switch it off
-    use_velocity_head=True, #False to switch it off
-    manning=0.013, #culvert manning's n
-    logging=True, #False to switch it off
-    label='culvert_log_file', 
-    verbose=False) 
-\end{verbatim}
-
-\subsection{User developed operators}
-
-Suppose we want to add water to our domain at a steady rate. we can create a class based on the \code{Operator} class which has a \code{__call__} function to do the required update to the centroid values of the stage variable.
-
-Here is some code to do this:
-\begin{verbatim}
-from anuga import Operator
-
-class My_operator(Operator):
-    """
-   An operator which adds water to the domain at a given rate (m/sec)
-    """
-
-    def __init__(self,
-                 domain,
-                 rate=0.0,
-                 description = None,
-                 label = None,
-                 logging = False,
-                 verbose = False):
-	  """Initialize the user defined operator
-	  """
-
-        Operator.__init__(self, domain, description, label, logging, verbose)
-       
-        self.rate = rate
-
-
-    def __call__(self):
-        """
-        Apply rate to all triangles 
-        """
-
-        timestep = self.domain.get_timestep()
-
-        self.stage_c[:] = self.stage_c[:]  + self.rate*timestep
-\end{verbatim}
-
-And then in your script you would asosciated \code{My_operator} with the \code{domain} with a call 
-\begin{verbatim}
-My_operator(domain,rate=1.0)
-\end{verbatim}
-
-
-
-
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-\section{Evolution}\index{evolution}
-\label{sec:evolution}
-
-\begin{methoddesc}{\emph{<domain>}.evolve}{yieldstep=None,
-                                           finaltime=None,
-                                           duration=None,
-                                           skip_initial_step=False}
-Module: \module{abstract_2d_finite_volumes.domain}
-
-This method is invoked once all the
-preliminaries have been completed, and causes the model to progress
-through successive steps in its evolution, storing results and
-outputting statistics whenever a user-specified period
-\code{yieldstep} is completed.  Generally during this period the
-model will evolve through several steps internally
-as the method forces the water speed to be calculated
-on successive new cells. 
-
-\code{yieldstep} is the interval in seconds between yields where results are
-stored, statistics written and the domain is inspected or possibly modified.
-If omitted an internal predefined \code{yieldstep} is used.  Internally, smaller
-timesteps may be taken.
-
-\code{duration} is the duration of the simulation in seconds.
-
-\code{finaltime} is the time in seconds where simulation should end. This is currently
-relative time, so it's the same as \code{duration}.  If both \code{duration} and
-\code{finaltime} are given an exception is thrown.
-
-\code{skip_initial_step} is a boolean flag that decides whether the first
-yield step is skipped or not. This is useful for example to avoid
-duplicate steps when multiple evolve processes are dove tailed.
-
-The code specified by the user in the block following the evolve statement is
-only executed once every \code{yieldstep} even though 
-\anuga typically will take many more internal steps behind the scenes.
-
-You can include \method{evolve} in a statement of the type:
-
-\begin{verbatim}
-for t in domain.evolve(yieldstep, finaltime):
-    <Do something with domain and t>
-\end{verbatim}
-\end{methoddesc}
-
-\subsection{Diagnostics}
-\label{sec:diagnostics}
-
-\begin{methoddesc}{\emph{<domain>}.statistics}{}
-Module: \module{abstract\_2d\_finite\_volumes.domain}
-
-Outputs statistics about the mesh within the \code{Domain}.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<domain>}.timestepping_statistics}{track_speeds=False, triangle_id=None}
-Module: \module{abstract_2d_finite_volumes.domain}
-
-Returns a string of the following type for each timestep:\\
-\code{Time = 0.9000, delta t in [0.00598964, 0.01177388], steps=12 (12)}
-
-Here the numbers in \code{steps=12 (12)} indicate the number of steps taken and
-the number of first-order steps, respectively.
-
-The optional keyword argument \code{track_speeds} will
-generate a histogram of speeds generated by each triangle if set to \code{True}. The
-speeds relate to the size of the timesteps used by \anuga and
-this diagnostics may help pinpoint problem areas where excessive speeds
-are generated.
-
-The optional keyword argument \code{triangle_id} can be used to specify a particular
-triangle rather than the one with the largest speed.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<domain>}.boundary_statistics}{quantities=None,
-                                                      tags=None}
-Module: \module{abstract_2d_finite_volumes.domain}
-
-Generates output about boundary forcing at each timestep.
-
-\code{quantities} names the quantities to be reported -- may be \code{None},
-a string or a list of strings.
-
-\code{tags} names the tags to be reported -- may be either None, a string or a list of strings.
-
-When \code{quantities = 'stage'} and \code{tags = ['top', 'bottom']}
-will return a string like:
-
-\begin{verbatim}
-Boundary values at time 0.5000:
-    top:
-        stage in [ -0.25821218,  -0.02499998]
-    bottom:
-        stage in [ -0.27098821,  -0.02499974]
-\end{verbatim}
-\end{methoddesc}
-
-\begin{methoddesc}{Q = \emph{<domain>}.get_quantity}{name,
-                                               location='vertices',
-                                               indices=None}
-Module: \module{abstract_2d_finite_volumes.domain}
-
-This function returns a Quantity object Q.
-Access to its values should be done through \code{Q.get_values()} documented on Page \pageref{pg:get values}.
-
-\code{name} is the name of the quantity to retrieve.
-
-\code{location} is
-
-\code{indices} is
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<domain>}.set_quantities_to_be_monitored}{quantity,
-                                             polygon=None,
-                                             time_interval=None}
-Module: \module{abstract\_2d\_finite\_volumes.domain}
-
-Selects quantities and derived quantities for which extrema attained at internal timesteps
-will be collected.
-
-\code{quantity} specifies the quantity or quantities to be monitored and must be either:
-\begin{itemize}
-  \item the name of a quantity or derived quantity such as 'stage-elevation',
-  \item a list of quantity names, or
-  \item \code{None}.
-\end{itemize}
-
-\code{polygon} can be used to monitor only triangles inside the polygon. Otherwise
-all triangles will be included.
-
-\code{time_interval} will restrict monitoring to time steps in the interval. Otherwise
-all timesteps will be included.
-
-Information can be tracked in the evolve loop by printing \code{quantity_statistics} and
-collected data will be stored in the SWW file.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<domain>}.quantity_statistics}{precision='\%.4f'}
-Module: \module{abstract_2d_finite_volumes.domain}
-
-Reports on extrema attained by selected quantities.
-
-Returns a string of the following type for each timestep:
-
-\begin{verbatim}
-Monitored quantities at time 1.0000:
-  stage-elevation:
-    values since time = 0.00 in [0.00000000, 0.30000000]
-    minimum attained at time = 0.00000000, location = (0.16666667, 0.33333333)
-    maximum attained at time = 0.00000000, location = (0.83333333, 0.16666667)
-  ymomentum:
-    values since time = 0.00 in [0.00000000, 0.06241221]
-    minimum attained at time = 0.00000000, location = (0.33333333, 0.16666667)
-    maximum attained at time = 0.22472667, location = (0.83333333, 0.66666667)
-  xmomentum:
-    values since time = 0.00 in [-0.06062178, 0.47886313]
-    minimum attained at time = 0.00000000, location = (0.16666667, 0.33333333)
-    maximum attained at time = 0.35103646, location = (0.83333333, 0.16666667)
-\end{verbatim}
-
-The quantities (and derived quantities) listed here must be selected at model
-initialisation time using the method \code{domain.set_quantities_to_be_monitored()}.
-
-The optional keyword argument \code{precision='\%.4f'} will
-determine the precision used for floating point values in the output.
-This diagnostics helps track extrema attained by the selected quantities
-at every internal timestep.
-
-These values are also stored in the SWW file for post-processing.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<quantity>}.get_values}{interpolation_points=None,
-                   location='vertices',
-                   indices=None,
-                   use_cache=False,
-                   verbose=False}
-\label{pg:get values}
-Module: \module{abstract_2d_finite_volumes.quantity}
-
-Extract values for quantity as a numeric array.
-
-\code{interpolation_points} is a list of (x, y) coordinates where the value is
-sought (using interpolation). If \code{interpolation_points} is given, values
-for \code{location} and \code{indices} are ignored.
-Assume either an absolute UTM coordinates or geospatial data object.
-    
-\code{location} specifies where values are to be stored.
-Permissible options are 'vertices', 'edges', 'centroids' or 'unique vertices'.
-
-The returned values will have the leading dimension equal to length of the \code{indices} list or
-N (all values) if \code{indices} is \code{None}.
-
-If \code{location} is 'centroids' the dimension of returned
-values will be a list or a numeric array of length N, N being
-the number of elements.
-      
-If \code{location} is 'vertices' or 'edges' the dimension of
-returned values will be of dimension \code{Nx3}.
-
-If \code{location} is 'unique vertices' the average value at
-each vertex will be returned and the dimension of returned values
-will be a 1d array of length "number of vertices" 
-      
-\code{indices} is the set of element ids that the operation applies to.
-
-The values will be stored in elements following their internal ordering.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<quantity>}.set_values}{numeric=None,
-                               quantity=None,
-                               function=None,
-                               geospatial_data=None,
-                               filename=None,
-                               attribute_name=None,
-                               alpha=None,
-                               location='vertices',
-                               polygon=None,
-                               indices=None,
-                               smooth=False,
-                               verbose=False,
-                               use_cache=False}
-Module: \module{abstract_2d_finite_volumes.quantity}
-
-Assign values to a quantity object.
-
-This method works the same way as \code{set_quantity()} except that it doesn't take
-a quantity name as the first argument since it is applied directly to the quantity.
-Use \code{set_values} is used to assign
-values to a new quantity that has been created but which is
-not part of the domain's predefined quantities.
-
-\code{location} is ??
-
-\code{indices} is ??
-
-The method \code{set_values()} is always called by \code{set_quantity()}
-behind the scenes.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<quantity>}.get_integral}{}
-Module: \module{abstract_2d_finite_volumes.quantity}
-
-Return the computed integral over the entire domain for the quantity.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<quantity>}.get_maximum_value}{indices=None}
-Module: \module{abstract_2d_finite_volumes.quantity}
-
-Return the maximum value of a quantity (on centroids).
-
-\code{indices} is the optional set of element \code{id}s that
-the operation applies to.
-
-We do not seek the maximum at vertices as each vertex can
-have multiple values -- one for each triangle sharing it.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<quantity>}.get_maximum_location}{indices=None}
-Module: \module{abstract_2d_finite_volumes.quantity}
-
-Return the location of the maximum value of a quantity (on centroids).
-
-\code{indices} is the optional set of element \code{id}s that
-the operation applies to.
-
-We do not seek the maximum at vertices as each vertex can
-have multiple values -- one for each triangle sharing it.
-
-If there are multiple cells with the same maximum value, the
-first cell encountered in the triangle array is returned.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<domain>}.get_wet_elements}{indices=None}
-Module: \module{shallow_water.shallow_water_domain}
-
-Returns the indices for elements where h $>$ minimum_allowed_height
-
-\code{indices} is the optional set of element \code{id}s that
-the operation applies to.
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<domain>}.get_maximum_inundation_elevation}{indices=None}
-Module: \module{shallow_water.shallow_water_domain}
-
-Return highest elevation where h $>$ 0.
-
-\code{indices} is the optional set of element \code{id}s that
-the operation applies to.
-
-Example to find maximum runup elevation:
-\begin{verbatim}
-z = domain.get_maximum_inundation_elevation()
-\end{verbatim}
-\end{methoddesc}
-
-\begin{methoddesc}{\emph{<domain>}.get_maximum_inundation_location}{indices=None}
-Module: \module{shallow_water.shallow_water_domain}
-
-Return location (x,y) of highest elevation where h $>$ 0.
-
-\code{indices} is the optional set of element \code{id}s that
-the operation applies to.
-
-Example to find maximum runup location:
-\begin{verbatim}
-x, y = domain.get_maximum_inundation_location()
-\end{verbatim}
-\end{methoddesc}
-
-
-\section{Queries of SWW model output files}
-After a model has been run, it is often useful to extract various information from the SWW
-output file (see Section \ref{sec:sww format}). This is typically more convenient than using the
-diagnostics described in Section \ref{sec:diagnostics} which rely on the model running -- something
-that can be very time consuming. The SWW files are easy and quick to read and offer information
-about the model results such as runup heights, time histories of selected quantities,
-flow through cross sections and much more.
-
-\begin{funcdesc}{elevation = get_maximum_inundation_elevation}{filename,
-                                     polygon=None,
-                                     time_interval=None,
-                                     verbose=False}
-Module: \module{shallow_water.data_manager}
-
-Return the highest elevation where depth is positive ($h > 0$).
-
-\code{filename} is the path to a NetCDF SWW file containing \anuga model output.
-
-\code{polygon} restricts the query to the points that lie within the polygon.
-
-\code {time_interval} restricts the query to within the time interval.
-
-Example to find maximum runup elevation:
-
-\begin{verbatim}
-max_runup = get_maximum_inundation_elevation(filename)
-\end{verbatim}
-
-If no inundation is found (within the \code{polygon} and \code{time_interval}, if specified)
-the return value is \code{None}. This indicates "No Runup" or "Everything is dry".
-\end{funcdesc}
-
-\begin{funcdesc}{(x, y) = get_maximum_inundation_location}{filename,
-                                    polygon=None,
-                                    time_interval=None,
-                                    verbose=False}
-Module: \module{shallow_water.data_manager}
-
-Return location (x,y) of the highest elevation where depth is positive ($h > 0$).
-
-\code{filename} is the path to a NetCDF SWW file containing \anuga model output.
-
-\code{polygon} restricts the query to the points that lie within the polygon.
-
-\code {time_interval} restricts the query to within the time interval.
-
-Example to find maximum runup location:
-
-\begin{verbatim}
-max_runup_location = get_maximum_inundation_location(filename)
-\end{verbatim}
-
-If no inundation is found (within the \code{polygon} and \code{time_interval}, if specified)
-the return value is \code{None}. This indicates "No Runup" or "Everything is dry".
-\end{funcdesc}
-
-\begin{funcdesc}{sww2timeseries}{swwfiles,
-                                 gauge_filename,
-                                 production_dirs,
-                                 report=None,
-                                 reportname=None,
-                                 plot_quantity=None,
-                                 generate_fig=False,
-                                 surface=None,
-                                 time_min=None,
-                                 time_max=None,
-                                 output_centroids=False,
-                                 time_thinning=1,
-                                 time_unit=None,
-                                 title_on=None,
-                                 use_cache=False,
-                                 verbose=False}
-Module: \module{abstract_2d_finite_volumes.util}
-
-Read a set of SWW files and plot the time series for the prescribed quantities
-at defined gauge locations and prescribed time range.
-
-\code{swwfiles} is a dictionary of SWW files with keys of \code{label_id}.
-
-\code{gauge_filename} is the path to a file containing gauge data.
-
-THIS NEEDS MORE WORK.  WORK ON FUNCTION __DOC__ STRING, IF NOTHING ELSE!
-\end{funcdesc}
-
-\begin{funcdesc}{(time, Q) = get_flow_through_cross_section}{filename, polyline, verbose=False}
-Module: \module{shallow_water.data_manager}
-
-Obtain flow ($m^3/s$) perpendicular to specified cross section.
-
-\code{filename} is the path to the SWW file.
-
-\code{output_centroids} set to true to sample at the centre of the triangle containing the point.
-This may be useful for debugging. (new in 1.2)
-
-\code{polyline} is the representation of the desired cross section -- it may contain
-multiple sections allowing for complex shapes. Assumes absolute UTM coordinates.
-
-Returns a tuple \code{time, Q} where:
-
-\code{time} is all the stored times in the SWW file.
-
-\code{Q} is a hydrograph of total flow across the given segments for all stored times.
-
-The normal flow is computed for each triangle intersected by the \code{polyline} and
-added up.  If multiple segments at different angles are specified the normal flows
-may partially cancel each other.
-
-Example to find flow through cross section:
-
-\begin{verbatim}
-cross_section = [[x, 0], [x, width]]
-time, Q = get_flow_through_cross_section(filename, cross_section)
-\end{verbatim}
-\end{funcdesc}
+
+
+
+
+
+
 
 
@@ -3692,7 +3759,6 @@
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
-\chapter{\anuga System Architecture}
-
-\section{File Formats}
+\chapter{\anuga File Formats}
+
 \label{sec:file formats}
 
@@ -3732,5 +3798,5 @@
 %SOMETHING MISSING HERE!?
 
-\subsection{Manually Created Files}
+\section{Manually Created Files}
 
 \begin{tabular}{ll}
@@ -3739,5 +3805,5 @@
 \end{tabular}
 
-\subsection{Automatically Created Files}
+\section{Automatically Created Files}
 
 \begin{tabular}{ll}
@@ -3752,5 +3818,5 @@
 \bigskip
 
-\subsection{SWW, STS and TMS Formats}
+\section{SWW, STS and TMS Formats}
 \label{sec:sww format}
 The SWW, STS and TMS formats are all NetCDF formats and are of key importance for \anuga.
@@ -3780,5 +3846,5 @@
 
 %FIXME (Ole): Should put in example with nonzero xllcorner, yllcorner
-\verbatiminputunderscore{examples/bedslopeexcerpt.cdl}
+\verbatiminputunderscore{../../anuga_core/examples/bedslopeexcerpt.cdl}
 
 The SWW format is used not only for output but also serves as input
@@ -3806,5 +3872,5 @@
 A TMS file is used to store time series data that is independent of position.
 
-\subsection{Mesh File Formats}
+\section{Mesh File Formats}
 
 A mesh file is a file that has a specific format suited to
@@ -3831,5 +3897,5 @@
 offset to be applied to $x$ and $y$ values -- e.g.\ to the vertices.
 
-\subsection{Formats for Storing Arbitrary Points and Attributes}
+\section{Formats for Storing Arbitrary Points and Attributes}
 
 A CSV/TXT/XYA file is used to store data representing
@@ -3859,5 +3925,5 @@
 attribute.
 
-\subsection{ArcView Formats}
+\section{ArcView Formats}
 
 Files of the three formats ASC, PRJ and ERS are all associated with
@@ -3882,11 +3948,11 @@
 to represent metadata for a DEM.
 
-\subsection{DEM Format}
+\section{DEM Format}
 
 A DEM file in \anuga is a NetCDF representation of regular DEM data.
 
-\subsection{Other Formats}
-
-\subsection{Basic File Conversions}
+
+
+\section{Basic File Conversions}
 \label{sec:basicfileconversions}
 
@@ -5574,5 +5640,5 @@
 \end{methoddesc}
 
-\subsection{Miscellaneous Functions}
+\section{Miscellaneous Functions}
 
 The functions here are not \code{Geospatial_data} object methods, but are used with them.
Index: /trunk/anuga_user_manual/source/copyright.tex
===================================================================
--- /trunk/anuga_user_manual/source/copyright.tex	(revision 9726)
+++ /trunk/anuga_user_manual/source/copyright.tex	(revision 9727)
@@ -73,5 +73,5 @@
 \begin{itemize}
 \item \anuga was developed by Stephen Roberts, Ole Nielsen, Duncan Gray and Jane Sexton. It is currently being developed and
- maintained by Nariman Habili and Stephen Roberts.
+ maintained by Stephen Roberts and Gareth Davies.
 \index{ANUGA!credits|textit}  
 \end{itemize}
@@ -87,5 +87,5 @@
 \item Ole Nielsen, James Hudson, John Jakeman, Rudy van Drie, Ted Rigby, 
       Petar Milevski, Joaquim Luis, Nils Goseberg, William Power,
-      Trevor Dhu, Linda Stals, Matt Hardy, Jack Kelly and Christopher
+      Trevor Dhu, Nariman Habili, Gareth Davies, Linda Stals, Matt Hardy, Jack Kelly and Christopher
       Zoppou who contributed to this project at various times.
 \item A stand alone visualiser (anuga\_viewer) based on Open-scene-graph was developed by Darran Edmundson and James Hudson.
Index: /trunk/anuga_user_manual/source/update_anuga_user_manual.py
===================================================================
--- /trunk/anuga_user_manual/source/update_anuga_user_manual.py	(revision 9726)
+++ /trunk/anuga_user_manual/source/update_anuga_user_manual.py	(revision 9727)
@@ -41,5 +41,5 @@
 from os.path import expanduser, split, join
 from anuga.utilities.system_tools import get_revision_number, get_pathname_from_package
-from anuga.config import major_revision
+from anuga import __version__ as major_revision
 from sys import argv
 
@@ -115,5 +115,5 @@
 
 system('cp *.pdf ../')
-system('mv *.pdf ../../anuga_core/user_manual')
+system('cp anuga_user_manual.pdf ../../anuga_core/doc')
 
 print 'Cleanup aux tex files'
Index: /trunk/anuga_user_manual/source/version.tex
===================================================================
--- /trunk/anuga_user_manual/source/version.tex	(revision 9726)
+++ /trunk/anuga_user_manual/source/version.tex	(revision 9727)
@@ -7,3 +7,3 @@
 % release version; this is used to define the
 % \version macro
-\release{1.3.1}
+\release{1.3.11}
