Changeset 9727


Ignore:
Timestamp:
Apr 28, 2015, 10:01:09 AM (9 years ago)
Author:
steve
Message:

Updating manual to include references to github

Location:
trunk/anuga_user_manual/source
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/anuga_user_manual/source/anuga_user_manual.tex

    r9665 r9727  
    129129cells, the topography and bathymetry, frictional resistance, initial
    130130values for water level (called \emph{stage}\index{stage} within \anuga),
    131 boundary conditions and forces such as rainfall, stream flows, windstress or pressure gradients if applicable.
     131boundary conditions and operators such as rainfall,
     132stream flows, windstress or pressure gradients if applicable.
    132133
    133134\anuga tracks the evolution of water depth and horizontal momentum
     
    160161inundation software system, describe what it can do and give step-by-step
    161162instructions for setting up and running hydrodynamic simulations.
    162 The stable release of \anuga and this manual are available on sourceforge at
    163 \url{http://sourceforge.net/projects/anuga}. A snapshot of work in progress is
    164 available through the \anuga software repository at
    165 \url{https://anuga.anu.edu.au/svn/anuga/trunk/anuga_core/source/anuga}
    166 where the more adventurous reader might like to go.
     163The stable releases of \anuga and this manual are available on github at
     164\url{http://github.com/GeoscienceAustralia/anuga_core/releases}.
     165
     166A snapshot of the current work in progress is at 
     167\url{http://github.com/GeoscienceAustralia/anuga_core}.
    167168
    168169This manual describes \anuga version \version. To check for later versions of this manual
    169 go to \url{https://anuga.anu.edu.au}.
     170go to \url{http://github.com/GeoscienceAustralia/anuga_core/doc/anuga_user_manual.pdf}.
    170171
    171172\section{Scope}
     
    178179
    179180The latest installation instructions may be found at:
    180 \url{http://anuga.anu.edu.au/raw-attachment/wiki/WikiStart/anuga_installation_guide-1.2.1.pdf}.
     181\url{http://github.com/GeoscienceAustralia/anuga_core/INSTALL.rst}.
    181182
    182183\section{Audience}
     
    184185Readers are assumed to be familiar with the Python Programming language and
    185186its object oriented approach.
    186 Python tutorials include
    187 \url{http://docs.python.org/tut} and \url{http://www.sthurlow.com/python}.
     187A good Python tutorial is
     188\url{https://docs.python.org/2/tutorial}.
    188189
    189190Readers also need to have a general understanding of scientific modelling,
     
    195196%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    196197
    197 \chapter{Background}
     198\section{Background}
    198199
    199200Modelling the effects on the built environment of natural hazards such
     
    217218method to accommodate discontinuities in the solution\footnote{
    218219While \anuga works with discontinuities in the conserved quantities stage,
    219 xmomentum and ymomentum, it does not allow discontinuities in the bed elevation.}.
     220xmomentum and ymomentum, it does not allow discontinuities in the bed elevation.
     221The next version 1.4 works with discontinuous bed elevation, and the current version
     222can use discontinuous bed if the flow algorithm is set to \texttt{DE0} or \texttt{DE1}}.
    220223
    221224To set up a particular scenario the user specifies the geometry
    222225(bathymetry and topography), the initial water level (stage),
    223 boundary conditions such as tide, and any forcing terms that may
    224 drive the system such as rainfall, abstraction of water, wind stress or atmospheric pressure
    225 gradients. Gravity and frictional resistance from the different
    226 terrains in the model are represented by predefined forcing terms.
    227 See section \ref{sec:forcing terms} for details on forcing terms available in \anuga.
     226boundary conditions such as tide, and any operators  that may
     227drive the system such as rainfall, abstraction of water,  erosion, culverts
     228See section \ref{sec:operators} for details of operators available in \anuga.
    228229
    229230The built-in mesh generator, called \code{graphical\_mesh\_generator},
     
    238239and can be readily adapted to changing requirements throughout its
    239240lifetime.  Computationally intensive components are written for
    240 efficiency in C routines working directly with Python numeric
    241 structures.  The animation tool developed for \anuga is based on
     241efficiency in C routines working directly with Python numpy
     242structures. 
     243
     244The visualisation tool developed for \anuga is based on
    242245OpenSceneGraph, an Open Source Software (OSS) component allowing high
    243246level interaction with sophisticated graphics primitives.
    244247See \cite{nielsen2005} for more background on \anuga.
    245248
    246 \chapter{Restrictions and limitations on \anuga}
     249
     250 
     251
     252 
     253
     254
     255\section{Restrictions and limitations on \anuga}
    247256\label{ch:limitations}
    248257
     
    270279\end{itemize}
    271280
     281\pagebreak
     282\section{Citing \anuga}
     283
     284When citing \anuga cite this manual:
     285
     286Bibtex entry:
     287\begin{verbatim}
     288@manual{anugamanual,
     289  title={ANUGA User Manual},
     290  author={Roberts, S and Nielsen, O. and Gray, D. and Sexton, J.},
     291  organization={Geoscience Australia},
     292  year={2015}
     293}
     294\end{verbatim}
     295
     296or this article
     297
     298Bibtex entry:
     299\begin{verbatim}
     300@article{nielsen2005hydrodynamic,
     301  title={Hydrodynamic modelling of coastal inundation},
     302  author={Nielsen, O and Roberts, S and Gray, D and McPherson, A and Hitchman, A},
     303  journal={MODSIM 2005 International Congress on Modelling and Simulation},
     304  pages={518--523},
     305  year={2005}
     306}
     307\end{verbatim}
     308
     309
     310
    272311%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    273312
     
    283322\section{A Simple Example}
    284323\label{sec:simpleexample}
    285 
    286 \subsection{Overview}
    287324
    288325What follows is a discussion of the structure and operation of a
     
    320357\end{verbatim}
    321358
    322 \subsection{Outline of the Program}
     359\section{Outline of the Program}
    323360
    324361In outline, \file{runup.py} performs the following steps:
     
    336373\end{enumerate}
    337374
    338 \subsection{The Code}
     375\section{The Code}
    339376
    340377For reference we include below the complete code listing for
     
    347384
    348385
    349 \subsection{Establishing the Domain}\index{domain, establishing}
     386\section{Establishing the Domain}\index{domain, establishing}
    350387
    351388The very first thing to do is import the various modules, of which the
     
    376413   \item a list \code{vertices} specifying the three vertices of each triangle, and
    377414   \item a dictionary \code{boundary} that stores the edges on
    378          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.
     415         the boundary and associates with each a symbolic tag.
     416         The edges are represented as pairs (i, j) where i refers
     417         to the triangle id and j to the edge id of that triangle.
    379418         Edge ids are enumerated from 0 to 2 based on the id of the vertex opposite.
    380419\end{itemize}
     
    411450However, this is not necessary, as the above is the default behaviour.
    412451
    413 \subsection{Initial Conditions}
     452\section{Initial Conditions}
    414453
    415454The next task is to specify a number of quantities that we wish to
     
    427466forms in which we can invoke \method{set\_quantity}.
    428467
    429 \subsubsection{Elevation}
     468\subsection{Elevation}
    430469
    431470The elevation, or height of the bed, is set using a function
     
    454493compatible. For example, using square root will not work.
    455494
    456 \subsubsection{Friction}
     495\subsection{Friction}
    457496
    458497The assignment of the friction quantity (a forcing term)
     
    467506at every mesh point.
    468507
    469 \subsubsection{Stage}
     508\subsection{Stage}
    470509
    471510The stage (the height of the water surface) is related to the
     
    509548details.
    510549
    511 \subsection{Boundary Conditions}\index{boundary conditions}
     550\section{Boundary Conditions}\index{boundary conditions}
    512551
    513552The boundary conditions are specified as follows:
     
    629668\end{verbatim}
    630669
    631 \subsection{Evolution}\index{evolution}
     670\section{Evolution}\index{evolution}
    632671
    633672The final statement:
     
    644683outputs.  Behind the scenes more time steps are generally taken.
    645684
    646 \subsection{Output}
     685\section{Output}
    647686
    648687The output is a NetCDF file with the extension \code{.sww}. It
     
    701740
    702741
    703 \section{A slightly more complex example}
     742%====================================================
     743\chapter{A slightly more complex example}
    704744\label{sec:channelexample}
    705745
    706 \subsection{Overview}
    707746
    708747The next example is about water-flow in a channel with varying boundary conditions and
     
    711750The example will be built up through three progressively more complex scripts.
    712751
    713 \subsection{Overview}
     752\section{Overview}
    714753
    715754As in the case of \file{runup.py}, the actions carried
     
    727766\end{enumerate}
    728767
    729 \subsection{The Code}
     768\section{The Code}
    730769
    731770Here is the code for the first version of the channel flow \file{channel1.py}:
     
    736775given above, discussing each major step of the code in turn.
    737776
    738 \subsection{Establishing the Mesh}\index{mesh, establishing}
     777\section{Establishing the Mesh}\index{mesh, establishing}
    739778
    740779In this example we use a similar simple structured triangular mesh
     
    796835Here is the code for the second version of the channel flow \file{channel2.py}:
    797836
    798 \verbatiminputunderscore{demos/channel2.py}
     837\verbatiminputunderscore{../../anuga_core/examples/channel2.py}
    799838
    800839This example differs from the first version in that a constant outflow boundary condition has
     
    830869used to model the breakdown of a sluice door when water exceeds a certain level.
    831870
    832 \subsection{Output}
     871\section{Output}
    833872
    834873The text output from this example looks like this:
     
    845884\end{verbatim}
    846885
    847 \subsection{Flow through more complex topographies}
     886\section{Flow through more complex topographies}
    848887
    849888Here is the code for the third version of the channel flow \file{channel3.py}:
     
    13511390\end{figure}
    13521391
    1353 %=================================
    1354 \chapter{Parallel Simulation}
    1355 
    1356 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.
    1357 
    1358 One way to alleviate this is to run your simulation in parallel and use more processors to spread the computational cost.
    1359 
    1360 
    1361 \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}.
    1362 
    1363 Such jobs are run using the command
    1364 
    1365 \begin{verbatim}
    1366 mpirun -np n python runscript.py
    1367 \end{verbatim}
    1368 where \code{n} is the total number of processors being used for the parallel run.
    1369 
    1370 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.
    1371 
    1372 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.
    1373 We would expect 70\% scalability, and so expect a speedup of $250 \times 0.7 = 175$.
    1374 
    1375 \section{The Code}
    1376 
    1377 Here is the code for \file{runParallelCairns.py} which is found in the \file{demos/cairns} directory:
    1378 
    1379 \verbatiminputunderscore{../../anuga_core/demos/cairns/runParallelCairns.py}
    1380 
    1381 
    1382 \section{Structure of the Code}
    1383 
    1384 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.
    1385 
    1386 
    1387 We first import a few procedures  need for the parallel code.
    1388 \begin{verbatim}
    1389 from anuga import distribute, myid, numprocs, finalize, barrier
    1390 \end{verbatim}
    1391 
    1392 \code{myid} returns the id of the current processor running the code.
    1393 
    1394 \code{numprocs} returns the total number of processors involved in this parallel jobs (the \code{n} in the original \code{mpirun} command.
    1395 
    1396 \code{finalize} is called at the end of a script to close down the parallel job.
    1397 
    1398 \code{distribute} is used to partition and setup the parallel domains.
    1399 
    1400 \code{barrier} is a command for processors to wait as all the other processor to catchup to this point.
    1401 
    1402 
    1403 
    1404 
    1405 The creation of the \code{domain} is only done on processor 0. Hence we have the structure:
    1406 
    1407 \begin{verbatim}
    1408 #-------------------------------------
    1409 # Do the domain creation on processor 0
    1410 #-------------------------------------
    1411 if myid == 0:
    1412     ....
    1413     domain = ...
    1414        
    1415 else:
    1416     domain = None
    1417 \end{verbatim}
    1418 
    1419 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).
    1420 
    1421 Once we have our code{domain} setup we partition it and send the partitions to each of the other processors, via the command
    1422 
    1423 \begin{verbatim}
    1424 #-------------------------------------
    1425 # Now produce parallel domain
    1426 #-------------------------------------
    1427 domain = distribute(domain)
    1428 \end{verbatim}
    1429 
    1430 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.
    1431 
    1432 It is important to apply the boundary conditions after the \code{distribute}
    1433 
    1434 
    1435 
    1436 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.
    1437 
    1438 We have a procedure to merge these individual sww files via the command
    1439 
    1440 \begin{verbatim}
    1441 domain.sww_merge()
    1442 \end{verbatim}
    1443 
    1444 And we close down the parallel job by issuing the command
    1445 
    1446 \begin{verbatim}
    1447 finalize()
    1448 \end{verbatim}
    1449 
    1450 
    1451 %=================================
    1452 \chapter{Checkpointing}
    1453 
    1454 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.
    1455 
    1456 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.
    1457 
    1458 \section{The Code}
    1459 
    1460 Here is the code for \file{runCheckpoint.py}:
    1461 
    1462 \verbatiminputunderscore{../../anuga_core/demos/checkpointing/runCheckpoint.py}
    1463 
    1464 \section{Structure of the Code}
    1465 
    1466 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.
    1467 
    1468 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.
    1469 
    1470 
    1471 \begin{verbatim}
    1472 try:   
    1473     from anuga import load_checkpoint_file
    1474        
    1475     domain = load_checkpoint_file(domain_name = domain_name, checkpoint_dir = checkpoint_dir)
    1476 
    1477 except:
    1478     # create the domain as usual
    1479 \end{verbatim}
    1480 
    1481 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'').
    1482 
    1483 When the domain is originally created, we need to setup checkpointing via the command:
    1484 \begin{verbatim}
    1485     if useCheckpointing:
    1486         domain.set_checkpointing(checkpoint_time = 5)
    1487 \end{verbatim}
    1488 
    1489 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.
    14901392
    14911393
     
    15021404will run the tests in verbose mode (produces output) and in parallel using 6 processors (if parallel version of anuga has been setup).
    15031405
    1504 A subset of these validation tests can be run using \code{run_auto_validation_tests.py} from the \module{anuga_core} directory.
     1406A subset of these validation tests can be run using \code{run_validations.py} from the \module{anuga_core} directory.
    15051407
    15061408
     
    15851487
    15861488
    1587 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    1588 
    1589 \chapter{\anuga Public Interface}
    1590 \label{ch:interface}
    1591 
    1592 This chapter gives an overview of the features of \anuga available
    1593 to the user at the public interface. These are grouped under the
    1594 following headings, which correspond to the outline of the examples
    1595 described in Chapter \ref{ch:getstarted}:
    1596 \begin{itemize}
    1597     \item Establishing the Mesh: Section \ref{sec:establishing the mesh}
    1598     \item Initialising the Domain: Section \ref{sec:initialising the domain}
    1599 %    \item Specifying the Quantities: Section \ref{sec:quantities}
    1600     \item Initial Conditions: Section \ref{sec:initial conditions}
    1601     \item Boundary Conditions: Section \ref{sec:boundary conditions}
    1602     \item Operators: Section \ref{sec:operators}
    1603     \item Evolution: Section \ref{sec:evolution}
    1604 \end{itemize}
    1605 
    1606 \section{Documentation}\index{Documentation}
    1607 
    1608 The listings here are intended merely to give the reader an idea of what
    1609 each feature is and how it can be used -- they do
    1610 not give full specifications; for these the reader
    1611 may consult the programmer's guide. The code for every function or class contains
    1612 a documentation string, or 'docstring', that specifies the precise
    1613 syntax for its use. This appears immediately after the line
    1614 introducing the code, between two sets of triple quotes.
    1615 
    1616 Python has a handy tool that lets you easily navigate this documentation,
    1617 called \code{pydoc}. In Linux, it runs as a server, which serves the
    1618 documentation up to your web browser:
    1619 
    1620 1. Open a terminal at your \code{anuga} folder
    1621 2. Start the python documentation server with \code{pydoc -p 6767}
    1622 3. Open a browser and type in \code{http://localhost:6767/}
    1623 
    1624 Now you have a real-time programmers' guide for \anuga, and an easy way to find
    1625 the functions you are interested in. Pydoctor and other Python doc generators
    1626 look nicer and have graphs, etc, but pydoc works straight out of the box.
    1627 
    1628 \section{Public vs Private Interface}\index{public vs private interface}
    1629 
    1630 To simplify the process of writing scripts, \anuga has a public API which packages
    1631 up the commonly-used functionality of \anuga in the one place. To use it, simply
    1632 import the anuga module like so:
    1633 
    1634 \begin{verbatim}
    1635 import anuga
    1636 \end{verbatim}
    1637 
    1638 You can now use the public API like so. Note the \code{anuga.} prefix:
    1639 
    1640 \begin{verbatim}
    1641 anuga.sww2dem('in.sww', 'out.asc')
    1642 \end{verbatim}
    1643 
    1644 If you wish to delve "under the hood" and modify the way \anuga runs at a more
    1645 advanced level, you need to specify the full location of the module like so:
    1646 
    1647 \begin{verbatim}
    1648 from anuga.fit_interpolate.interpolate import Interpolation_interface
    1649 \end{verbatim}
    1650 
    1651 All modules are in the folder \file{anuga} or one of its subfolders, and the
    1652 location of each module is described relative to \file{anuga}. Rather
    1653 than using pathnames, whose syntax depends on the operating system,
    1654 we use the format adopted for importing the function or class for
    1655 use in Python code. For example, suppose we wish to specify that the
    1656 function \function{create\_mesh\_from\_regions} is in a module called
    1657 \module{mesh\_interface} in a subfolder of \module{anuga} called
    1658 \code{pmesh}. In Linux or Unix syntax, the pathname of the file
    1659 containing the function, relative to \file{anuga}, would be:
    1660 
    1661 \begin{verbatim}
    1662 pmesh/mesh_interface.py
    1663 \end{verbatim}
    1664 
    1665 \label{sec:mesh interface}
    1666 while in Windows syntax it would be:
    1667 
    1668 \begin{verbatim}
    1669 pmesh\mesh_interface.py
    1670 \end{verbatim}
    1671 
    1672 Rather than using either of these forms, in this chapter we specify
    1673 the location simply as \code{anuga.pmesh.mesh_interface}, in keeping with
    1674 the usage in the Python statement for importing the function,
    1675 namely:
    1676 
    1677 \begin{verbatim}
    1678 from anuga.pmesh.mesh_interface import create_mesh_from_regions
    1679 \end{verbatim}
    1680 
    1681 
    1682 The following parameters are common to many functions and classes
    1683 and are omitted from the descriptions given below:
    1684 
    1685 %\begin{tabular}{ll}
    1686 \begin{tabular}{p{2.0cm} p{14.0cm}}
    1687   \emph{use\_cache} & Specifies whether caching is to be used for improved performance. See Section \ref{sec:caching} for details\\
    1688                     & on the underlying caching functionality\\                     
    1689   \emph{verbose}    & If \code{True}, provides detailed terminal output to the user\\
    1690 \end{tabular}
    1691 
    1692 %=======================================================================
    1693 \section{Mesh Generation}\index{Mesh!generation}
     1489%=================================
     1490\chapter{Mesh Generation}\index{Mesh!generation}
    16941491\label{sec:establishing the mesh}
    16951492Before discussing the part of the interface relating to mesh
     
    17071504  \label{fig:simplemesh}
    17081505\end{figure}
    1709 
    1710 \clearpage
    17111506
    17121507The variables \code{points}, \code{triangles} and \code{boundary}
     
    17741569
    17751570
    1776 %===========================================
    1777 \subsection{High Level Mesh Generation Functions}
     1571%------------------------------------------------------------------
     1572\section{High Level Mesh Generation Functions}
    17781573% \refmodindex[pmesh.meshinterface]{pmesh.mesh\_interface}
    17791574\label{sec:meshgeneration}
     
    19121707\end{funcdesc}
    19131708
    1914 \subsection{Advanced mesh generation}
     1709%-----------------------------------------------------
     1710\section{Advanced mesh generation}
    19151711
    19161712For more control over the creation of the mesh outline, use the
     
    19381734\end{classdesc}
    19391735
    1940 \subsubsection{Key Methods of Class Mesh}
     1736\subsection{Key Methods of Class Mesh}
    19411737
    19421738\begin{methoddesc}{\emph{<mesh>}.add_hole}{x, y, geo_reference=None}
     
    21031899\end{methoddesc}
    21041900
    2105 %==============================================
    2106 \section{Initialising the Domain}\index{Initialising the Domain}
     1901
     1902%=================================
     1903\chapter{The Domain}\index{Initialising the Domain}
    21071904\label{sec:initialising the domain}
    21081905
     
    21261923Module: \refmodule{abstract_2d_finite_volumes.domain}
    21271924
    2128 This class is used to create an instance of a structure used to
    2129 store and manipulate data associated with a mesh. The mesh is
     1925This class isused to
     1926store and manipulate data associated with a mesh and to setup the particulars for
     1927the evolution of the specific problem. The mesh is
    21301928specified either by assigning the name of a mesh file to
    21311929\code{source} or by specifying the points, triangle and boundary of the
     
    21331931\end{classdesc}
    21341932
    2135 \subsection{Key Methods of Domain}
     1933%----------------------------------------------
     1934\section{Key Methods of Domain}
    21361935
    21371936\begin{methoddesc}{\emph{<domain>}.set_name}{name}
     
    23652164\end{methoddesc}
    23662165
    2367 
    2368 \section{Initial Conditions}\index{Initial Conditions}
     2166%-----------------------------------------
     2167\section{Evolution}\index{evolution}
     2168\label{sec:evolution}
     2169
     2170\begin{methoddesc}{\emph{<domain>}.evolve}{yieldstep=None,
     2171                                           finaltime=None,
     2172                                           duration=None,
     2173                                           skip_initial_step=False}
     2174Module: \module{abstract_2d_finite_volumes.domain}
     2175
     2176This method is invoked once all the
     2177preliminaries have been completed, and causes the model to progress
     2178through successive steps in its evolution, storing results and
     2179outputting statistics whenever a user-specified period
     2180\code{yieldstep} is completed.  Generally during this period the
     2181model will evolve through several steps internally
     2182as the method forces the water speed to be calculated
     2183on successive new cells.
     2184
     2185\code{yieldstep} is the interval in seconds between yields where results are
     2186stored, statistics written and the domain is inspected or possibly modified.
     2187If omitted an internal predefined \code{yieldstep} is used.  Internally, smaller
     2188timesteps may be taken.
     2189
     2190\code{duration} is the duration of the simulation in seconds.
     2191
     2192\code{finaltime} is the time in seconds where simulation should end. This is currently
     2193relative time, so it's the same as \code{duration}.  If both \code{duration} and
     2194\code{finaltime} are given an exception is thrown.
     2195
     2196\code{skip_initial_step} is a boolean flag that decides whether the first
     2197yield step is skipped or not. This is useful for example to avoid
     2198duplicate steps when multiple evolve processes are dove tailed.
     2199
     2200The code specified by the user in the block following the evolve statement is
     2201only executed once every \code{yieldstep} even though
     2202\anuga typically will take many more internal steps behind the scenes.
     2203
     2204You can include \method{evolve} in a statement of the type:
     2205
     2206\begin{verbatim}
     2207for t in domain.evolve(yieldstep, finaltime):
     2208    <Do something with domain and t>
     2209\end{verbatim}
     2210\end{methoddesc}
     2211
     2212%-----------------------------------------
     2213\section{Diagnostics}
     2214\label{sec:diagnostics}
     2215
     2216\begin{methoddesc}{\emph{<domain>}.statistics}{}
     2217Module: \module{abstract\_2d\_finite\_volumes.domain}
     2218
     2219Outputs statistics about the mesh within the \code{Domain}.
     2220\end{methoddesc}
     2221
     2222\begin{methoddesc}{\emph{<domain>}.timestepping_statistics}{track_speeds=False, triangle_id=None}
     2223Module: \module{abstract_2d_finite_volumes.domain}
     2224
     2225Returns a string of the following type for each timestep:\\
     2226\code{Time = 0.9000, delta t in [0.00598964, 0.01177388], steps=12 (12)}
     2227
     2228Here the numbers in \code{steps=12 (12)} indicate the number of steps taken and
     2229the number of first-order steps, respectively.
     2230
     2231The optional keyword argument \code{track_speeds} will
     2232generate a histogram of speeds generated by each triangle if set to \code{True}. The
     2233speeds relate to the size of the timesteps used by \anuga and
     2234this diagnostics may help pinpoint problem areas where excessive speeds
     2235are generated.
     2236
     2237The optional keyword argument \code{triangle_id} can be used to specify a particular
     2238triangle rather than the one with the largest speed.
     2239\end{methoddesc}
     2240
     2241\begin{methoddesc}{\emph{<domain>}.boundary_statistics}{quantities=None,
     2242                                                      tags=None}
     2243Module: \module{abstract_2d_finite_volumes.domain}
     2244
     2245Generates output about boundary forcing at each timestep.
     2246
     2247\code{quantities} names the quantities to be reported -- may be \code{None},
     2248a string or a list of strings.
     2249
     2250\code{tags} names the tags to be reported -- may be either None, a string or a list of strings.
     2251
     2252When \code{quantities = 'stage'} and \code{tags = ['top', 'bottom']}
     2253will return a string like:
     2254
     2255\begin{verbatim}
     2256Boundary values at time 0.5000:
     2257    top:
     2258        stage in [ -0.25821218,  -0.02499998]
     2259    bottom:
     2260        stage in [ -0.27098821,  -0.02499974]
     2261\end{verbatim}
     2262\end{methoddesc}
     2263
     2264\begin{methoddesc}{Q = \emph{<domain>}.get_quantity}{name,
     2265                                               location='vertices',
     2266                                               indices=None}
     2267Module: \module{abstract_2d_finite_volumes.domain}
     2268
     2269This function returns a Quantity object Q.
     2270Access to its values should be done through \code{Q.get_values()} documented on Page \pageref{pg:get values}.
     2271
     2272\code{name} is the name of the quantity to retrieve.
     2273
     2274\code{location} is
     2275
     2276\code{indices} is
     2277\end{methoddesc}
     2278
     2279\begin{methoddesc}{\emph{<domain>}.set_quantities_to_be_monitored}{quantity,
     2280                                             polygon=None,
     2281                                             time_interval=None}
     2282Module: \module{abstract\_2d\_finite\_volumes.domain}
     2283
     2284Selects quantities and derived quantities for which extrema attained at internal timesteps
     2285will be collected.
     2286
     2287\code{quantity} specifies the quantity or quantities to be monitored and must be either:
     2288\begin{itemize}
     2289  \item the name of a quantity or derived quantity such as 'stage-elevation',
     2290  \item a list of quantity names, or
     2291  \item \code{None}.
     2292\end{itemize}
     2293
     2294\code{polygon} can be used to monitor only triangles inside the polygon. Otherwise
     2295all triangles will be included.
     2296
     2297\code{time_interval} will restrict monitoring to time steps in the interval. Otherwise
     2298all timesteps will be included.
     2299
     2300Information can be tracked in the evolve loop by printing \code{quantity_statistics} and
     2301collected data will be stored in the SWW file.
     2302\end{methoddesc}
     2303
     2304\begin{methoddesc}{\emph{<domain>}.quantity_statistics}{precision='\%.4f'}
     2305Module: \module{abstract_2d_finite_volumes.domain}
     2306
     2307Reports on extrema attained by selected quantities.
     2308
     2309Returns a string of the following type for each timestep:
     2310
     2311\begin{verbatim}
     2312Monitored quantities at time 1.0000:
     2313  stage-elevation:
     2314    values since time = 0.00 in [0.00000000, 0.30000000]
     2315    minimum attained at time = 0.00000000, location = (0.16666667, 0.33333333)
     2316    maximum attained at time = 0.00000000, location = (0.83333333, 0.16666667)
     2317  ymomentum:
     2318    values since time = 0.00 in [0.00000000, 0.06241221]
     2319    minimum attained at time = 0.00000000, location = (0.33333333, 0.16666667)
     2320    maximum attained at time = 0.22472667, location = (0.83333333, 0.66666667)
     2321  xmomentum:
     2322    values since time = 0.00 in [-0.06062178, 0.47886313]
     2323    minimum attained at time = 0.00000000, location = (0.16666667, 0.33333333)
     2324    maximum attained at time = 0.35103646, location = (0.83333333, 0.16666667)
     2325\end{verbatim}
     2326
     2327The quantities (and derived quantities) listed here must be selected at model
     2328initialisation time using the method \code{domain.set_quantities_to_be_monitored()}.
     2329
     2330The optional keyword argument \code{precision='\%.4f'} will
     2331determine the precision used for floating point values in the output.
     2332This diagnostics helps track extrema attained by the selected quantities
     2333at every internal timestep.
     2334
     2335These values are also stored in the SWW file for post-processing.
     2336\end{methoddesc}
     2337
     2338\begin{methoddesc}{\emph{<quantity>}.get_values}{interpolation_points=None,
     2339                   location='vertices',
     2340                   indices=None,
     2341                   use_cache=False,
     2342                   verbose=False}
     2343\label{pg:get values}
     2344Module: \module{abstract_2d_finite_volumes.quantity}
     2345
     2346Extract values for quantity as a numeric array.
     2347
     2348\code{interpolation_points} is a list of (x, y) coordinates where the value is
     2349sought (using interpolation). If \code{interpolation_points} is given, values
     2350for \code{location} and \code{indices} are ignored.
     2351Assume either an absolute UTM coordinates or geospatial data object.
     2352   
     2353\code{location} specifies where values are to be stored.
     2354Permissible options are 'vertices', 'edges', 'centroids' or 'unique vertices'.
     2355
     2356The returned values will have the leading dimension equal to length of the \code{indices} list or
     2357N (all values) if \code{indices} is \code{None}.
     2358
     2359If \code{location} is 'centroids' the dimension of returned
     2360values will be a list or a numeric array of length N, N being
     2361the number of elements.
     2362     
     2363If \code{location} is 'vertices' or 'edges' the dimension of
     2364returned values will be of dimension \code{Nx3}.
     2365
     2366If \code{location} is 'unique vertices' the average value at
     2367each vertex will be returned and the dimension of returned values
     2368will be a 1d array of length "number of vertices"
     2369     
     2370\code{indices} is the set of element ids that the operation applies to.
     2371
     2372The values will be stored in elements following their internal ordering.
     2373\end{methoddesc}
     2374
     2375\begin{methoddesc}{\emph{<quantity>}.set_values}{numeric=None,
     2376                               quantity=None,
     2377                               function=None,
     2378                               geospatial_data=None,
     2379                               filename=None,
     2380                               attribute_name=None,
     2381                               alpha=None,
     2382                               location='vertices',
     2383                               polygon=None,
     2384                               indices=None,
     2385                               smooth=False,
     2386                               verbose=False,
     2387                               use_cache=False}
     2388Module: \module{abstract_2d_finite_volumes.quantity}
     2389
     2390Assign values to a quantity object.
     2391
     2392This method works the same way as \code{set_quantity()} except that it doesn't take
     2393a quantity name as the first argument since it is applied directly to the quantity.
     2394Use \code{set_values} is used to assign
     2395values to a new quantity that has been created but which is
     2396not part of the domain's predefined quantities.
     2397
     2398\code{location} is ??
     2399
     2400\code{indices} is ??
     2401
     2402The method \code{set_values()} is always called by \code{set_quantity()}
     2403behind the scenes.
     2404\end{methoddesc}
     2405
     2406\begin{methoddesc}{\emph{<quantity>}.get_integral}{}
     2407Module: \module{abstract_2d_finite_volumes.quantity}
     2408
     2409Return the computed integral over the entire domain for the quantity.
     2410\end{methoddesc}
     2411
     2412\begin{methoddesc}{\emph{<quantity>}.get_maximum_value}{indices=None}
     2413Module: \module{abstract_2d_finite_volumes.quantity}
     2414
     2415Return the maximum value of a quantity (on centroids).
     2416
     2417\code{indices} is the optional set of element \code{id}s that
     2418the operation applies to.
     2419
     2420We do not seek the maximum at vertices as each vertex can
     2421have multiple values -- one for each triangle sharing it.
     2422\end{methoddesc}
     2423
     2424\begin{methoddesc}{\emph{<quantity>}.get_maximum_location}{indices=None}
     2425Module: \module{abstract_2d_finite_volumes.quantity}
     2426
     2427Return the location of the maximum value of a quantity (on centroids).
     2428
     2429\code{indices} is the optional set of element \code{id}s that
     2430the operation applies to.
     2431
     2432We do not seek the maximum at vertices as each vertex can
     2433have multiple values -- one for each triangle sharing it.
     2434
     2435If there are multiple cells with the same maximum value, the
     2436first cell encountered in the triangle array is returned.
     2437\end{methoddesc}
     2438
     2439\begin{methoddesc}{\emph{<domain>}.get_wet_elements}{indices=None}
     2440Module: \module{shallow_water.shallow_water_domain}
     2441
     2442Returns the indices for elements where h $>$ minimum_allowed_height
     2443
     2444\code{indices} is the optional set of element \code{id}s that
     2445the operation applies to.
     2446\end{methoddesc}
     2447
     2448\begin{methoddesc}{\emph{<domain>}.get_maximum_inundation_elevation}{indices=None}
     2449Module: \module{shallow_water.shallow_water_domain}
     2450
     2451Return highest elevation where h $>$ 0.
     2452
     2453\code{indices} is the optional set of element \code{id}s that
     2454the operation applies to.
     2455
     2456Example to find maximum runup elevation:
     2457\begin{verbatim}
     2458z = domain.get_maximum_inundation_elevation()
     2459\end{verbatim}
     2460\end{methoddesc}
     2461
     2462\begin{methoddesc}{\emph{<domain>}.get_maximum_inundation_location}{indices=None}
     2463Module: \module{shallow_water.shallow_water_domain}
     2464
     2465Return location (x,y) of highest elevation where h $>$ 0.
     2466
     2467\code{indices} is the optional set of element \code{id}s that
     2468the operation applies to.
     2469
     2470Example to find maximum runup location:
     2471\begin{verbatim}
     2472x, y = domain.get_maximum_inundation_location()
     2473\end{verbatim}
     2474\end{methoddesc}
     2475
     2476
     2477
     2478
     2479
     2480%=================================
     2481\chapter{Initial Conditions}\index{Initial Conditions}
    23692482\label{sec:initial conditions}
     2483
    23702484In standard usage of partial differential equations, initial conditions
    23712485refers to the values associated to the system variables (the conserved
     
    25042618\end{methoddesc}
    25052619
    2506 \begin{funcdesc}{slump_tsunami}{length, depth, slope, width=None, thickness=None,
     2620\begin{funcdesc}{\emph{<callable_object>} =slump_tsunami}{length, depth, slope, width=None, thickness=None,
    25072621                                radius=None, dphi=0.48, x0=0.0, y0=0.0, alpha=0.0,
    25082622                                gravity=9.8, gamma=1.85,
     
    26312745
    26322746
    2633 \section{Boundary Conditions}\index{boundary conditions}
     2747
     2748%=================================
     2749
     2750\chapter{Boundary Conditions}\index{boundary conditions}
    26342751\label{sec:boundary conditions}
    26352752
     
    26592776\end{methoddesc}
    26602777
    2661 \subsection{Predefined boundary conditions}
     2778%----------------------------------------------------------
     2779\section{Predefined boundary conditions}
    26622780
    26632781\begin{classdesc}{Reflective_boundary}{domain=None}
     
    28062924\end{classdesc}
    28072925
    2808 \subsection{User-defined boundary conditions}
     2926\section{User-defined boundary conditions}
    28092927\label{sec:roll your own}
    28102928
     
    28212939Please refer to the source code for the existing boundary conditions
    28222940for examples of how to implement boundary conditions.
     2941
     2942
     2943%=================================
     2944\chapter{Operators and Structures}\index{operators}\index{structures}
     2945\label{sec:operators}
     2946
     2947A way to effect the evolution is via fractional step \code{operators}.  The idea is that at each inner timestep we can use
     2948an \code{operator} to change the centroid values of quantities such as the stage and the xmomentum and the ymomentum. T
     2949
     2950If 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.
     2951The moral, play with \code{elevation} at your peril if using the old algorithms.
     2952
     2953The newer algorithms (check that you are using discontinuous elevation with the domain member function \code{<domain>.get_using_discontinuous_elevation()}
     2954
     2955Currently predefined operators:
     2956
     2957
     2958\begin{classdesc}{Set_elevation_operator}{domain,
     2959                 elevation=None,
     2960                 indices=None,
     2961                 polygon=None,
     2962                 center=None,
     2963                 radius=None,
     2964                 description = None,
     2965                 label = None,
     2966                 logging = False,
     2967                 verbose = False}
     2968Module: \module{anuga.operators.set_elevation_operator}
     2969
     2970Set the elevation in a region (careful to maintain continuity of elevation)
     2971
     2972
     2973\code{domain} is the domain being evolved.
     2974
     2975\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.
     2976
     2977\code{indices} is a list of triangle indices where the function \code{elevation} will be applied.
     2978
     2979\code{polygon} is a polygon. The function \code{elevation} will be applied to triangles with in the polygon.
     2980
     2981\code{center} is the center of a circle where the function \code{elevation} will be applied. (\code{radius} needs to be set).
     2982
     2983\code{radius} is the radius of a circle where the function \code{elevation} will be applied. (\code{center} needs to be set).
     2984
     2985\code{description} is a description of this operator.
     2986
     2987\code{label} is the label used when reporting on the operator.
     2988
     2989\code{logging} is the boolean flag to set logging to a file.
     2990
     2991\code{verbose} is the boolean flag to setup extended output from the operator.
     2992
     2993Example:
     2994
     2995\begin{verbatim}
     2996op0 = anuga.Set_elevation_operator(domain, elevation = lambda t : 0.01*t)
     2997\end{verbatim}
     2998would setup an operator which raises the elevation over the whole domain 1 cm per second.
     2999
     3000
     3001While
     3002\begin{verbatim}
     3003P = [[x0, y0], [x1, y0], [x1, y1], [x0, y1]]
     3004
     3005op0 = anuga.Set_elevation_operator(domain, \
     3006                       elevation = lambda t : 0.01*t, polygon=P)
     3007\end{verbatim}
     3008would setup an operator which raises the elevation over the polygon \code{P} 1 cm per second.
     3009
     3010\end{classdesc}
     3011
     3012%--------------------------------------------
     3013\section{Culvert operators}
     3014
     3015A 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.
     3016
     3017Only 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).
     3018
     3019Usage
     3020
     3021In \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.
     3022
     3023The 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).
     3024
     3025The 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.
     3026
     3027Culvert hydraulics information can be logged by switching logging on or off.
     3028
     3029To 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.
     3030\begin{verbatim}
     3031#------------------------------------- -----------------
     3032# Box Culvert modelled along two lines
     3033#------------------------------------- -----------------
     3034losses = {'inlet':0.5, 'outlet':1.0, 'bend':0.0, 'grate':0.0, 'pier': 0.0, 'other': 0.0}
     3035el0 = numpy.array([[297750.2,6181379.1], [297753.5,6181380.9]])
     3036el1 = numpy.array([[297731.8,6181416.4], [297735.1,6181418.3]]) 
     3037culvert = anuga.Boyd_box_operator(domain,
     3038    losses=losses,
     3039    width=1.5,
     3040    height=1.5,
     3041    end_points=[ep0, ep1],
     3042    use_momentum_jet=True, #False to switch it off
     3043    use_velocity_head=True, #False to switch it off
     3044    manning=0.013, #culvert manning's n
     3045    logging=True, #False to switch it off
     3046    label='culvert_log_file',
     3047    verbose=False) 
     3048\end{verbatim}
     3049
     3050\begin{verbatim}
     3051#------------------------------------- ------------------
     3052# Pipe Culvert modelled with a single point
     3053#------------------------------------- ------------------
     3054losses = {'inlet':0.5, 'outlet':1.0, 'bend':0.0, 'grate':0.0, 'pier': 0.0, 'other': 0.0}
     3055ep0 = numpy.array([296653.0,6180014.9])
     3056ep1 = numpy.array([296642.5,6180036.3])     
     3057culvert = anuga.Boyd_pipe_operator(domain,
     3058    losses=losses,
     3059    diameter=1.5,
     3060    end_points=[ep0, ep1],
     3061    use_momentum_jet=True, #False to switch it off
     3062    use_velocity_head=True, #False to switch it off
     3063    manning=0.013, #culvert manning's n
     3064    logging=True, #False to switch it off
     3065    label='culvert_log_file',
     3066    verbose=False)
     3067\end{verbatim}
     3068
     3069\section{User developed operators}
     3070
     3071Suppose 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.
     3072
     3073Here is some code to do this:
     3074\begin{verbatim}
     3075from anuga import Operator
     3076
     3077class My_operator(Operator):
     3078    """
     3079   An operator which adds water to the domain at a given rate (m/sec)
     3080    """
     3081
     3082    def __init__(self,
     3083                 domain,
     3084                 rate=0.0,
     3085                 description = None,
     3086                 label = None,
     3087                 logging = False,
     3088                 verbose = False):
     3089          """Initialize the user defined operator
     3090          """
     3091
     3092        Operator.__init__(self, domain, description, label, logging, verbose)
     3093       
     3094        self.rate = rate
     3095
     3096
     3097    def __call__(self):
     3098        """
     3099        Apply rate to all triangles
     3100        """
     3101
     3102        timestep = self.domain.get_timestep()
     3103
     3104        self.stage_c[:] = self.stage_c[:]  + self.rate*timestep
     3105\end{verbatim}
     3106
     3107And then in your script you would asosciated \code{My_operator} with the \code{domain} with a call
     3108\begin{verbatim}
     3109My_operator(domain,rate=1.0)
     3110\end{verbatim}
     3111
     3112
     3113%=================================
     3114\chapter{Querying SWW model output files}
     3115
     3116After a model has been run, it is often useful to extract various information from the SWW
     3117output file (see Section \ref{sec:sww format}). This is typically more convenient than using the
     3118diagnostics described in Section \ref{sec:diagnostics} which rely on the model running -- something
     3119that can be very time consuming. The SWW files are easy and quick to read and offer information
     3120about the model results such as runup heights, time histories of selected quantities,
     3121flow through cross sections and much more.
     3122
     3123\begin{funcdesc}{elevation = get_maximum_inundation_elevation}{filename,
     3124                                     polygon=None,
     3125                                     time_interval=None,
     3126                                     verbose=False}
     3127Module: \module{shallow_water.data_manager}
     3128
     3129Return the highest elevation where depth is positive ($h > 0$).
     3130
     3131\code{filename} is the path to a NetCDF SWW file containing \anuga model output.
     3132
     3133\code{polygon} restricts the query to the points that lie within the polygon.
     3134
     3135\code {time_interval} restricts the query to within the time interval.
     3136
     3137Example to find maximum runup elevation:
     3138
     3139\begin{verbatim}
     3140max_runup = get_maximum_inundation_elevation(filename)
     3141\end{verbatim}
     3142
     3143If no inundation is found (within the \code{polygon} and \code{time_interval}, if specified)
     3144the return value is \code{None}. This indicates "No Runup" or "Everything is dry".
     3145\end{funcdesc}
     3146
     3147\begin{funcdesc}{(x, y) = get_maximum_inundation_location}{filename,
     3148                                    polygon=None,
     3149                                    time_interval=None,
     3150                                    verbose=False}
     3151Module: \module{shallow_water.data_manager}
     3152
     3153Return location (x,y) of the highest elevation where depth is positive ($h > 0$).
     3154
     3155\code{filename} is the path to a NetCDF SWW file containing \anuga model output.
     3156
     3157\code{polygon} restricts the query to the points that lie within the polygon.
     3158
     3159\code {time_interval} restricts the query to within the time interval.
     3160
     3161Example to find maximum runup location:
     3162
     3163\begin{verbatim}
     3164max_runup_location = get_maximum_inundation_location(filename)
     3165\end{verbatim}
     3166
     3167If no inundation is found (within the \code{polygon} and \code{time_interval}, if specified)
     3168the return value is \code{None}. This indicates "No Runup" or "Everything is dry".
     3169\end{funcdesc}
     3170
     3171\begin{funcdesc}{sww2timeseries}{swwfiles,
     3172                                 gauge_filename,
     3173                                 production_dirs,
     3174                                 report=None,
     3175                                 reportname=None,
     3176                                 plot_quantity=None,
     3177                                 generate_fig=False,
     3178                                 surface=None,
     3179                                 time_min=None,
     3180                                 time_max=None,
     3181                                 output_centroids=False,
     3182                                 time_thinning=1,
     3183                                 time_unit=None,
     3184                                 title_on=None,
     3185                                 use_cache=False,
     3186                                 verbose=False}
     3187Module: \module{abstract_2d_finite_volumes.util}
     3188
     3189Read a set of SWW files and plot the time series for the prescribed quantities
     3190at defined gauge locations and prescribed time range.
     3191
     3192\code{swwfiles} is a dictionary of SWW files with keys of \code{label_id}.
     3193
     3194\code{gauge_filename} is the path to a file containing gauge data.
     3195
     3196THIS NEEDS MORE WORK.  WORK ON FUNCTION __DOC__ STRING, IF NOTHING ELSE!
     3197\end{funcdesc}
     3198
     3199\begin{funcdesc}{(time, Q) = get_flow_through_cross_section}{filename, polyline, verbose=False}
     3200Module: \module{shallow_water.data_manager}
     3201
     3202Obtain flow ($m^3/s$) perpendicular to specified cross section.
     3203
     3204\code{filename} is the path to the SWW file.
     3205
     3206\code{output_centroids} set to true to sample at the centre of the triangle containing the point.
     3207This may be useful for debugging. (new in 1.2)
     3208
     3209\code{polyline} is the representation of the desired cross section -- it may contain
     3210multiple sections allowing for complex shapes. Assumes absolute UTM coordinates.
     3211
     3212Returns a tuple \code{time, Q} where:
     3213
     3214\code{time} is all the stored times in the SWW file.
     3215
     3216\code{Q} is a hydrograph of total flow across the given segments for all stored times.
     3217
     3218The normal flow is computed for each triangle intersected by the \code{polyline} and
     3219added up.  If multiple segments at different angles are specified the normal flows
     3220may partially cancel each other.
     3221
     3222Example to find flow through cross section:
     3223
     3224\begin{verbatim}
     3225cross_section = [[x, 0], [x, width]]
     3226time, Q = get_flow_through_cross_section(filename, cross_section)
     3227\end{verbatim}
     3228\end{funcdesc}
     3229
     3230
     3231%=================================
     3232\chapter{Parallel Simulation}
     3233
     3234The 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.
     3235
     3236One way to alleviate this is to run your simulation in parallel and use more processors to spread the computational cost.
     3237
     3238
     3239\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}.
     3240
     3241Such jobs are run using the command
     3242
     3243\begin{verbatim}
     3244mpirun -np n python runscript.py
     3245\end{verbatim}
     3246where \code{n} is the total number of processors being used for the parallel run.
     3247
     3248Essentially 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.
     3249
     3250For 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.
     3251We would expect 70\% scalability, and so expect a speedup of $250 \times 0.7 = 175$.
     3252
     3253\section{The Code}
     3254
     3255Here is the code for \file{runParallelCairns.py} which is found in the \file{examples/cairns} directory:
     3256
     3257\verbatiminputunderscore{../../anuga_core/examples/cairns/runParallelCairns.py}
     3258
     3259
     3260\section{Structure of the Code}
     3261
     3262The 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.
     3263
     3264
     3265We first import a few procedures  need for the parallel code.
     3266\begin{verbatim}
     3267from anuga import distribute, myid, numprocs, finalize, barrier
     3268\end{verbatim}
     3269
     3270\code{myid} returns the id of the current processor running the code.
     3271
     3272\code{numprocs} returns the total number of processors involved in this parallel jobs (the \code{n} in the original \code{mpirun} command.
     3273
     3274\code{finalize} is called at the end of a script to close down the parallel job.
     3275
     3276\code{distribute} is used to partition and setup the parallel domains.
     3277
     3278\code{barrier} is a command for processors to wait as all the other processor to catchup to this point.
     3279
     3280
     3281
     3282
     3283The creation of the \code{domain} is only done on processor 0. Hence we have the structure:
     3284
     3285\begin{verbatim}
     3286#-------------------------------------
     3287# Do the domain creation on processor 0
     3288#-------------------------------------
     3289if myid == 0:
     3290    ....
     3291    domain = ...
     3292       
     3293else:
     3294    domain = None
     3295\end{verbatim}
     3296
     3297We 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).
     3298
     3299Once we have our code{domain} setup we partition it and send the partitions to each of the other processors, via the command
     3300
     3301\begin{verbatim}
     3302#-------------------------------------
     3303# Now produce parallel domain
     3304#-------------------------------------
     3305domain = distribute(domain)
     3306\end{verbatim}
     3307
     3308This 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.
     3309
     3310It is important to apply the boundary conditions after the \code{distribute}
     3311
     3312
     3313
     3314As 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.
     3315
     3316We have a procedure to merge these individual sww files via the command
     3317
     3318\begin{verbatim}
     3319domain.sww_merge()
     3320\end{verbatim}
     3321
     3322And we close down the parallel job by issuing the command
     3323
     3324\begin{verbatim}
     3325finalize()
     3326\end{verbatim}
     3327
     3328
     3329%=================================
     3330\chapter{Checkpointing}
     3331
     3332When 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.
     3333
     3334The 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.
     3335
     3336\section{The Code}
     3337
     3338Here is the code for \file{runCheckpoint.py}:
     3339
     3340\verbatiminputunderscore{../../anuga_core/examples/checkpointing/runCheckpoint.py}
     3341
     3342\section{Structure of the Code}
     3343
     3344As usual the code needs to import the required modules, and setup parameters and in this case procedures for setting up the stage and elevation.
     3345
     3346Then 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.
     3347
     3348
     3349\begin{verbatim}
     3350try:   
     3351    from anuga import load_checkpoint_file
     3352       
     3353    domain = load_checkpoint_file(domain_name = domain_name,
     3354                                  checkpoint_dir = checkpoint_dir)
     3355except:
     3356    # create the domain as usual
     3357\end{verbatim}
     3358
     3359The \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'').
     3360
     3361When the domain is originally created, we need to setup checkpointing via the command:
     3362\begin{verbatim}
     3363    if useCheckpointing:
     3364        domain.set_checkpointing(checkpoint_time = 5)
     3365\end{verbatim}
     3366
     3367Here 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.
     3368
     3369
     3370
     3371
     3372%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     3373
     3374\chapter{\anuga Public Interface}
     3375\label{ch:interface}
     3376
     3377This chapter gives an overview of the features of \anuga available
     3378to the user at the public interface. These are grouped under the
     3379following headings, which correspond to the outline of the examples
     3380described in Chapter \ref{ch:getstarted}:
     3381\begin{itemize}
     3382    \item Establishing the Mesh: Section \ref{sec:establishing the mesh}
     3383    \item Initialising the Domain: Section \ref{sec:initialising the domain}
     3384%    \item Specifying the Quantities: Section \ref{sec:quantities}
     3385    \item Initial Conditions: Section \ref{sec:initial conditions}
     3386    \item Boundary Conditions: Section \ref{sec:boundary conditions}
     3387    \item Operators: Section \ref{sec:operators}
     3388    \item Evolution: Section \ref{sec:evolution}
     3389\end{itemize}
     3390
     3391\section{Documentation}\index{Documentation}
     3392
     3393The listings here are intended merely to give the reader an idea of what
     3394each feature is and how it can be used -- they do
     3395not give full specifications; for these the reader
     3396may consult the programmer's guide. The code for every function or class contains
     3397a documentation string, or 'docstring', that specifies the precise
     3398syntax for its use. This appears immediately after the line
     3399introducing the code, between two sets of triple quotes.
     3400
     3401Python has a handy tool that lets you easily navigate this documentation,
     3402called \code{pydoc}. In Linux, it runs as a server, which serves the
     3403documentation up to your web browser:
     3404
     34051. Open a terminal at your \code{anuga} folder
     34062. Start the python documentation server with \code{pydoc -p 6767}
     34073. Open a browser and type in \code{http://localhost:6767/}
     3408
     3409Now you have a real-time programmers' guide for \anuga, and an easy way to find
     3410the functions you are interested in. Pydoctor and other Python doc generators
     3411look nicer and have graphs, etc, but pydoc works straight out of the box.
     3412
     3413\section{Public vs Private Interface}\index{public vs private interface}
     3414
     3415To simplify the process of writing scripts, \anuga has a public API which packages
     3416up the commonly-used functionality of \anuga in the one place. To use it, simply
     3417import the anuga module like so:
     3418
     3419\begin{verbatim}
     3420import anuga
     3421\end{verbatim}
     3422
     3423You can now use the public API like so. Note the \code{anuga.} prefix:
     3424
     3425\begin{verbatim}
     3426anuga.sww2dem('in.sww', 'out.asc')
     3427\end{verbatim}
     3428
     3429If you wish to delve "under the hood" and modify the way \anuga runs at a more
     3430advanced level, you need to specify the full location of the module like so:
     3431
     3432\begin{verbatim}
     3433from anuga.fit_interpolate.interpolate import Interpolation_interface
     3434\end{verbatim}
     3435
     3436All modules are in the folder \file{anuga} or one of its subfolders, and the
     3437location of each module is described relative to \file{anuga}. Rather
     3438than using pathnames, whose syntax depends on the operating system,
     3439we use the format adopted for importing the function or class for
     3440use in Python code. For example, suppose we wish to specify that the
     3441function \function{create\_mesh\_from\_regions} is in a module called
     3442\module{mesh\_interface} in a subfolder of \module{anuga} called
     3443\code{pmesh}. In Linux or Unix syntax, the pathname of the file
     3444containing the function, relative to \file{anuga}, would be:
     3445
     3446\begin{verbatim}
     3447pmesh/mesh_interface.py
     3448\end{verbatim}
     3449
     3450\label{sec:mesh interface}
     3451while in Windows syntax it would be:
     3452
     3453\begin{verbatim}
     3454pmesh\mesh_interface.py
     3455\end{verbatim}
     3456
     3457Rather than using either of these forms, in this chapter we specify
     3458the location simply as \code{anuga.pmesh.mesh_interface}, in keeping with
     3459the usage in the Python statement for importing the function,
     3460namely:
     3461
     3462\begin{verbatim}
     3463from anuga.pmesh.mesh_interface import create_mesh_from_regions
     3464\end{verbatim}
     3465
     3466
     3467The following parameters are common to many functions and classes
     3468and are omitted from the descriptions given below:
     3469
     3470%\begin{tabular}{ll}
     3471\begin{tabular}{p{2.0cm} p{14.0cm}}
     3472  \emph{use\_cache} & Specifies whether caching is to be used for improved performance. See Section \ref{sec:caching} for details\\
     3473                    & on the underlying caching functionality\\                     
     3474  \emph{verbose}    & If \code{True}, provides detailed terminal output to the user\\
     3475\end{tabular}
     3476
     3477
     3478
    28233479
    28243480
     
    30713727
    30723728%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    3073 \section{Operators}\index{operators}
    3074 \label{sec:operators}
    3075 
    3076 A way to effect the evolution is via fractional step \code{operators}.  The idea is that at each inner timestep we can use
    3077 an \code{operator} to change the centroid values of quantities such as the stage and the xmomentum and the ymomentum. T
    3078 
    3079 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.
    3080 The moral, play with \code{elevation} at your peril if using the old algorithms.
    3081 
    3082 The newer algorithms (check that you are using discontinuous elevation with the domain member function \code{<domain>.get_using_discontinuous_elevation()}
    3083 
    3084 Currently predefined operators:
    3085 
    3086 
    3087 \begin{classdesc}{Set_elevation_operator}{domain,
    3088                  elevation=None,
    3089                  indices=None,
    3090                  polygon=None,
    3091                  center=None,
    3092                  radius=None,
    3093                  description = None,
    3094                  label = None,
    3095                  logging = False,
    3096                  verbose = False}
    3097 Module: \module{anuga.operators.set_elevation_operator}
    3098 
    3099 Set the elevation in a region (careful to maintain continuity of elevation)
    3100 
    3101 
    3102 \code{domain} is the domain being evolved.
    3103 
    3104 \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.
    3105 
    3106 \code{indices} is a list of triangle indices where the function \code{elevation} will be applied.
    3107 
    3108 \code{polygon} is a polygon. The function \code{elevation} will be applied to triangles with in the polygon.
    3109 
    3110 \code{center} is the center of a circle where the function \code{elevation} will be applied. (\code{radius} needs to be set).
    3111 
    3112 \code{radius} is the radius of a circle where the function \code{elevation} will be applied. (\code{center} needs to be set).
    3113 
    3114 \code{description} is a description of this operator.
    3115 
    3116 \code{label} is the label used when reporting on the operator.
    3117 
    3118 \code{logging} is the boolean flag to set logging to a file.
    3119 
    3120 \code{verbose} is the boolean flag to setup extended output from the operator.
    3121 
    3122 Example:
    3123 
    3124 \begin{verbatim}
    3125 op0 = anuga.Set_elevation_operator(domain, elevation = lambda t : 0.01*t)
    3126 \end{verbatim}
    3127 would setup an operator which raises the elevation over the whole domain 1 cm per second.
    3128 
    3129 
    3130 While
    3131 \begin{verbatim}
    3132 P = [[x0, y0], [x1, y0], [x1, y1], [x0, y1]]
    3133 
    3134 op0 = anuga.Set_elevation_operator(domain, \
    3135                        elevation = lambda t : 0.01*t, polygon=P)
    3136 \end{verbatim}
    3137 would setup an operator which raises the elevation over the polygon \code{P} 1 cm per second.
    3138 
    3139 \end{classdesc}
    3140 
    3141 %--------------------------------------------
    3142 \subsection{Culvert operators}
    3143 
    3144 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.
    3145 
    3146 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).
    3147 
    3148 Usage
    3149 
    3150 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.
    3151 
    3152 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).
    3153 
    3154 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.
    3155 
    3156 Culvert hydraulics information can be logged by switching logging on or off.
    3157 
    3158 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.
    3159 \begin{verbatim}
    3160 #------------------------------------- -----------------
    3161 # Box Culvert modelled along two lines
    3162 #------------------------------------- -----------------
    3163 losses = {'inlet':0.5, 'outlet':1.0, 'bend':0.0, 'grate':0.0, 'pier': 0.0, 'other': 0.0}
    3164 el0 = numpy.array([[297750.2,6181379.1], [297753.5,6181380.9]])
    3165 el1 = numpy.array([[297731.8,6181416.4], [297735.1,6181418.3]]) 
    3166 culvert = anuga.Boyd_box_operator(domain,
    3167     losses=losses,
    3168     width=1.5,
    3169     height=1.5,
    3170     end_points=[ep0, ep1],
    3171     use_momentum_jet=True, #False to switch it off
    3172     use_velocity_head=True, #False to switch it off
    3173     manning=0.013, #culvert manning's n
    3174     logging=True, #False to switch it off
    3175     label='culvert_log_file',
    3176     verbose=False) 
    3177 \end{verbatim}
    3178 
    3179 \begin{verbatim}
    3180 #------------------------------------- ------------------
    3181 # Pipe Culvert modelled with a single point
    3182 #------------------------------------- ------------------
    3183 losses = {'inlet':0.5, 'outlet':1.0, 'bend':0.0, 'grate':0.0, 'pier': 0.0, 'other': 0.0}
    3184 ep0 = numpy.array([296653.0,6180014.9])
    3185 ep1 = numpy.array([296642.5,6180036.3])     
    3186 culvert = anuga.Boyd_pipe_operator(domain,
    3187     losses=losses,
    3188     diameter=1.5,
    3189     end_points=[ep0, ep1],
    3190     use_momentum_jet=True, #False to switch it off
    3191     use_velocity_head=True, #False to switch it off
    3192     manning=0.013, #culvert manning's n
    3193     logging=True, #False to switch it off
    3194     label='culvert_log_file',
    3195     verbose=False)
    3196 \end{verbatim}
    3197 
    3198 \subsection{User developed operators}
    3199 
    3200 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.
    3201 
    3202 Here is some code to do this:
    3203 \begin{verbatim}
    3204 from anuga import Operator
    3205 
    3206 class My_operator(Operator):
    3207     """
    3208    An operator which adds water to the domain at a given rate (m/sec)
    3209     """
    3210 
    3211     def __init__(self,
    3212                  domain,
    3213                  rate=0.0,
    3214                  description = None,
    3215                  label = None,
    3216                  logging = False,
    3217                  verbose = False):
    3218           """Initialize the user defined operator
    3219           """
    3220 
    3221         Operator.__init__(self, domain, description, label, logging, verbose)
    3222        
    3223         self.rate = rate
    3224 
    3225 
    3226     def __call__(self):
    3227         """
    3228         Apply rate to all triangles
    3229         """
    3230 
    3231         timestep = self.domain.get_timestep()
    3232 
    3233         self.stage_c[:] = self.stage_c[:]  + self.rate*timestep
    3234 \end{verbatim}
    3235 
    3236 And then in your script you would asosciated \code{My_operator} with the \code{domain} with a call
    3237 \begin{verbatim}
    3238 My_operator(domain,rate=1.0)
    3239 \end{verbatim}
    3240 
    3241 
    3242 
    3243 
    3244 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    3245 \section{Evolution}\index{evolution}
    3246 \label{sec:evolution}
    3247 
    3248 \begin{methoddesc}{\emph{<domain>}.evolve}{yieldstep=None,
    3249                                            finaltime=None,
    3250                                            duration=None,
    3251                                            skip_initial_step=False}
    3252 Module: \module{abstract_2d_finite_volumes.domain}
    3253 
    3254 This method is invoked once all the
    3255 preliminaries have been completed, and causes the model to progress
    3256 through successive steps in its evolution, storing results and
    3257 outputting statistics whenever a user-specified period
    3258 \code{yieldstep} is completed.  Generally during this period the
    3259 model will evolve through several steps internally
    3260 as the method forces the water speed to be calculated
    3261 on successive new cells.
    3262 
    3263 \code{yieldstep} is the interval in seconds between yields where results are
    3264 stored, statistics written and the domain is inspected or possibly modified.
    3265 If omitted an internal predefined \code{yieldstep} is used.  Internally, smaller
    3266 timesteps may be taken.
    3267 
    3268 \code{duration} is the duration of the simulation in seconds.
    3269 
    3270 \code{finaltime} is the time in seconds where simulation should end. This is currently
    3271 relative time, so it's the same as \code{duration}.  If both \code{duration} and
    3272 \code{finaltime} are given an exception is thrown.
    3273 
    3274 \code{skip_initial_step} is a boolean flag that decides whether the first
    3275 yield step is skipped or not. This is useful for example to avoid
    3276 duplicate steps when multiple evolve processes are dove tailed.
    3277 
    3278 The code specified by the user in the block following the evolve statement is
    3279 only executed once every \code{yieldstep} even though
    3280 \anuga typically will take many more internal steps behind the scenes.
    3281 
    3282 You can include \method{evolve} in a statement of the type:
    3283 
    3284 \begin{verbatim}
    3285 for t in domain.evolve(yieldstep, finaltime):
    3286     <Do something with domain and t>
    3287 \end{verbatim}
    3288 \end{methoddesc}
    3289 
    3290 \subsection{Diagnostics}
    3291 \label{sec:diagnostics}
    3292 
    3293 \begin{methoddesc}{\emph{<domain>}.statistics}{}
    3294 Module: \module{abstract\_2d\_finite\_volumes.domain}
    3295 
    3296 Outputs statistics about the mesh within the \code{Domain}.
    3297 \end{methoddesc}
    3298 
    3299 \begin{methoddesc}{\emph{<domain>}.timestepping_statistics}{track_speeds=False, triangle_id=None}
    3300 Module: \module{abstract_2d_finite_volumes.domain}
    3301 
    3302 Returns a string of the following type for each timestep:\\
    3303 \code{Time = 0.9000, delta t in [0.00598964, 0.01177388], steps=12 (12)}
    3304 
    3305 Here the numbers in \code{steps=12 (12)} indicate the number of steps taken and
    3306 the number of first-order steps, respectively.
    3307 
    3308 The optional keyword argument \code{track_speeds} will
    3309 generate a histogram of speeds generated by each triangle if set to \code{True}. The
    3310 speeds relate to the size of the timesteps used by \anuga and
    3311 this diagnostics may help pinpoint problem areas where excessive speeds
    3312 are generated.
    3313 
    3314 The optional keyword argument \code{triangle_id} can be used to specify a particular
    3315 triangle rather than the one with the largest speed.
    3316 \end{methoddesc}
    3317 
    3318 \begin{methoddesc}{\emph{<domain>}.boundary_statistics}{quantities=None,
    3319                                                       tags=None}
    3320 Module: \module{abstract_2d_finite_volumes.domain}
    3321 
    3322 Generates output about boundary forcing at each timestep.
    3323 
    3324 \code{quantities} names the quantities to be reported -- may be \code{None},
    3325 a string or a list of strings.
    3326 
    3327 \code{tags} names the tags to be reported -- may be either None, a string or a list of strings.
    3328 
    3329 When \code{quantities = 'stage'} and \code{tags = ['top', 'bottom']}
    3330 will return a string like:
    3331 
    3332 \begin{verbatim}
    3333 Boundary values at time 0.5000:
    3334     top:
    3335         stage in [ -0.25821218,  -0.02499998]
    3336     bottom:
    3337         stage in [ -0.27098821,  -0.02499974]
    3338 \end{verbatim}
    3339 \end{methoddesc}
    3340 
    3341 \begin{methoddesc}{Q = \emph{<domain>}.get_quantity}{name,
    3342                                                location='vertices',
    3343                                                indices=None}
    3344 Module: \module{abstract_2d_finite_volumes.domain}
    3345 
    3346 This function returns a Quantity object Q.
    3347 Access to its values should be done through \code{Q.get_values()} documented on Page \pageref{pg:get values}.
    3348 
    3349 \code{name} is the name of the quantity to retrieve.
    3350 
    3351 \code{location} is
    3352 
    3353 \code{indices} is
    3354 \end{methoddesc}
    3355 
    3356 \begin{methoddesc}{\emph{<domain>}.set_quantities_to_be_monitored}{quantity,
    3357                                              polygon=None,
    3358                                              time_interval=None}
    3359 Module: \module{abstract\_2d\_finite\_volumes.domain}
    3360 
    3361 Selects quantities and derived quantities for which extrema attained at internal timesteps
    3362 will be collected.
    3363 
    3364 \code{quantity} specifies the quantity or quantities to be monitored and must be either:
    3365 \begin{itemize}
    3366   \item the name of a quantity or derived quantity such as 'stage-elevation',
    3367   \item a list of quantity names, or
    3368   \item \code{None}.
    3369 \end{itemize}
    3370 
    3371 \code{polygon} can be used to monitor only triangles inside the polygon. Otherwise
    3372 all triangles will be included.
    3373 
    3374 \code{time_interval} will restrict monitoring to time steps in the interval. Otherwise
    3375 all timesteps will be included.
    3376 
    3377 Information can be tracked in the evolve loop by printing \code{quantity_statistics} and
    3378 collected data will be stored in the SWW file.
    3379 \end{methoddesc}
    3380 
    3381 \begin{methoddesc}{\emph{<domain>}.quantity_statistics}{precision='\%.4f'}
    3382 Module: \module{abstract_2d_finite_volumes.domain}
    3383 
    3384 Reports on extrema attained by selected quantities.
    3385 
    3386 Returns a string of the following type for each timestep:
    3387 
    3388 \begin{verbatim}
    3389 Monitored quantities at time 1.0000:
    3390   stage-elevation:
    3391     values since time = 0.00 in [0.00000000, 0.30000000]
    3392     minimum attained at time = 0.00000000, location = (0.16666667, 0.33333333)
    3393     maximum attained at time = 0.00000000, location = (0.83333333, 0.16666667)
    3394   ymomentum:
    3395     values since time = 0.00 in [0.00000000, 0.06241221]
    3396     minimum attained at time = 0.00000000, location = (0.33333333, 0.16666667)
    3397     maximum attained at time = 0.22472667, location = (0.83333333, 0.66666667)
    3398   xmomentum:
    3399     values since time = 0.00 in [-0.06062178, 0.47886313]
    3400     minimum attained at time = 0.00000000, location = (0.16666667, 0.33333333)
    3401     maximum attained at time = 0.35103646, location = (0.83333333, 0.16666667)
    3402 \end{verbatim}
    3403 
    3404 The quantities (and derived quantities) listed here must be selected at model
    3405 initialisation time using the method \code{domain.set_quantities_to_be_monitored()}.
    3406 
    3407 The optional keyword argument \code{precision='\%.4f'} will
    3408 determine the precision used for floating point values in the output.
    3409 This diagnostics helps track extrema attained by the selected quantities
    3410 at every internal timestep.
    3411 
    3412 These values are also stored in the SWW file for post-processing.
    3413 \end{methoddesc}
    3414 
    3415 \begin{methoddesc}{\emph{<quantity>}.get_values}{interpolation_points=None,
    3416                    location='vertices',
    3417                    indices=None,
    3418                    use_cache=False,
    3419                    verbose=False}
    3420 \label{pg:get values}
    3421 Module: \module{abstract_2d_finite_volumes.quantity}
    3422 
    3423 Extract values for quantity as a numeric array.
    3424 
    3425 \code{interpolation_points} is a list of (x, y) coordinates where the value is
    3426 sought (using interpolation). If \code{interpolation_points} is given, values
    3427 for \code{location} and \code{indices} are ignored.
    3428 Assume either an absolute UTM coordinates or geospatial data object.
    3429    
    3430 \code{location} specifies where values are to be stored.
    3431 Permissible options are 'vertices', 'edges', 'centroids' or 'unique vertices'.
    3432 
    3433 The returned values will have the leading dimension equal to length of the \code{indices} list or
    3434 N (all values) if \code{indices} is \code{None}.
    3435 
    3436 If \code{location} is 'centroids' the dimension of returned
    3437 values will be a list or a numeric array of length N, N being
    3438 the number of elements.
    3439      
    3440 If \code{location} is 'vertices' or 'edges' the dimension of
    3441 returned values will be of dimension \code{Nx3}.
    3442 
    3443 If \code{location} is 'unique vertices' the average value at
    3444 each vertex will be returned and the dimension of returned values
    3445 will be a 1d array of length "number of vertices"
    3446      
    3447 \code{indices} is the set of element ids that the operation applies to.
    3448 
    3449 The values will be stored in elements following their internal ordering.
    3450 \end{methoddesc}
    3451 
    3452 \begin{methoddesc}{\emph{<quantity>}.set_values}{numeric=None,
    3453                                quantity=None,
    3454                                function=None,
    3455                                geospatial_data=None,
    3456                                filename=None,
    3457                                attribute_name=None,
    3458                                alpha=None,
    3459                                location='vertices',
    3460                                polygon=None,
    3461                                indices=None,
    3462                                smooth=False,
    3463                                verbose=False,
    3464                                use_cache=False}
    3465 Module: \module{abstract_2d_finite_volumes.quantity}
    3466 
    3467 Assign values to a quantity object.
    3468 
    3469 This method works the same way as \code{set_quantity()} except that it doesn't take
    3470 a quantity name as the first argument since it is applied directly to the quantity.
    3471 Use \code{set_values} is used to assign
    3472 values to a new quantity that has been created but which is
    3473 not part of the domain's predefined quantities.
    3474 
    3475 \code{location} is ??
    3476 
    3477 \code{indices} is ??
    3478 
    3479 The method \code{set_values()} is always called by \code{set_quantity()}
    3480 behind the scenes.
    3481 \end{methoddesc}
    3482 
    3483 \begin{methoddesc}{\emph{<quantity>}.get_integral}{}
    3484 Module: \module{abstract_2d_finite_volumes.quantity}
    3485 
    3486 Return the computed integral over the entire domain for the quantity.
    3487 \end{methoddesc}
    3488 
    3489 \begin{methoddesc}{\emph{<quantity>}.get_maximum_value}{indices=None}
    3490 Module: \module{abstract_2d_finite_volumes.quantity}
    3491 
    3492 Return the maximum value of a quantity (on centroids).
    3493 
    3494 \code{indices} is the optional set of element \code{id}s that
    3495 the operation applies to.
    3496 
    3497 We do not seek the maximum at vertices as each vertex can
    3498 have multiple values -- one for each triangle sharing it.
    3499 \end{methoddesc}
    3500 
    3501 \begin{methoddesc}{\emph{<quantity>}.get_maximum_location}{indices=None}
    3502 Module: \module{abstract_2d_finite_volumes.quantity}
    3503 
    3504 Return the location of the maximum value of a quantity (on centroids).
    3505 
    3506 \code{indices} is the optional set of element \code{id}s that
    3507 the operation applies to.
    3508 
    3509 We do not seek the maximum at vertices as each vertex can
    3510 have multiple values -- one for each triangle sharing it.
    3511 
    3512 If there are multiple cells with the same maximum value, the
    3513 first cell encountered in the triangle array is returned.
    3514 \end{methoddesc}
    3515 
    3516 \begin{methoddesc}{\emph{<domain>}.get_wet_elements}{indices=None}
    3517 Module: \module{shallow_water.shallow_water_domain}
    3518 
    3519 Returns the indices for elements where h $>$ minimum_allowed_height
    3520 
    3521 \code{indices} is the optional set of element \code{id}s that
    3522 the operation applies to.
    3523 \end{methoddesc}
    3524 
    3525 \begin{methoddesc}{\emph{<domain>}.get_maximum_inundation_elevation}{indices=None}
    3526 Module: \module{shallow_water.shallow_water_domain}
    3527 
    3528 Return highest elevation where h $>$ 0.
    3529 
    3530 \code{indices} is the optional set of element \code{id}s that
    3531 the operation applies to.
    3532 
    3533 Example to find maximum runup elevation:
    3534 \begin{verbatim}
    3535 z = domain.get_maximum_inundation_elevation()
    3536 \end{verbatim}
    3537 \end{methoddesc}
    3538 
    3539 \begin{methoddesc}{\emph{<domain>}.get_maximum_inundation_location}{indices=None}
    3540 Module: \module{shallow_water.shallow_water_domain}
    3541 
    3542 Return location (x,y) of highest elevation where h $>$ 0.
    3543 
    3544 \code{indices} is the optional set of element \code{id}s that
    3545 the operation applies to.
    3546 
    3547 Example to find maximum runup location:
    3548 \begin{verbatim}
    3549 x, y = domain.get_maximum_inundation_location()
    3550 \end{verbatim}
    3551 \end{methoddesc}
    3552 
    3553 
    3554 \section{Queries of SWW model output files}
    3555 After a model has been run, it is often useful to extract various information from the SWW
    3556 output file (see Section \ref{sec:sww format}). This is typically more convenient than using the
    3557 diagnostics described in Section \ref{sec:diagnostics} which rely on the model running -- something
    3558 that can be very time consuming. The SWW files are easy and quick to read and offer information
    3559 about the model results such as runup heights, time histories of selected quantities,
    3560 flow through cross sections and much more.
    3561 
    3562 \begin{funcdesc}{elevation = get_maximum_inundation_elevation}{filename,
    3563                                      polygon=None,
    3564                                      time_interval=None,
    3565                                      verbose=False}
    3566 Module: \module{shallow_water.data_manager}
    3567 
    3568 Return the highest elevation where depth is positive ($h > 0$).
    3569 
    3570 \code{filename} is the path to a NetCDF SWW file containing \anuga model output.
    3571 
    3572 \code{polygon} restricts the query to the points that lie within the polygon.
    3573 
    3574 \code {time_interval} restricts the query to within the time interval.
    3575 
    3576 Example to find maximum runup elevation:
    3577 
    3578 \begin{verbatim}
    3579 max_runup = get_maximum_inundation_elevation(filename)
    3580 \end{verbatim}
    3581 
    3582 If no inundation is found (within the \code{polygon} and \code{time_interval}, if specified)
    3583 the return value is \code{None}. This indicates "No Runup" or "Everything is dry".
    3584 \end{funcdesc}
    3585 
    3586 \begin{funcdesc}{(x, y) = get_maximum_inundation_location}{filename,
    3587                                     polygon=None,
    3588                                     time_interval=None,
    3589                                     verbose=False}
    3590 Module: \module{shallow_water.data_manager}
    3591 
    3592 Return location (x,y) of the highest elevation where depth is positive ($h > 0$).
    3593 
    3594 \code{filename} is the path to a NetCDF SWW file containing \anuga model output.
    3595 
    3596 \code{polygon} restricts the query to the points that lie within the polygon.
    3597 
    3598 \code {time_interval} restricts the query to within the time interval.
    3599 
    3600 Example to find maximum runup location:
    3601 
    3602 \begin{verbatim}
    3603 max_runup_location = get_maximum_inundation_location(filename)
    3604 \end{verbatim}
    3605 
    3606 If no inundation is found (within the \code{polygon} and \code{time_interval}, if specified)
    3607 the return value is \code{None}. This indicates "No Runup" or "Everything is dry".
    3608 \end{funcdesc}
    3609 
    3610 \begin{funcdesc}{sww2timeseries}{swwfiles,
    3611                                  gauge_filename,
    3612                                  production_dirs,
    3613                                  report=None,
    3614                                  reportname=None,
    3615                                  plot_quantity=None,
    3616                                  generate_fig=False,
    3617                                  surface=None,
    3618                                  time_min=None,
    3619                                  time_max=None,
    3620                                  output_centroids=False,
    3621                                  time_thinning=1,
    3622                                  time_unit=None,
    3623                                  title_on=None,
    3624                                  use_cache=False,
    3625                                  verbose=False}
    3626 Module: \module{abstract_2d_finite_volumes.util}
    3627 
    3628 Read a set of SWW files and plot the time series for the prescribed quantities
    3629 at defined gauge locations and prescribed time range.
    3630 
    3631 \code{swwfiles} is a dictionary of SWW files with keys of \code{label_id}.
    3632 
    3633 \code{gauge_filename} is the path to a file containing gauge data.
    3634 
    3635 THIS NEEDS MORE WORK.  WORK ON FUNCTION __DOC__ STRING, IF NOTHING ELSE!
    3636 \end{funcdesc}
    3637 
    3638 \begin{funcdesc}{(time, Q) = get_flow_through_cross_section}{filename, polyline, verbose=False}
    3639 Module: \module{shallow_water.data_manager}
    3640 
    3641 Obtain flow ($m^3/s$) perpendicular to specified cross section.
    3642 
    3643 \code{filename} is the path to the SWW file.
    3644 
    3645 \code{output_centroids} set to true to sample at the centre of the triangle containing the point.
    3646 This may be useful for debugging. (new in 1.2)
    3647 
    3648 \code{polyline} is the representation of the desired cross section -- it may contain
    3649 multiple sections allowing for complex shapes. Assumes absolute UTM coordinates.
    3650 
    3651 Returns a tuple \code{time, Q} where:
    3652 
    3653 \code{time} is all the stored times in the SWW file.
    3654 
    3655 \code{Q} is a hydrograph of total flow across the given segments for all stored times.
    3656 
    3657 The normal flow is computed for each triangle intersected by the \code{polyline} and
    3658 added up.  If multiple segments at different angles are specified the normal flows
    3659 may partially cancel each other.
    3660 
    3661 Example to find flow through cross section:
    3662 
    3663 \begin{verbatim}
    3664 cross_section = [[x, 0], [x, width]]
    3665 time, Q = get_flow_through_cross_section(filename, cross_section)
    3666 \end{verbatim}
    3667 \end{funcdesc}
     3729
     3730
     3731
     3732
     3733
     3734
    36683735
    36693736
     
    36923759%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    36933760
    3694 \chapter{\anuga System Architecture}
    3695 
    3696 \section{File Formats}
     3761\chapter{\anuga File Formats}
     3762
    36973763\label{sec:file formats}
    36983764
     
    37323798%SOMETHING MISSING HERE!?
    37333799
    3734 \subsection{Manually Created Files}
     3800\section{Manually Created Files}
    37353801
    37363802\begin{tabular}{ll}
     
    37393805\end{tabular}
    37403806
    3741 \subsection{Automatically Created Files}
     3807\section{Automatically Created Files}
    37423808
    37433809\begin{tabular}{ll}
     
    37523818\bigskip
    37533819
    3754 \subsection{SWW, STS and TMS Formats}
     3820\section{SWW, STS and TMS Formats}
    37553821\label{sec:sww format}
    37563822The SWW, STS and TMS formats are all NetCDF formats and are of key importance for \anuga.
     
    37803846
    37813847%FIXME (Ole): Should put in example with nonzero xllcorner, yllcorner
    3782 \verbatiminputunderscore{examples/bedslopeexcerpt.cdl}
     3848\verbatiminputunderscore{../../anuga_core/examples/bedslopeexcerpt.cdl}
    37833849
    37843850The SWW format is used not only for output but also serves as input
     
    38063872A TMS file is used to store time series data that is independent of position.
    38073873
    3808 \subsection{Mesh File Formats}
     3874\section{Mesh File Formats}
    38093875
    38103876A mesh file is a file that has a specific format suited to
     
    38313897offset to be applied to $x$ and $y$ values -- e.g.\ to the vertices.
    38323898
    3833 \subsection{Formats for Storing Arbitrary Points and Attributes}
     3899\section{Formats for Storing Arbitrary Points and Attributes}
    38343900
    38353901A CSV/TXT/XYA file is used to store data representing
     
    38593925attribute.
    38603926
    3861 \subsection{ArcView Formats}
     3927\section{ArcView Formats}
    38623928
    38633929Files of the three formats ASC, PRJ and ERS are all associated with
     
    38823948to represent metadata for a DEM.
    38833949
    3884 \subsection{DEM Format}
     3950\section{DEM Format}
    38853951
    38863952A DEM file in \anuga is a NetCDF representation of regular DEM data.
    38873953
    3888 \subsection{Other Formats}
    3889 
    3890 \subsection{Basic File Conversions}
     3954
     3955
     3956\section{Basic File Conversions}
    38913957\label{sec:basicfileconversions}
    38923958
     
    55745640\end{methoddesc}
    55755641
    5576 \subsection{Miscellaneous Functions}
     5642\section{Miscellaneous Functions}
    55775643
    55785644The functions here are not \code{Geospatial_data} object methods, but are used with them.
  • trunk/anuga_user_manual/source/copyright.tex

    r9090 r9727  
    7373\begin{itemize}
    7474\item \anuga was developed by Stephen Roberts, Ole Nielsen, Duncan Gray and Jane Sexton. It is currently being developed and
    75  maintained by Nariman Habili and Stephen Roberts.
     75 maintained by Stephen Roberts and Gareth Davies.
    7676\index{ANUGA!credits|textit} 
    7777\end{itemize}
     
    8787\item Ole Nielsen, James Hudson, John Jakeman, Rudy van Drie, Ted Rigby,
    8888      Petar Milevski, Joaquim Luis, Nils Goseberg, William Power,
    89       Trevor Dhu, Linda Stals, Matt Hardy, Jack Kelly and Christopher
     89      Trevor Dhu, Nariman Habili, Gareth Davies, Linda Stals, Matt Hardy, Jack Kelly and Christopher
    9090      Zoppou who contributed to this project at various times.
    9191\item A stand alone visualiser (anuga\_viewer) based on Open-scene-graph was developed by Darran Edmundson and James Hudson.
  • trunk/anuga_user_manual/source/update_anuga_user_manual.py

    r9665 r9727  
    4141from os.path import expanduser, split, join
    4242from anuga.utilities.system_tools import get_revision_number, get_pathname_from_package
    43 from anuga.config import major_revision
     43from anuga import __version__ as major_revision
    4444from sys import argv
    4545
     
    115115
    116116system('cp *.pdf ../')
    117 system('mv *.pdf ../../anuga_core/user_manual')
     117system('cp anuga_user_manual.pdf ../../anuga_core/doc')
    118118
    119119print 'Cleanup aux tex files'
  • trunk/anuga_user_manual/source/version.tex

    r9269 r9727  
    77% release version; this is used to define the
    88% \version macro
    9 \release{1.3.1}
     9\release{1.3.11}
Note: See TracChangeset for help on using the changeset viewer.