source: documentation/user_manual/anuga_user_manual.tex @ 3383

Last change on this file since 3383 was 3275, checked in by sexton, 19 years ago

updates on sydney demo for user manual

  • Property svn:keywords set to LastChangedDate LastChangedRevision LastChangedBy HeadURL Id
File size: 106.6 KB
Line 
1% Complete documentation on the extended LaTeX markup used for Python
2% documentation is available in ``Documenting Python'', which is part
3% of the standard documentation for Python.  It may be found online
4% at:
5%
6%     http://www.python.org/doc/current/doc/doc.html
7
8
9%labels
10%Sections and subsections \label{sec: }
11%Chapters \label{ch: }
12%Equations \label{eq: }
13%Figures \label{fig: }
14
15
16
17\documentclass{manual}
18
19\usepackage{graphicx}
20\usepackage{datetime}
21
22\input{definitions}
23
24\title{\anuga User Manual}
25\author{Geoscience Australia and the Australian National University}
26
27% Please at least include a long-lived email address;
28% the rest is at your discretion.
29\authoraddress{Geoscience Australia \\
30  Email: \email{ole.nielsen@ga.gov.au}
31}
32
33%Draft date
34
35% update before release!
36% Use an explicit date so that reformatting
37% doesn't cause a new date to be used.  Setting
38% the date to \today can be used during draft
39% stages to make it easier to handle versions.
40
41
42\longdate       % Make date format long using datetime.sty
43%\settimeformat{xxivtime} % 24 hour Format
44\settimeformat{oclock} % Verbose
45\date{\today, \ \currenttime}
46%\hyphenation{set\_datadir}
47
48\ifhtml
49\date{\today} % latex2html does not know about datetime
50\fi
51
52
53
54
55
56\release{1.0}   % release version; this is used to define the
57                % \version macro
58
59\makeindex          % tell \index to actually write the .idx file
60\makemodindex       % If this contains a lot of module sections.
61
62\setcounter{tocdepth}{3}
63\setcounter{secnumdepth}{3}
64
65
66\begin{document}
67\maketitle
68
69
70% This makes the contents more accessible from the front page of the HTML.
71\ifhtml
72\chapter*{Front Matter\label{front}}
73\fi
74
75%Subversion keywords:
76%
77%$LastChangedDate: 2006-07-04 04:34:21 +0000 (Tue, 04 Jul 2006) $
78%$LastChangedRevision: 3275 $
79%$LastChangedBy: ole $
80
81\input{copyright}
82
83
84\begin{abstract}
85\label{def:anuga}
86
87\noindent \anuga\index{\anuga} is a hydrodynamic modelling tool that
88allows users to model realistic flow problems in complex geometries.
89Examples include dam breaks or the effects of natural hazards such
90as riverine flooding, storm surges and tsunami.
91
92The user must specify a study area represented by a mesh of triangular
93cells, the topography and bathymetry, frictional resistance, initial
94values for water level (called \emph{stage}\index{stage} within \anuga),
95boundary
96conditions and forces such as windstress or pressure gradients if
97applicable.
98
99\anuga tracks the evolution of water depth and horizontal momentum
100within each cell over time by solving the shallow water wave equation
101governing equation using a finite-volume method.
102
103\anuga cannot model details of breaking waves, flow under ceilings such
104as pipes, turbulence and vortices, vertical convection or viscous
105flows.
106
107\anuga also incorporates a mesh generator, called \code{pmesh}, that
108allows the user to set up the geometry of the problem interactively as
109well as tools for interpolation and surface fitting, and a number of
110auxiliary tools for visualising and interrogating the model output.
111
112Most \anuga components are written in the object-oriented programming
113language Python and most users will interact with \anuga by writing
114small Python programs based on the \anuga library
115functions. Computationally intensive components are written for
116efficiency in C routines working directly with the Numerical Python
117structures.
118
119
120\end{abstract}
121
122\tableofcontents
123
124
125\chapter{Introduction}
126
127
128\section{Purpose}
129
130The purpose of this user manual is to introduce the new user to the
131inundation software, describe what it can do and give step-by-step
132instructions for setting up and running hydrodynamic simulations.
133
134\section{Scope}
135
136This manual covers only what is needed to operate the software after
137installation and configuration. It does not includes instructions
138for installing the software or detailed API documentation, both of
139which will be covered in separate publications and by documentation
140in the source code.
141
142\section{Audience}
143
144Readers are assumed to be familiar with the operating environment
145and have a general understanding of the subject matter, as well as
146enough programming experience to adapt the code to different
147requirements and to understand the basic terminology of
148object-oriented programming.
149
150\pagebreak
151\chapter{Background}
152
153
154Modelling the effects on the built environment of natural hazards such
155as riverine flooding, storm surges and tsunami is critical for
156understanding their economic and social impact on our urban
157communities.  Geoscience Australia and the Australian National
158University are developing a hydrodynamic inundation modelling tool
159called \anuga to help simulate the impact of these hazards.
160
161The core of \anuga is the fluid dynamics module, called pyvolution,
162which is based on a finite-volume method for solving the shallow water
163wave equation.  The study area is represented by a mesh of triangular
164cells.  By solving the governing equation within each cell, water
165depth and horizontal momentum are tracked over time.
166
167A major capability of pyvolution is that it can model the process of
168wetting and drying as water enters and leaves an area.  This means
169that it is suitable for simulating water flow onto a beach or dry land
170and around structures such as buildings.  Pyvolution is also capable
171of modelling hydraulic jumps due to the ability of the finite-volume
172method to accommodate discontinuities in the solution.
173
174To set up a particular scenario the user specifies the geometry
175(bathymetry and topography), the initial water level (stage),
176boundary conditions such as tide, and any forcing terms that may
177drive the system such as wind stress or atmospheric pressure
178gradients. Gravity and frictional resistance from the different
179terrains in the model are represented by predefined forcing terms.
180
181A mesh generator, called pmesh, allows the user to set up the geometry
182of the problem interactively and to identify boundary segments and
183regions using symbolic tags.  These tags may then be used to set the
184actual boundary conditions and attributes for different regions
185(e.g. the Manning friction coefficient) for each simulation.
186
187Most \anuga components are written in the object-oriented programming
188language Python.  Software written in Python can be produced quickly
189and can be readily adapted to changing requirements throughout its
190lifetime.  Computationally intensive components are written for
191efficiency in C routines working directly with the Numerical Python
192structures.  The animation tool developed for \anuga is based on
193OpenSceneGraph, an Open Source Software (OSS) component allowing high
194level interaction with sophisticated graphics primitives.
195See \cite{nielsen2005} for more background on \anuga.
196
197\chapter{Restrictions and limitations on \anuga}
198
199Although a powerful and flexible tool for hydrodynamic modelling, \anuga has a
200number of limitations that any potential user need to be aware of. They are
201
202\begin{itemize}
203  \item The mathematical model is the 2D shallow water wave equation.
204  As such it cannot resolve vertical convection and consequently not breaking
205  waves or 3D turbulence (e.g.\ vorticity).
206  \item The finite volume is a very robust and flexible numerical technique,
207  but it is not the fastest method around. If the geometry is sufficiently
208  simple and if there is no need for wetting or drying, a finite-difference
209  method may be able to solve the problem faster than \anuga.
210  %\item Mesh resolutions near coastlines with steep gradients need to be...   
211  \item Frictional resistance is implemented using Manning's formula, but
212  \anuga has not yet been validated in regard to bottom roughness
213\end{itemize}
214
215
216
217\chapter{Getting Started}
218\label{ch:getstarted}
219
220This section is designed to assist the reader to get started with
221\anuga by working through some examples. Two examples are discussed;
222the first is a simple example to illustrate many of the ideas, and
223the second is a more realistic example.
224
225\section{A Simple Example}
226\label{sec:simpleexample}
227
228\subsection{Overview}
229
230What follows is a discussion of the structure and operation of a
231script called \file{runup.py}.
232
233This example carries out the solution of the shallow-water wave
234equation in the simple case of a configuration comprising a flat
235bed, sloping at a fixed angle in one direction and having a
236constant depth across each line in the perpendicular direction.
237
238The example demonstrates the basic ideas involved in setting up a
239complex scenario. In general the user specifies the geometry
240(bathymetry and topography), the initial water level, boundary
241conditions such as tide, and any forcing terms that may drive the
242system such as wind stress or atmospheric pressure gradients.
243Frictional resistance from the different terrains in the model is
244represented by predefined forcing terms. In this example, the
245boundary is reflective on three sides and a time dependent wave on
246one side.
247
248The present example represents a simple scenario and does not
249include any forcing terms, nor is the data taken from a file as it
250would typically be.
251
252The conserved quantities involved in the
253problem are stage (absolute height of water surface),
254$x$-momentum and $y$-momentum. Other quantities
255involved in the computation are the friction and elevation.
256
257Water depth can be obtained through the equation
258
259\begin{tabular}{rcrcl}
260  \code{depth} &=& \code{stage} &$-$& \code{elevation}
261\end{tabular}
262
263
264\subsection{Outline of the Program}
265
266In outline, \file{runup.py} performs the following steps:
267
268\begin{enumerate}
269
270   \item Sets up a triangular mesh.
271
272   \item Sets certain parameters governing the mode of
273operation of the model-specifying, for instance, where to store the model output.
274
275   \item Inputs various quantities describing physical measurements, such
276as the elevation, to be specified at each mesh point (vertex).
277
278   \item Sets up the boundary conditions.
279
280   \item Carries out the evolution of the model through a series of time
281steps and outputs the results, providing a results file that can
282be visualised.
283
284\end{enumerate}
285
286\subsection{The Code}
287
288%FIXME: we are using the \code function here.
289%This should be used wherever possible
290For reference we include below the complete code listing for
291\file{runup.py}. Subsequent paragraphs provide a
292`commentary' that describes each step of the program and explains it
293significance.
294
295\verbatiminput{examples/runup.py}
296%\verbatiminput{examples/bedslope.py}
297
298\subsection{Establishing the Mesh}\index{mesh, establishing}
299
300The first task is to set up the triangular mesh to be used for the
301scenario. This is carried out through the statement:
302
303{\small \begin{verbatim}
304    points, vertices, boundary = rectangular(10, 10)
305\end{verbatim}}
306
307The function \function{rectangular} is imported from a module
308\module{mesh\_factory} defined elsewhere. (\anuga also contains
309several other schemes that can be used for setting up meshes, but we
310shall not discuss these.) The above assignment sets up a $10 \times
31110$ rectangular mesh, triangulated in a regular way. The assignment
312
313{\small \begin{verbatim}
314    points, vertices, boundary = rectangular(m, n)
315\end{verbatim}}
316
317returns:
318
319\begin{itemize}
320
321   \item a list \code{points} giving the coordinates of each mesh point,
322
323   \item a list \code{vertices} specifying the three vertices of each triangle, and
324
325   \item a dictionary \code{boundary} that stores the edges on
326   the boundary and associates each with one of the symbolic tags \code{`left'}, \code{`right'},
327   \code{`top'} or \code{`bottom'}.
328
329\end{itemize}
330
331(For more details on symbolic tags, see page
332\pageref{ref:tagdescription}.)
333
334An example of a general unstructured mesh and the associated data
335structures \code{points}, \code{vertices} and \code{boundary} is
336given in Section \ref{sec:meshexample}.
337
338
339
340
341\subsection{Initialising the Domain}
342
343These variables are then used to set up a data structure
344\code{domain}, through the assignment:
345
346{\small \begin{verbatim}
347    domain = Domain(points, vertices, boundary)
348\end{verbatim}}
349
350This creates an instance of the \class{Domain} class, which
351represents the domain of the simulation. Specific options are set at
352this point, including the basename for the output file and the
353directory to be used for data:
354
355{\small \begin{verbatim}
356    domain.set_name('bedslope')
357\end{verbatim}}
358
359{\small \begin{verbatim}
360    domain.set_datadir('.')
361\end{verbatim}}
362
363In addition, the following statement now specifies that the
364quantities \code{stage}, \code{xmomentum} and \code{ymomentum} are
365to be stored:
366
367{\small \begin{verbatim}
368    domain.set_quantities_to_be_stored(['stage', 'xmomentum',
369    'ymomentum'])
370\end{verbatim}}
371
372
373\subsection{Initial Conditions}
374
375The next task is to specify a number of quantities that we wish to
376set for each mesh point. The class \class{Domain} has a method
377\method{set\_quantity}, used to specify these quantities. It is a
378flexible method that allows the user to set quantities in a variety
379of ways---using constants, functions, numeric arrays, expressions
380involving other quantities, or arbitrary data points with associated
381values, all of which can be passed as arguments. All quantities can
382be initialised using \method{set\_quantity}. For a conserved
383quantity (such as \code{stage, xmomentum, ymomentum}) this is called
384an \emph{initial condition}. However, other quantities that aren't
385updated by the equation are also assigned values using the same
386interface. The code in the present example demonstrates a number of
387forms in which we can invoke \method{set\_quantity}.
388
389
390\subsubsection{Elevation}
391
392The elevation, or height of the bed, is set using a function,
393defined through the statements below, which is specific to this
394example and specifies a particularly simple initial configuration
395for demonstration purposes:
396
397{\small \begin{verbatim}
398    def f(x,y):
399        return -x/2
400\end{verbatim}}
401
402This simply associates an elevation with each point \code{(x, y)} of
403the plane.  It specifies that the bed slopes linearly in the
404\code{x} direction, with slope $-\frac{1}{2}$,  and is constant in
405the \code{y} direction.
406
407Once the function \function{f} is specified, the quantity
408\code{elevation} is assigned through the simple statement:
409
410{\small \begin{verbatim}
411    domain.set_quantity('elevation', f)
412\end{verbatim}}
413
414
415\subsubsection{Friction}
416
417The assignment of the friction quantity (a forcing term)
418demonstrates another way we can use \method{set\_quantity} to set
419quantities---namely, assign them to a constant numerical value:
420
421{\small \begin{verbatim}
422    domain.set_quantity('friction', 0.1)
423\end{verbatim}}
424
425This specifies that the Manning friction coefficient is set to 0.1
426at every mesh point.
427
428\subsubsection{Stage}
429
430The stage (the height of the water surface) is related to the
431elevation and the depth at any time by the equation
432
433{\small \begin{verbatim}
434    stage = elevation + depth
435\end{verbatim}}
436
437
438For this example, we simply assign a constant value to \code{stage},
439using the statement
440
441{\small \begin{verbatim}
442    domain.set_quantity('stage', -.4)
443\end{verbatim}}
444
445which specifies that the surface level is set to a height of $-0.4$,
446i.e. 0.4 units (m) below the zero level.
447
448Although it is not necessary for this example, it may be useful to
449digress here and mention a variant to this requirement, which allows
450us to illustrate another way to use \method{set\_quantity}---namely,
451incorporating an expression involving other quantities. Suppose,
452instead of setting a constant value for the stage, we wished to
453specify a constant value for the \emph{depth}. For such a case we
454need to specify that \code{stage} is everywhere obtained by adding
455that value to the value already specified for \code{elevation}. We
456would do this by means of the statements:
457
458{\small \begin{verbatim}
459    h = 0.05 # Constant depth
460    domain.set_quantity('stage', expression = 'elevation + %f' %h)
461\end{verbatim}}
462
463That is, the value of \code{stage} is set to $\code{h} = 0.05$ plus
464the value of \code{elevation} already defined.
465
466The reader will probably appreciate that this capability to
467incorporate expressions into statements using \method{set\_quantity}
468greatly expands its power.) See Section \ref{sec:Initial Conditions} for more
469details.
470
471\subsection{Boundary Conditions}\index{boundary conditions}
472
473The boundary conditions are specified as follows:
474
475{\small \begin{verbatim}
476    Br = Reflective_boundary(domain)
477
478    Bt = Transmissive_boundary(domain)
479
480    Bd = Dirichlet_boundary([0.2,0.,0.])
481
482    Bw = Time_boundary(domain=domain,
483                f=lambda t: [(0.1*sin(t*2*pi)-0.3), 0.0, 0.0])
484\end{verbatim}}
485
486The effect of these statements is to set up a selection of different
487alternative boundary conditions and store them in variables that can be
488assigned as needed. Each boundary condition specifies the
489behaviour at a boundary in terms of the behaviour in neighbouring
490elements. The boundary conditions introduced here may be briefly described as
491follows:
492
493\begin{itemize}
494    \item \textbf{Reflective boundary}\label{def:reflective boundary} Returns same \code{stage} as
495      as present in its neighbour volume but momentum vector
496      reversed 180 degrees (reflected).
497      Specific to the shallow water equation as it works with the
498      momentum quantities assumed to be the second and third conserved
499      quantities. A reflective boundary condition models a solid wall.
500    \item \textbf{Transmissive boundary}\label{def:transmissive boundary} Returns same conserved quantities as
501      those present in its neighbour volume. This is one way of modelling
502      outflow from a domain, but it should be used with caution if flow is
503      not steady state as replication of momentum at the boundary
504      may cause occasional spurious effects. If this occurs,
505      consider using e.g. a Dirichlet boundary condition.
506    \item \textbf{Dirichlet boundary}\label{def:dirichlet boundary} Specifies
507      constant values for stage, $x$-momentum and $y$-momentum at the boundary.
508    \item \textbf{Time boundary}\label{def:time boundary} Like a Dirichlet
509      boundary but with behaviour varying with time.
510\end{itemize}
511
512\label{ref:tagdescription}Before describing how these boundary
513conditions are assigned, we recall that a mesh is specified using
514three variables \code{points}, \code{vertices} and \code{boundary}.
515In the code we are discussing, these three variables are returned by
516the function \code{rectangular}; however, the example given in
517Section \ref{sec:realdataexample} illustrates another way of
518assigning the values, by means of the function
519\code{create\_mesh\_from\_regions}.
520
521These variables store the data determining the mesh as follows. (You
522may find that the example given in Section \ref{sec:meshexample}
523helps to clarify the following discussion, even though that example
524is a \emph{non-rectangular} mesh.)
525
526\begin{itemize}
527\item The variable \code{points} stores a list of 2-tuples giving the
528coordinates of the mesh points.
529
530\item The variable \code{vertices} stores a list of 3-tuples of
531numbers, representing vertices of triangles in the mesh. In this
532list, the triangle whose vertices are \code{points[i]},
533\code{points[j]}, \code{points[k]} is represented by the 3-tuple
534\code{(i, j, k)}.
535
536\item The variable \code{boundary} is a Python dictionary that
537not only stores the edges that make up the boundary but also assigns
538symbolic tags to these edges to distinguish different parts of the
539boundary. An edge with endpoints \code{points[i]} and
540\code{points[j]} is represented by the 2-tuple \code{(i, j)}. The
541keys for the dictionary are the 2-tuples \code{(i, j)} corresponding
542to boundary edges in the mesh, and the values are the tags are used
543to label them. In the present example, the value \code{boundary[(i,
544j)]} assigned to \code{(i, j)]} is one of the four tags
545\code{`left'}, \code{`right'}, \code{`top'} or \code{`bottom'},
546depending on whether the boundary edge represented by \code{(i, j)}
547occurs at the left, right, top or bottom of the rectangle bounding
548the mesh. The function \code{rectangular} automatically assigns
549these tags to the boundary edges when it generates the mesh.
550\end{itemize}
551
552The tags provide the means to assign different boundary conditions
553to an edge depending on which part of the boundary it belongs to.
554(In Section \ref{sec:realdataexample} we describe an example that
555uses different boundary tags---in general, the possible tags are not
556limited to `left', `right', `top' and `bottom', but can be specified
557by the user.)
558
559Using the boundary objects described above, we assign a boundary
560condition to each part of the boundary by means of a statement like
561
562{\small \begin{verbatim}
563    domain.set_boundary({'left': Br, 'right': Bw, 'top': Br, 'bottom': Br})
564\end{verbatim}}
565
566This statement stipulates that, in the current example, the right
567boundary varies with time, as defined by the lambda function, while the other
568boundaries are all reflective.
569
570The reader may wish to experiment by varying the choice of boundary
571types for one or more of the boundaries. (In the case of \code{Bd}
572and \code{Bw}, the three arguments in each case represent the
573\code{stage}, $x$-momentum and $y$-momentum, respectively.)
574
575{\small \begin{verbatim}
576    Bw = Time_boundary(domain=domain,
577                       f=lambda t: [(0.1*sin(t*2*pi)-0.3), 0.0, 0.0])
578\end{verbatim}}
579
580
581
582\subsection{Evolution}\index{evolution}
583
584The final statement \nopagebreak[3]
585{\small \begin{verbatim}
586    for t in domain.evolve(yieldstep = 0.1, duration = 4.0):
587        print domain.timestepping_statistics()
588\end{verbatim}}
589
590causes the configuration of the domain to `evolve', over a series of
591steps indicated by the values of \code{yieldstep} and
592\code{duration}, which can be altered as required.  The value of
593\code{yieldstep} controls the time interval between successive model
594outputs.  Behind the scenes more time steps are generally taken.
595
596
597\subsection{Output}
598
599The output is a NetCDF file with the extension \code{.sww}. It
600contains stage and momentum information and can be used with the
601\code{swollen} (see Section \ref{sec:swollen}) visualisation package
602to generate a visual display. See Section \ref{sec:file formats}
603(page \pageref{sec:file formats}) for more on NetCDF and other file
604formats.
605
606The following is a listing of the screen output seen by the user
607when this example is run:
608
609\verbatiminput{examples/runupoutput.txt}
610
611
612\section{How to Run the Code}
613
614The code can be run in various ways:
615
616\begin{itemize}
617  \item{from a Windows or Unix command line} as in\ \code{python runup.py}
618  \item{within the Python IDLE environment}
619  \item{within emacs}
620  \item{within Windows, by double-clicking the \code{runup.py}
621  file.}
622\end{itemize}
623
624
625\section{Exploring the Model Output}
626
627The following figures are screenshots from the \anuga visualisation
628tool \code{swollen}. Figure \ref{fig:runupstart} shows the domain
629with water surface as specified by the initial condition, $t=0$.
630Figure \ref{fig:bedslope2} shows later snapshots for $t=2.3$ and
631$t=4$ where the system has been evolved and the wave is encroaching
632on the previously dry bed.  All figures are screenshots from an
633interactive animation tool called Swollen which is part of \anuga.
634Swollen is described in more detail is Section \ref{sec:swollen}.
635
636\begin{figure}[hbt]
637
638  \centerline{\includegraphics[width=75mm, height=75mm]
639    {examples/runupstart.eps}}
640
641  \caption{Bedslope example viewed with Swollen}
642  \label{fig:runupstart}
643\end{figure}
644
645
646\begin{figure}[hbt]
647
648  \centerline{
649    \includegraphics[width=75mm, height=75mm]{examples/runupduring.eps}
650    \includegraphics[width=75mm, height=75mm]{examples/runupend.eps}
651   }
652
653  \caption{Bedslope example viewed with Swollen}
654  \label{fig:bedslope2}
655\end{figure}
656
657
658
659
660\clearpage
661
662%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
663
664\section{An Example with Real Data}
665\label{sec:realdataexample} The following discussion builds on the
666concepts introduced through the \file{runup.py} example and
667introduces a second example, \file{run\_sydney.py}.  This refers to
668a real-life scenario, in which the domain of interest surrounds the
669Sydney region, and predominantly covers Sydney Harbour. A
670hypothetical tsunami wave is generated by a submarine mass failure
671situated on the edge of the continental shelf.
672
673\subsection{Overview}
674As in the case of \file{runup.py}, the actions carried
675out by the program can be organised according to this outline:
676
677\begin{enumerate}
678
679   \item Set up a triangular mesh.
680
681   \item Set certain parameters governing the mode of
682operation of the model---specifying, for instance, where to store the
683model output.
684
685   \item Input various quantities describing physical measurements, such
686as the elevation, to be specified at each mesh point (vertex).
687
688   \item Set up the boundary conditions.
689
690   \item Carry out the evolution of the model through a series of time
691steps and output the results, providing a results file that can be
692visualised.
693
694\end{enumerate}
695
696
697
698\subsection{The Code}
699
700Here is the code for \file{run\_sydney\_smf.py}:
701
702\verbatiminput{examples/runsydney.py}
703
704In discussing the details of this example, we follow the outline
705given above, discussing each major step of the code in turn.
706
707\subsection{Establishing the Mesh}\index{mesh, establishing}
708
709One obvious way that the present example differs from
710\file{runup.py} is in the use of a more complex method to
711create the mesh. Instead of imposing a mesh structure on a
712rectangular grid, the technique used for this example involves
713building mesh structures inside polygons specified by the user,
714using a mesh-generator referred to as \code{pmesh}.
715
716In its simplest form, \code{pmesh} creates the mesh within a single
717polygon whose vertices are at geographical locations specified by
718the user. The user specifies the \emph{resolution}---that is, the
719maximal area of a triangle used for triangulation---and a triangular
720mesh is created inside the polygon using a mesh generation engine.
721On any given platform, the same mesh will be returned. Figure
722\ref{fig:pentagon} shows a simple example of this, in which the
723triangulation is carried out within a pentagon.
724
725
726\begin{figure}[hbt]
727
728  \caption{Mesh points are created inside the polygon}
729  \label{fig:pentagon}
730\end{figure}
731
732Boundary tags are not restricted to \code{`left'}, \code{`right'},
733\code{`bottom'} and \code{`top'}, as in the case of
734\file{runup.py}. Instead the user specifies a list of
735tags appropriate to the configuration being modelled.
736
737In addition, \code{pmesh} provides a way to adapt to geographic or
738other features in the landscape, whose presence may require an
739increase in resolution. This is done by allowing the user to specify
740a number of \emph{interior polygons}, each with a specified
741resolution, see Figure \ref{fig:interior meshes}. It is also
742possible to specify one or more `holes'---that is, areas bounded by
743polygons in which no triangulation is required.
744
745\begin{figure}[hbt]
746
747
748
749  \caption{Interior meshes with individual resolution}
750  \label{fig:interior meshes}
751\end{figure}
752
753In its general form, \code{pmesh} takes for its input a bounding
754polygon and (optionally) a list of interior polygons. The user
755specifies resolutions, both for the bounding polygon and for each of
756the interior polygons. Given this data, \code{pmesh} first creates a
757triangular mesh with varying resolution.
758
759The function used to implement this process is
760\function{create\_mesh\_from\_regions}. Its arguments include the
761bounding polygon and its resolution, a list of boundary tags, and a
762list of pairs \code{[polygon, resolution]}, specifying the interior
763polygons and their resolutions.
764
765In practice, the details of the polygons used are read from a
766separate file \file{project.py}. Here is a complete listing of
767\file{project.py}:
768
769\verbatiminput{examples/project.py}
770
771The resulting mesh is output to a \emph{mesh file}\index{mesh
772file}\label{def:mesh file}. This term is used to describe a file of
773a specific format used to store the data specifying a mesh. (There
774are in fact two possible formats for such a file: it can either be a
775binary file, with extension \code{.msh}, or an ASCII file, with
776extension \code{.tsh}. In the present case, the binary file format
777\code{.msh} is used. See Section \ref{sec:file formats} (page
778\pageref{sec:file formats}) for more on file formats.)
779
780The statements
781
782{\small \begin{verbatim}
783    interior_res = 5000
784    interior_regions = [[project.northern_polygon, interior_res],
785                        [project.harbour_polygon, interior_res],
786                        [project.southern_polygon, interior_res],
787                        [project.manly_polygon, interior_res],
788                        [project.top_polygon, interior_res]]
789\end{verbatim}}
790
791are used to read in the specific polygons \code{project.harbour\_polygon\_2} and
792\code{botanybay\_polygon\_2} from \file{project.py} and assign a
793common resolution of to each. The statement
794
795{\small \begin{verbatim}
796create_mesh_from_regions(project.demopoly,
797                         boundary_tags= {'oceannorth': [0],
798                                         'onshorenorth1': [1],
799                                         'onshorenorthwest1': [2],
800                                         'onshorenorthwest2': [3],
801                                         'onshorenorth2': [4],
802                                         'onshorewest': [5],
803                                         'onshoresouth': [6],
804                                         'oceansouth': [7],
805                                         'oceaneast': [8]},
806                         maximum_triangle_area=100000,
807                         filename=meshname,           
808                         interior_regions=interior_regions)   
809\end{verbatim}}
810
811is then used to create the mesh, taking the bounding polygon to be
812the polygon \code{diffpolygonall} specified in \file{project.py}.
813The argument \code{boundary\_tags} assigns a dictionary, whose keys
814are the names of the boundary tags used for the bounding
815polygon---\code{`bottom'}, \code{`right0'}, \code{`right1'},
816\code{`right2'}, \code{`top'}, \code{`left1'}, \code{`left2'} and
817\code{`left3'}--- and whose values identify the indices of the
818segments associated with each of these tags. (The value associated
819with each boundary tag is a one-element list.)
820
821
822\subsection{Initialising the Domain}
823
824As with \file{runup.py}, once we have created the mesh, the next
825step is to create the data structure \code{domain}. We did this for
826\file{runup.py} by inputting lists of points and triangles and
827specifying the boundary tags directly. However, in the present case,
828we use a method that works directly with the mesh file
829\code{meshname}, as follows:
830
831
832{\small \begin{verbatim}
833    domain = Domain(meshname, use_cache=True, verbose=True)
834\end{verbatim}}
835
836Providing a filename instead of the lists used in \file{runup.py}
837above causes \code{Domain} to convert a mesh file \code{meshname}
838into an instance of \code{Domain}, allowing us to use methods like
839\method{set\_quantity} to set quantities and to apply other
840operations.
841
842%(In principle, the
843%second argument of \function{pmesh\_to\_domain\_instance} can be any
844%subclass of \class{Domain}, but for applications involving the
845%shallow-water wave equation, the second argument of
846%\function{pmesh\_to\_domain\_instance} can always be set simply to
847%\class{Domain}.)
848
849The following statements specify a basename and data directory, and
850identify quantities to be stored. For the first two, values are
851taken from \file{project.py}.
852
853{\small \begin{verbatim}
854    domain.set_name(project.basename)
855    domain.set_datadir(project.outputdir)
856    domain.set_quantities_to_be_stored(['stage', 'xmomentum',
857        'ymomentum'])
858\end{verbatim}}
859
860
861\subsection{Initial Conditions}
862Quantities for \file{runsydney.py} are set
863using similar methods to those in \file{runup.py}. However,
864in this case, many of the values are read from the auxiliary file
865\file{project.py} or, in the case of \code{elevation}, from an
866ancillary points file.
867
868
869
870\subsubsection{Stage}
871
872For the scenario we are modelling in this case, we use a callable
873object \code{tsunami\_source}, assigned by means of a function
874\function{slump\_tsunami}. This is similar to how we set elevation in
875\file{runup.py} using a function---however, in this case the
876function is both more complex and more interesting.
877
878The function returns the water displacement for all \code{x} and
879\code{y} in the domain. The water displacement is a double Gaussian
880function that depends on the characteristics of the slump (length,
881thickness, slope, etc), its location (origin) and the depth at that
882location.
883
884
885\subsubsection{Friction}
886
887We assign the friction exactly as we did for \file{runup.py}:
888
889{\small \begin{verbatim}
890    domain.set_quantity('friction', 0.0)
891\end{verbatim}}
892
893
894\subsubsection{Elevation}
895
896The elevation is specified by reading data from a file:
897
898{\small \begin{verbatim}
899    domain.set_quantity('elevation',
900                        filename = project.combineddemname + '.pts',
901                        use_cache = True,
902                        verbose = True)
903\end{verbatim}}
904
905However, before this step can be executed, some preliminary steps
906are needed to prepare the file from which the data is taken. Two
907source files are used for this data---their names are specified in
908the file \file{project.py}, in the variables \code{coarsedemname}
909and \code{finedemname}. They contain `coarse' and `fine' data,
910respectively---that is, data sampled at widely spaced points over a
911large region and data sampled at closely spaced points over a
912smaller subregion. The data in these files is combined through the
913statement
914
915{\small \begin{verbatim}
916combine_rectangular_points_files(project.finedemname + '.pts',
917                                 project.coarsedemname + '.pts',
918                                 project.combineddemname + '.pts')
919\end{verbatim}}
920
921The effect of this is simply to combine the datasets by eliminating
922any coarse data associated with points inside the smaller region
923common to both datasets. The name to be assigned to the resulting
924dataset is also derived from the name stored in the variable
925\code{combinedname} in the file \file{project.py}.
926
927\subsection{Boundary Conditions}\index{boundary conditions}
928
929Setting boundaries follows a similar pattern to the one used for
930\file{runup.py}, except that in this case we need to associate a
931boundary type with each of the
932boundary tag names introduced when we established the mesh. In place of the four
933boundary types introduced for \file{runup.py}, we use the reflective
934boundary for each of the
935eight tagged segments defined by \code{create_mesh_from_regions}:
936
937{\small \begin{verbatim}
938Bd = Dirichlet_boundary([0.0,0.0,0.0])
939domain.set_boundary( {'oceannorth': Bd,
940                      'onshorenorth1': Bd,
941                      'onshorenorthwest1': Bd,
942                      'onshorenorthwest2': Bd,
943                      'onshorenorth2': Bd,
944                      'onshorewest': Bd,
945                      'onshoresouth': Bd,
946                      'oceansouth': Bd,
947                      'oceaneast': Bd} )
948\end{verbatim}}
949
950\subsection{Evolution}
951
952With the basics established, the running of the `evolve' step is
953very similar to the corresponding step in \file{runup.py}. Here,
954the simulation is run for five hours (18000 seconds) with
955the output stored every two minutes (120 seconds).
956
957{\small \begin{verbatim}
958    import time t0 = time.time()
959
960    for t in domain.evolve(yieldstep = 120, duration = 18000):
961        print domain.timestepping_statistics()
962        print domain.boundary_statistics(tags = 'bottom')
963
964    print 'That took %.2f seconds' %(time.time()
965\end{verbatim}}
966
967%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
968%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
969
970\chapter{\anuga Public Interface}
971\label{ch:interface}
972
973This chapter gives an overview of the features of \anuga available
974to the user at the public interface. These are grouped under the
975following headings, which correspond to the outline of the examples
976described in Chapter \ref{ch:getstarted}:
977
978\begin{itemize}
979    \item Establishing the Mesh
980    \item Initialising the Domain
981    \item Specifying the Quantities
982    \item Initial Conditions
983    \item Boundary Conditions
984    \item Forcing Functions
985    \item Evolution
986\end{itemize}
987
988The listings are intended merely to give the reader an idea of what
989each feature is, where to find it and how it can be used---they do
990not give full specifications; for these the reader
991may consult the code. The code for every function or class contains
992a documentation string, or `docstring', that specifies the precise
993syntax for its use. This appears immediately after the line
994introducing the code, between two sets of triple quotes.
995
996Each listing also describes the location of the module in which
997the code for the feature being described can be found. All modules
998are in the folder \file{inundation} or one of its subfolders, and the
999location of each module is described relative to \file{inundation}. Rather
1000than using pathnames, whose syntax depends on the operating system,
1001we use the format adopted for importing the function or class for
1002use in Python code. For example, suppose we wish to specify that the
1003function \function{create\_mesh\_from\_regions} is in a module called
1004\module{mesh\_interface} in a subfolder of \module{inundation} called
1005\code{pmesh}. In Linux or Unix syntax, the pathname of the file
1006containing the function, relative to \file{inundation}, would be
1007
1008\begin{center}
1009%    \code{pmesh/mesh\_interface.py}
1010    \code{pmesh}$\slash$\code{mesh\_interface.py}
1011\end{center}
1012
1013while in Windows syntax it would be
1014
1015\begin{center}
1016    \code{pmesh}$\backslash$\code{mesh\_interface.py}
1017\end{center}
1018
1019Rather than using either of these forms, in this chapter we specify
1020the location simply as \code{pmesh.mesh\_interface}, in keeping with
1021the usage in the Python statement for importing the function,
1022namely:
1023\begin{center}
1024    \code{from pmesh.mesh\_interface import create\_mesh\_from\_regions}
1025\end{center}
1026
1027Each listing details the full set of parameters for the class or
1028function; however, the description is generally limited to the most
1029important parameters and the reader is again referred to the code
1030for more details.
1031
1032The following parameters are common to many functions and classes
1033and are omitted from the descriptions given below:
1034
1035%\begin{center}
1036\begin{tabular}{ll}  %\hline
1037%\textbf{Name } & \textbf{Description}\\
1038%\hline
1039\emph{use\_cache} & Specifies whether caching is to be used for improved performance. See Section \ref{sec:caching} for details on the underlying caching functionality\\
1040\emph{verbose} & If \code{True}, provides detailed terminal output
1041to the user\\  % \hline
1042\end{tabular}
1043%\end{center}
1044
1045\section{Mesh Generation}
1046
1047Before discussing the part of the interface relating to mesh
1048generation, we begin with a description of a simple example of a
1049mesh and use it to describe how mesh data is stored.
1050
1051\label{sec:meshexample} Figure \ref{fig:simplemesh} represents a
1052very simple mesh comprising just 11 points and 10 triangles.
1053
1054
1055\begin{figure}[h]
1056  \begin{center}
1057    \includegraphics[width=90mm, height=90mm]{triangularmesh.eps}
1058  \end{center}
1059
1060  \caption{A simple mesh}
1061  \label{fig:simplemesh}
1062\end{figure}
1063
1064
1065The variables \code{points}, \code{vertices} and \code{boundary}
1066represent the data displayed in Figure \ref{fig:simplemesh} as
1067follows. The list \code{points} stores the coordinates of the
1068points, and may be displayed schematically as in Table
1069\ref{tab:points}.
1070
1071
1072\begin{table}
1073  \begin{center}
1074    \begin{tabular}[t]{|c|cc|} \hline
1075      index & \code{x} & \code{y}\\  \hline
1076      0 & 1 & 1\\
1077      1 & 4 & 2\\
1078      2 & 8 & 1\\
1079      3 & 1 & 3\\
1080      4 & 5 & 5\\
1081      5 & 8 & 6\\
1082      6 & 11 & 5\\
1083      7 & 3 & 6\\
1084      8 & 1 & 8\\
1085      9 & 4 & 9\\
1086      10 & 10 & 7\\  \hline
1087    \end{tabular}
1088  \end{center}
1089
1090  \caption{Point coordinates for mesh in
1091    Figure \protect \ref{fig:simplemesh}}
1092  \label{tab:points}
1093\end{table}
1094
1095The list \code{vertices} specifies the triangles that make up the
1096mesh. It does this by specifying, for each triangle, the indices
1097(the numbers shown in the first column above) that correspond to the
1098three points at its vertices, taken in an anti-clockwise order
1099around the triangle. Thus, in the example shown in Figure
1100\ref{fig:simplemesh}, the variable \code{vertices} contains the
1101entries shown in Table \ref{tab:vertices}. The starting point is
1102arbitrary so triangle $(0,1,3)$ is considered the same as $(1,3,0)$
1103and $(3,0,1)$.
1104
1105
1106\begin{table}
1107  \begin{center}
1108    \begin{tabular}{|c|ccc|} \hline
1109      index & \multicolumn{3}{c|}{\code{vertices}}\\ \hline
1110      0 & 0 & 1 & 3\\
1111      1 & 1 & 2 & 4\\
1112      2 & 2 & 5 & 4\\
1113      3 & 2 & 6 & 5\\
1114      4 & 4 & 5 & 9\\
1115      5 & 4 & 9 & 7\\
1116      6 & 3 & 4 & 7\\
1117      7 & 7 & 9 & 8\\
1118      8 & 1 & 4 & 3\\
1119      9 & 5 & 10 & 9\\  \hline
1120    \end{tabular}
1121  \end{center}
1122
1123  \caption{Vertices for mesh in Figure \protect \ref{fig:simplemesh}}
1124  \label{tab:vertices}
1125\end{table}
1126
1127Finally, the variable \code{boundary} identifies the boundary
1128triangles and associates a tag with each.
1129
1130\refmodindex[pmesh.meshinterface]{pmesh.mesh\_interface}\label{sec:meshgeneration}
1131
1132\begin{funcdesc}  {create\_mesh\_from\_regions}{bounding_polygon,
1133                             boundary_tags,
1134                             maximum_triangle_area,
1135                             filename=None,
1136                             interior_regions=None,
1137                             poly_geo_reference=None,
1138                             mesh_geo_reference=None,
1139                             minimum_triangle_angle=28.0}
1140Module: \module{pmesh.mesh\_interface}
1141
1142This function allows a user to initiate the automatic creation of a
1143mesh inside a specified polygon (input \code{bounding_polygon}).
1144Among the parameters that can be set are the \emph{resolution}
1145(maximal area for any triangle in the mesh) and the minimal angle
1146allowable in any triangle. The user can specify a number of internal
1147polygons within each of which a separate mesh is to be created,
1148generally with a smaller resolution. \code{interior_regions} is a
1149paired list containing the interior polygon and its resolution.
1150Additionally, the user specifies a list of boundary tags, one for
1151each edge of the bounding polygon.
1152
1153\textbf{WARNING}. Note that the dictionary structure used for the
1154parameter \code{boundary\_tags} is different from that used for the
1155variable \code{boundary} that occurs in the specification of a mesh.
1156In the case of \code{boundary}, the tags are the \emph{values} of
1157the dictionary, whereas in the case of \code{boundary_tags}, the
1158tags are the \emph{keys} and the \emph{value} corresponding to a
1159particular tag is a list of numbers identifying boundary edges
1160labelled with that tag. Because of this, it is theoretically
1161possible to assign the same edge to more than one tag. However, an
1162attempt to do this will cause an error.
1163\end{funcdesc}
1164
1165
1166
1167
1168\begin{classdesc}  {Mesh}{userSegments=None,
1169                 userVertices=None,
1170                 holes=None,
1171                 regions=None,
1172                 geo_reference=None}
1173Module: \module{pmesh.mesh}
1174
1175A class used to build a mesh outline and generate a two-dimensional
1176triangular mesh. The mesh outline is used to describe features on the
1177mesh, such as the mesh boundary. Many of this classes methods are used
1178to build a mesh outline, such as \code{add\_vertices} and
1179\code{add\_region\_from\_polygon}.
1180
1181\end{classdesc}
1182
1183
1184\subsection{Key Methods of Class Mesh}
1185
1186
1187\begin{methoddesc} {add\_hole}{x,y, geo_reference=None}
1188Module: \module{pmesh.mesh},  Class: \class{Mesh}
1189
1190This method is used to build the mesh outline.  It defines a hole,
1191when the boundary of the hole has already been defined, by selecting a
1192point within the boundary.
1193
1194\end{methoddesc}
1195
1196
1197\begin{methoddesc}  {add\_hole\_from\_polygon}{self, polygon, tags=None,
1198    geo_reference=None}
1199Module: \module{pmesh.mesh},  Class: \class{Mesh}
1200
1201This method is used to add a `hole' within a region ---that is, to
1202define a interior region where the triangular mesh will not be
1203generated---to a \class{Mesh} instance. The region boundary is described by
1204the polygon passed in.  Additionally, the user specifies a list of
1205boundary tags, one for each edge of the bounding polygon.
1206\end{methoddesc}
1207
1208
1209\begin{methoddesc} {add\_region}{x,y, geo_reference=None}
1210Module: \module{pmesh.mesh},  Class: \class{Mesh}
1211
1212This method is used to build the mesh outline.  It defines a region,
1213when the boundary of the region has already been defined, by selecting
1214a point within the boundary.  A region instance is returned.  This can
1215be used to set the resolution.
1216
1217\end{methoddesc}
1218
1219\begin{methoddesc}  {add\_region\_from\_polygon}{self, polygon, tags=None,
1220                                max_triangle_area=None, geo_reference=None}
1221Module: \module{pmesh.mesh},  Class: \class{Mesh}
1222
1223This method is used to build the mesh outline.  It adds a region to a
1224\class{Mesh} instance.  Regions are commonly used to describe an area
1225with an increased density of triangles, by setting
1226\code{max_triangle_area}.  The
1227region boundary is described by the input \code{polygon}.  Additionally, the
1228user specifies a list of boundary tags, one for each edge of the
1229bounding polygon.
1230
1231\end{methoddesc}
1232
1233
1234
1235
1236
1237\begin{methoddesc} {add\_vertices}{point_data}
1238
1239Add user vertices. The point_data can be a list of (x,y) values, a numeric
1240array or a geospatial_data instance.
1241       
1242\end{methoddesc}
1243
1244\begin{methoddesc}  {export\_mesh_file}{self,ofile}
1245Module: \module{pmesh.mesh},  Class: \class{Mesh}
1246
1247This method is used to save the mesh to a file. \code{ofile} is the
1248name of the mesh file to be written, including the extension.  Use
1249the extension \code{.msh} for the file to be in NetCDF format and
1250\code{.tsh} for the file to be ASCII format.
1251\end{methoddesc}
1252
1253\begin{methoddesc}  {generate\_mesh}{self,
1254                      maximum_triangle_area=None,
1255                      minimum_triangle_angle=28.0,
1256                      verbose=False}
1257Module: \module{pmesh.mesh},  Class: \class{Mesh}
1258
1259This method is used to generate the triangular mesh.  The  maximal
1260area of any triangle in the mesh can be specified, which is used to
1261control the triangle density, along with the
1262minimum angle in any triangle.
1263\end{methoddesc}
1264
1265
1266
1267\begin{methoddesc}  {import_ungenerate_file}{self,ofile, tag=None}
1268Module: \module{pmesh.mesh},  Class: \class{Mesh}
1269
1270This method is used to import a polygon file in the ungenerate
1271format, which is used by arcGIS. The polygons from the file are converted to
1272vertices and segments. \code{ofile} is the name of the polygon file.
1273\code{tag} is the tag given to all the polygon's segments.
1274
1275This function can be used to import building footprints.
1276\end{methoddesc}
1277
1278%%%%%%
1279\section{Initialising the Domain}
1280
1281%Include description of the class Domain and the module domain.
1282
1283%FIXME (Ole): This is also defined in a later chapter
1284%\declaremodule{standard}{pyvolution.domain}
1285
1286\begin{classdesc} {Domain} {source=None,
1287                 triangles=None,
1288                 boundary=None,
1289                 conserved_quantities=None,
1290                 other_quantities=None,
1291                 tagged_elements=None,
1292                 geo_reference=None,
1293                 use_inscribed_circle=False,
1294                 mesh_filename=None,
1295                 use_cache=False,
1296                 verbose=False,
1297                 full_send_dict=None,
1298                 ghost_recv_dict=None,
1299                 processor=0,
1300                 numproc=1}
1301Module: \refmodule{pyvolution.domain}
1302
1303This class is used to create an instance of a data structure used to
1304store and manipulate data associated with a mesh. The mesh is
1305specified either by assigning the name of a mesh file to
1306\code{source} or by
1307\end{classdesc}
1308
1309\begin{funcdesc}  {pmesh\_to\_domain\_instance}{file_name, DomainClass, use_cache = False, verbose = False}
1310Module: \module{pyvolution.pmesh2domain}
1311
1312Once the initial mesh file has been created, this function is
1313applied to convert it to a domain object---that is, to a member of
1314the special Python class \class{Domain} (or a subclass of
1315\class{Domain}), which provides access to properties and methods
1316that allow quantities to be set and other operations to be carried
1317out.
1318
1319\code{file\_name} is the name of the mesh file to be converted,
1320including the extension. \code{DomainClass} is the class to be
1321returned, which must be a subclass of \class{Domain} having the same
1322interface as \class{Domain}---in practice, it can usually be set
1323simply to \class{Domain}.
1324
1325This is now superseded by Domain(mesh_filename).
1326\end{funcdesc}
1327
1328
1329\subsection{Key Methods of Domain}
1330
1331\begin{methoddesc} {set\_name}{name}
1332    Module: \refmodule{pyvolution.domain}, page \pageref{mod:pyvolution.domain}  %\code{pyvolution.domain}
1333
1334    Assigns the name \code{name} to the domain.
1335\end{methoddesc}
1336
1337\begin{methoddesc} {get\_name}{}
1338    Module: \module{pyvolution.domain}
1339
1340    Returns the name assigned to the domain by \code{set\_name}. If no name has been
1341    assigned, returns \code{`domain'}.
1342\end{methoddesc}
1343
1344\begin{methoddesc} {set\_datadir}{name}
1345    Module: \module{pyvolution.domain}
1346
1347    Specifies the directory used for SWW files, assigning it to the pathname \code{name}. The default value, before
1348    \code{set\_datadir} has been run, is the value \code{default\_datadir}
1349    specified in \code{config.py}.
1350
1351    Since different operating systems use different formats for specifying pathnames,
1352    it is necessary to specify path separators using the Python code \code{os.sep}, rather than
1353    the operating-specific ones such as `$\slash$' or `$\backslash$'.
1354    For this to work you will need to include the statement \code{import os}
1355    in your code, before the first appearance of \code{set\_datadir}.
1356
1357    For example, to set the data directory to a subdirectory
1358    \code{data} of the directory \code{project}, you could use
1359    the statements:
1360
1361    {\small \begin{verbatim}
1362        import os
1363        domain.set_datadir{'project' + os.sep + 'data'}
1364    \end{verbatim}}
1365\end{methoddesc}
1366
1367\begin{methoddesc} {get\_datadir}{}
1368    Module: \module{pyvolution.domain}
1369
1370    Returns the data directory set by \code{set\_datadir} or,
1371    if \code{set\_datadir} has not
1372    been run, returns the value \code{default\_datadir} specified in
1373    \code{config.py}.
1374\end{methoddesc}
1375
1376\begin{methoddesc} {set\_time}{time=0.0}
1377    Module: \module{pyvolution.domain}
1378
1379    Sets the initial time, in seconds, for the simulation. The
1380    default is 0.0.
1381\end{methoddesc}
1382
1383\begin{methoddesc} {set\_default\_order}{n}
1384    Sets the default (spatial) order to the value specified by
1385    \code{n}, which must be either 1 or 2. (Assigning any other value
1386    to \code{n} will cause an error.)
1387\end{methoddesc}
1388
1389
1390%%%%%%
1391\section{Initial Conditions}
1392\label{sec:Initial Conditions}
1393In standard usage of partial differential equations, initial conditions
1394refers to the values associated to the system variables (the conserved
1395quantities here) for \code{time = 0}. In setting up a scenario script
1396as described in Sections \ref{sec:simpleexample} and \ref{sec:realdataexample},
1397\code{set_quantity} is used to define the initial conditions of variables
1398other than the conserved quantities, such as friction. Here, we use the terminology
1399of initial conditions to refer to initial values for variables which need
1400prescription to solve the shallow water wave equation. Further, it must be noted
1401that \code{set_quantity} does not necessarily have to be used in the initial
1402condition setting; it can be used at any time throughout the simulation.
1403
1404\begin{methoddesc}{set\_quantity}{name,
1405    numeric = None,
1406    quantity = None,
1407    function = None,
1408    geospatial_data = None,
1409    filename = None,
1410    attribute_name = None,
1411    alpha = None,
1412    location = 'vertices',
1413    indices = None,
1414    verbose = False,
1415    use_cache = False}
1416  Module: \module{pyvolution.domain}
1417  (see also \module{pyvolution.quantity.set\_values})
1418
1419This function is used to assign values to individual quantities for a
1420domain. It is very flexible and can be used with many data types: a
1421statement of the form \code{domain.set\_quantity(name, x)} can be used
1422to define a quantity having the name \code{name}, where the other
1423argument \code{x} can be any of the following:
1424
1425\begin{itemize}
1426\item a number, in which case all vertices in the mesh gets that for
1427the quantity in question.
1428\item a list of numbers or a Numeric array ordered the same way as the mesh vertices.
1429\item a function (e.g.\ see the samples introduced in Chapter 2)
1430\item an expression composed of other quantities and numbers, arrays, lists (for
1431example, a linear combination of quantities, such as
1432\code{domain.set\_quantity('stage','elevation'+x))}
1433\item the name of a file from which the data can be read. In this case, the optional argument attribute\_name will select which attribute to use from the file. If left out, set\_quantity will pick one. This is useful in cases where there is only one attribute.
1434\item a geospatial dataset (See ?????). Optional argument attribute\_name applies here as with files.
1435\end{itemize}
1436
1437
1438Exactly one of the arguments
1439  numeric, quantity, function, points, filename
1440must be present.
1441
1442
1443Set quantity will look at the type of the second argument (\code{numeric}) and
1444determine what action to take.
1445
1446Values can also be set using the appropriate keyword arguments.
1447If x is a function, for example, \code{domain.set\_quantity(name, x)}, \code{domain.set\_quantity(name, numeric=x)}, and \code{domain.set\_quantity(name, function=x)}
1448are all equivalent.
1449
1450
1451Other optional arguments are
1452\begin{itemize}
1453\item \code{indices} which is a list of ids of triangles to which set\_quantity should apply its assignment of values.
1454\item \code{location} determines which part of the triangles to assign to. Options are 'vertices' (default), 'edges', and 'centroids'.
1455\end{itemize}
1456
1457
1458\end{methoddesc}
1459
1460
1461
1462
1463
1464
1465
1466%%%
1467\anuga provides a number of predefined initial conditions to be used
1468with \code{set\_quantity}.
1469
1470\begin{funcdesc}{slump_tsunami}{length, depth, slope, width=None, thickness=None,
1471                x0=0.0, y0=0.0, alpha=0.0,
1472                gravity=9.8, gamma=1.85,
1473                massco=1, dragco=1, frictionco=0, psi=0,
1474                dx=None, kappa=3.0, kappad=0.8, zsmall=0.01,
1475                domain=None,
1476                verbose=False}
1477Module: \module{pyvolution.smf}
1478
1479This function returns a callable object representing an initial water
1480displacement generated by a submarine sediment failure. These failures can take the form of
1481a submarine slump or slide. In the case of a slide, use \code{slide_tsunami} instead.
1482
1483The arguments include as a minimum, the slump or slide length, the water depth to the centre of sediment
1484mass, and the bathymetric slope. Other slump or slide parameters can be included if they are known.
1485\end{funcdesc}
1486
1487
1488%%%
1489\begin{funcdesc}{file\_function}{filename,
1490    domain = None,
1491    quantities = None,
1492    interpolation_points = None,
1493    verbose = False,
1494    use_cache = False}
1495Module: \module{pyvolution.util}
1496
1497Reads the time history of spatial data for
1498specified interpolation points from a NetCDF file (\code{filename})
1499and returns
1500a callable object. \code{filename} could be a \code{sww} file.
1501Returns interpolated values based on the input
1502file using the underlying \code{interpolation\_function}.
1503
1504\code{quantities} is either the name of a single quantity to be
1505interpolated or a list of such quantity names. In the second case, the resulting
1506function will return a tuple of values---one for each quantity.
1507
1508\code{interpolation\_points} is a list of absolute UTM coordinates
1509for points at which values are sought.
1510
1511The model time stored within the file function can be accessed using
1512the method \code{f.get\_time()}
1513\end{funcdesc}
1514
1515%%%
1516\begin{classdesc}{Interpolation\_function}{self,
1517    time,
1518    quantities,
1519    quantity_names = None,
1520    vertex_coordinates = None,
1521    triangles = None,
1522    interpolation_points = None,
1523    verbose = False}
1524Module: \module{pyvolution.least\_squares}
1525
1526Given a time series (i.e. a series of values associated with
1527different times), whose values are either just numbers or a set of
1528numbers defined at the vertices of a triangular mesh (such as those
1529stored in SWW files), \code{Interpolation\_function} is used to
1530create a callable object that interpolates a value for an arbitrary
1531time \code{t} within the model limits and possibly a point \code{(x,
1532y)} within a mesh region.
1533
1534The actual time series at which data is available is specified by
1535means of an array \code{time} of monotonically increasing times. The
1536quantities containing the values to be interpolated are specified in
1537an array---or dictionary of arrays (used in conjunction with the
1538optional argument \code{quantity\_names}) --- called
1539\code{quantities}. The optional arguments \code{vertex\_coordinates}
1540and \code{triangles} represent the spatial mesh associated with the
1541quantity arrays. If omitted the function created by
1542\code{Interpolation\_function} will be a function of \code{t} only.
1543
1544Since, in practice, values need to be computed at specified points,
1545the syntax allows the user to specify, once and for all, a list
1546\code{interpolation\_points} of points at which values are required.
1547In this case, the function may be called using the form \code{f(t,
1548id)}, where \code{id} is an index for the list
1549\code{interpolation\_points}.
1550
1551\end{classdesc}
1552
1553%%%
1554%\begin{funcdesc}{set\_region}{functions}
1555%[Low priority. Will be merged into set\_quantity]
1556
1557%Module:\module{pyvolution.domain}
1558%\end{funcdesc}
1559
1560
1561
1562%%%%%%
1563\section{Boundary Conditions}\index{boundary conditions}
1564
1565\anuga provides a large number of predefined boundary conditions,
1566represented by objects such as \code{Reflective\_boundary(domain)} and
1567\code{Dirichlet\_boundary([0.2, 0.0, 0.0])}, described in the examples
1568in Chapter 2. Alternatively, you may prefer to ``roll your own'',
1569following the method explained in Section \ref{sec:roll your own}.
1570
1571These boundary objects may be used with the function \code{set\_boundary} described below
1572to assign boundary conditions according to the tags used to label boundary segments.
1573
1574\begin{methoddesc}{set\_boundary}{boundary_map}
1575Module: \module{pyvolution.domain}
1576
1577This function allows you to assign a boundary object (corresponding to a
1578pre-defined or user-specified boundary condition) to every boundary segment that
1579has been assigned a particular tag.
1580
1581This is done by specifying a dictionary \code{boundary\_map}, whose values are the boundary objects
1582and whose keys are the symbolic tags.
1583
1584\end{methoddesc}
1585
1586\begin{methoddesc} {get\_boundary\_tags}{}
1587Module: \module{pyvolution.mesh}
1588
1589Returns a list of the available boundary tags.
1590\end{methoddesc}
1591
1592%%%
1593\subsection{Predefined boundary conditions}
1594
1595\begin{classdesc}{Reflective\_boundary}{Boundary}
1596Module: \module{pyvolution.shallow\_water}
1597
1598Reflective boundary returns same conserved quantities as those present in
1599the neighbouring volume but reflected.
1600
1601This class is specific to the shallow water equation as it works with the
1602momentum quantities assumed to be the second and third conserved quantities.
1603\end{classdesc}
1604
1605%%%
1606\begin{classdesc}{Transmissive\_boundary}{domain = None}
1607Module: \module{pyvolution.generic\_boundary\_conditions}
1608
1609A transmissive boundary returns the same conserved quantities as
1610those present in the neighbouring volume.
1611
1612The underlying domain must be specified when the boundary is instantiated.
1613\end{classdesc}
1614
1615%%%
1616\begin{classdesc}{Dirichlet\_boundary}{conserved_quantities=None}
1617Module: \module{pyvolution.generic\_boundary\_conditions}
1618
1619A Dirichlet boundary returns constant values for each of conserved
1620quantities. In the example of \code{Dirichlet\_boundary([0.2, 0.0, 0.0])},
1621the \code{stage} value at the boundary is 0.2 and the \code{xmomentum} and
1622\code{ymomentum} at the boundary are set to 0.0. The list must contain
1623a value for each conserved quantity.
1624\end{classdesc}
1625
1626%%%
1627\begin{classdesc}{Time\_boundary}{domain = None, f = None}
1628Module: \module{pyvolution.generic\_boundary\_conditions}
1629
1630A time-dependent boundary returns values for the conserved
1631quantities as a function \code{f(t)} of time. The user must specify
1632the domain to get access to the model time.
1633\end{classdesc}
1634
1635%%%
1636\begin{classdesc}{File\_boundary}{Boundary}
1637Module: \module{pyvolution.generic\_boundary\_conditions}
1638
1639This method may be used if the user wishes to apply a SWW file or
1640a time series file to a boundary segment or segments.
1641The boundary values are obtained from a file and interpolated to the
1642appropriate segments for each conserved quantity.
1643\end{classdesc}
1644
1645
1646
1647%%%
1648\begin{classdesc}{Transmissive\_Momentum\_Set\_Stage\_boundary}{Boundary}
1649Module: \module{pyvolution.shallow\_water}
1650
1651This boundary returns same momentum conserved quantities as
1652those present in its neighbour volume but sets stage as in a Time\_boundary.
1653The underlying domain must be specified when boundary is instantiated
1654
1655This type of boundary is useful when stage is known at the boundary as a
1656function of time, but momenta (or speeds) aren't.
1657
1658This class is specific to the shallow water equation as it works with the
1659momentum quantities assumed to be the second and third conserved quantities.
1660\end{classdesc}
1661
1662
1663
1664\subsection{User-defined boundary conditions}
1665\label{sec:roll your own}
1666[How to roll your own] TBA
1667
1668
1669
1670
1671
1672\section{Forcing Functions}
1673
1674\anuga provides a number of predefined forcing functions to be used with .....
1675
1676%\begin{itemize}
1677
1678
1679%  \item \indexedcode{}
1680%  [function, arguments]
1681
1682%  \item \indexedcode{}
1683
1684%\end{itemize}
1685
1686
1687
1688\section{Evolution}\index{evolution}
1689
1690  \begin{methoddesc}{evolve}{yieldstep = None, finaltime = None, duration = None, skip_initial_step = False}
1691
1692  Module: \module{pyvolution.domain}
1693
1694  This function (a method of \class{domain}) is invoked once all the
1695  preliminaries have been completed, and causes the model to progress
1696  through successive steps in its evolution, storing results and
1697  outputting statistics whenever a user-specified period
1698  \code{yieldstep} is completed (generally during this period the
1699  model will evolve through several steps internally
1700  as the method forces the water speed to be calculated
1701  on successive new cells). The user
1702  specifies the total time period over which the evolution is to take
1703  place, by specifying values (in seconds) for either \code{duration}
1704  or \code{finaltime}, as well as the interval in seconds after which
1705  results are to be stored and statistics output.
1706
1707  You can include \method{evolve} in a statement of the type:
1708
1709  {\small \begin{verbatim}
1710      for t in domain.evolve(yieldstep, finaltime):
1711          <Do something with domain and t>
1712  \end{verbatim}}
1713
1714  \end{methoddesc}
1715
1716
1717
1718\subsection{Diagnostics}
1719
1720
1721  \begin{funcdesc}{statistics}{}
1722  Module: \module{pyvolution.domain}
1723
1724  \end{funcdesc}
1725
1726  \begin{funcdesc}{timestepping\_statistics}{}
1727  Module: \module{pyvolution.domain}
1728
1729  Returns a string of the following type for each
1730  timestep:
1731
1732  \code{Time = 0.9000, delta t in [0.00598964, 0.01177388], steps=12
1733  (12)}
1734
1735  Here the numbers in \code{steps=12 (12)} indicate the number of steps taken and
1736  the number of first-order steps, respectively.
1737  \end{funcdesc}
1738
1739
1740  \begin{funcdesc}{boundary\_statistics}{quantities = None, tags = None}
1741  Module: \module{pyvolution.domain}
1742
1743  Returns a string of the following type when \code{quantities = 'stage'} and \code{tags = ['top', 'bottom']}:
1744
1745  {\small \begin{verbatim}
1746 Boundary values at time 0.5000:
1747    top:
1748        stage in [ -0.25821218,  -0.02499998]
1749    bottom:
1750        stage in [ -0.27098821,  -0.02499974]
1751  \end{verbatim}}
1752
1753  \end{funcdesc}
1754
1755
1756  \begin{funcdesc}{get\_quantity}{name, location='vertices', indices = None}
1757  Module: \module{pyvolution.domain}
1758  Allow access to individual quantities and their methods
1759
1760  \end{funcdesc}
1761
1762
1763  \begin{funcdesc}{get\_values}{location='vertices', indices = None}
1764  Module: \module{pyvolution.quantity}
1765
1766  Extract values for quantity as an array
1767
1768  \end{funcdesc}
1769
1770
1771  \begin{funcdesc}{get\_integral}{}
1772  Module: \module{pyvolution.quantity}
1773
1774  Return computed integral over entire domain for this quantity
1775
1776  \end{funcdesc}
1777
1778
1779\section{Other}
1780
1781  \begin{funcdesc}{domain.create\_quantity\_from\_expression}{???}
1782
1783  Handy for creating derived quantities on-the-fly.
1784  See \file{Analytical\_solution\_circular\_hydraulic\_jump.py} for an example of use.
1785  \end{funcdesc}
1786
1787
1788
1789
1790
1791%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1792%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1793
1794\chapter{\anuga System Architecture}
1795
1796From pyvolution/documentation
1797
1798\section{File Formats}
1799\label{sec:file formats}
1800
1801\anuga makes use of a number of different file formats. The
1802following table lists all these formats, which are described in more
1803detail in the paragraphs below.
1804
1805\bigskip
1806
1807\begin{center}
1808
1809\begin{tabular}{|ll|}  \hline
1810
1811\textbf{Extension} & \textbf{Description} \\
1812\hline\hline
1813
1814\code{.sww} & NetCDF format for storing model output
1815\code{f(t,x,y)}\\
1816
1817\code{.tms} & NetCDF format for storing time series \code{f(t)}\\
1818
1819\code{.xya} & ASCII format for storing arbitrary points and
1820associated attributes\\
1821
1822\code{.pts} & NetCDF format for storing arbitrary points and
1823associated attributes\\
1824
1825\code{.asc} & ASCII format of regular DEMs as output from ArcView\\
1826
1827\code{.prj} & Associated ArcView file giving more metadata for
1828\code{.asc} format\\
1829
1830\code{.ers} & ERMapper header format of regular DEMs for ArcView\\
1831
1832\code{.dem} & NetCDF representation of regular DEM data\\
1833
1834\code{.tsh} & ASCII format for storing meshes and associated
1835boundary and region info\\
1836
1837\code{.msh} & NetCDF format for storing meshes and associated
1838boundary and region info\\
1839
1840\code{.nc} & Native ferret NetCDF format\\
1841
1842\code{.geo} & Houdinis ASCII geometry format (?) \\  \par \hline
1843%\caption{File formats used by \anuga}
1844\end{tabular}
1845
1846
1847\end{center}
1848
1849The above table shows the file extensions used to identify the
1850formats of files. However, typically, in referring to a format we
1851capitalise the extension and omit the initial full stop---thus, we
1852refer, for example, to `SWW files' or `PRJ files'.
1853
1854\bigskip
1855
1856A typical dataflow can be described as follows:
1857
1858\subsection{Manually Created Files}
1859
1860\begin{tabular}{ll}
1861ASC, PRJ & Digital elevation models (gridded)\\
1862NC & Model outputs for use as boundary conditions (e.g. from MOST)
1863\end{tabular}
1864
1865\subsection{Automatically Created Files}
1866
1867\begin{tabular}{ll}
1868ASC, PRJ  $\rightarrow$  DEM  $\rightarrow$  PTS & Convert
1869DEMs to native \code{.pts} file\\
1870
1871NC $\rightarrow$ SWW & Convert MOST boundary files to
1872boundary \code{.sww}\\
1873
1874PTS + TSH $\rightarrow$ TSH with elevation & Least squares fit\\
1875
1876TSH $\rightarrow$ SWW & Convert TSH to \code{.sww}-viewable using
1877Swollen\\
1878
1879TSH + Boundary SWW $\rightarrow$ SWW & Simulation using
1880\code{pyvolution}\\
1881
1882Polygonal mesh outline $\rightarrow$ & TSH or MSH
1883\end{tabular}
1884
1885
1886
1887
1888\bigskip
1889
1890\subsection{SWW and TMS Formats}
1891
1892The SWW and TMS formats are both NetCDF formats, and are of key
1893importance for \anuga.
1894
1895An SWW file is used for storing \anuga output and therefore pertains
1896to a set of points and a set of times at which a model is evaluated.
1897It contains, in addition to dimension information, the following
1898variables:
1899
1900\begin{itemize}
1901    \item \code{x} and \code{y}: coordinates of the points, represented as Numeric arrays
1902    \item \code{elevation}, a Numeric array storing bed-elevations
1903    \item \code{volumes}, a list specifying the points at the vertices of each of the
1904    triangles
1905    % Refer here to the example to be provided in describing the simple example
1906    \item \code{time}, a Numeric array containing times for model
1907    evaluation
1908\end{itemize}
1909
1910
1911The contents of an SWW file may be viewed using the visualisation
1912tool \code{swollen}, which creates an on-screen geometric
1913representation. See section \ref{sec:swollen} (page
1914\pageref{sec:swollen}) in Appendix \ref{ch:supportingtools} for more
1915on \code{swollen}.
1916
1917Alternatively, there are tools, such as \code{ncdump}, that allow
1918you to convert an NetCDF file into a readable format such as the
1919Class Definition Language (CDL). The following is an excerpt from a
1920CDL representation of the output file \file{bedslope.sww} generated
1921from running the simple example \file{runup.py} of
1922Chapter \ref{ch:getstarted}:
1923
1924\verbatiminput{examples/bedslopeexcerpt.cdl}
1925
1926The SWW format is used not only for output but also serves as input
1927for functions such as \function{file\_boundary} and
1928\function{file\_function}, described in Chapter \ref{ch:interface}.
1929
1930A TMS file is used to store time series data that is independent of
1931position.
1932
1933
1934\subsection{Mesh File Formats}
1935
1936A mesh file is a file that has a specific format suited to
1937triangular meshes and their outlines. A mesh file can have one of
1938two formats: it can be either a TSH file, which is an ASCII file, or
1939an MSH file, which is a NetCDF file. A mesh file can be generated
1940from the function \function{create\_mesh\_from\_regions} (see
1941Section \ref{sec:meshgeneration}) and used to initialise a domain.
1942
1943A mesh file can define the outline of the mesh---the vertices and
1944line segments that enclose the region in which the mesh is
1945created---and the triangular mesh itself, which is specified by
1946listing the triangles and their vertices, and the segments, which
1947are those sides of the triangles that are associated with boundary
1948conditions.
1949
1950In addition, a mesh file may contain `holes' and/or `regions'. A
1951hole represents an area where no mesh is to be created, while a
1952region is a labelled area used for defining properties of a mesh,
1953such as friction values.  A hole or region is specified by a point
1954and bounded by a number of segments that enclose that point.
1955
1956A mesh file can also contain a georeference, which describes an
1957offset to be applied to $x$ and $y$ values---eg to the vertices.
1958
1959
1960\subsection{Formats for Storing Arbitrary Points and Attributes}
1961
1962An XYA file is used to store data representing arbitrary numerical
1963attributes associated with a set of points.
1964
1965The format for an XYA file is:\\
1966%\begin{verbatim}
1967
1968            first line:     \code{[attribute names]}\\
1969            other lines:  \code{x y [attributes]}\\
1970
1971            for example:\\
1972            \code{elevation, friction}\\
1973            \code{0.6, 0.7, 4.9, 0.3}\\
1974            \code{1.9, 2.8, 5, 0.3}\\
1975            \code{2.7, 2.4, 5.2, 0.3}
1976
1977        The first two columns are always implicitly assumed to be $x$, $y$ coordinates.
1978        Use the same delimiter for the attribute names and the data.
1979
1980        An XYA file can optionally end with a description of the georeference:
1981
1982            \code{\#geo reference}\\
1983            \code{56}\\
1984            \code{466600.0}\\
1985            \code{8644444.0}
1986
1987        Here the first number specifies the UTM zone (in this case zone 56)  and other numbers specify the
1988        easting and northing
1989        coordinates (in this case (466600.0, 8644444.0)) of the lower left corner.
1990
1991A PTS file is a NetCDF representation of the data held in an XYA
1992file. If the data is associated with a set of $N$ points, then the
1993data is stored using an $N \times 2$ Numeric array of float
1994variables for the points and an $N \times 1$ Numeric array for each
1995attribute.
1996
1997%\end{verbatim}
1998
1999\subsection{ArcView Formats}
2000
2001Files of the three formats ASC, PRJ and ERS are all associated with
2002data from ArcView.
2003
2004An ASC file is an ASCII representation of DEM output from ArcView.
2005It contains a header with the following format:
2006
2007\begin{tabular}{l l}
2008\code{ncols}      &   \code{753}\\
2009\code{nrows}      &   \code{766}\\
2010\code{xllcorner}  &   \code{314036.58727982}\\
2011\code{yllcorner}  & \code{6224951.2960092}\\
2012\code{cellsize}   & \code{100}\\
2013\code{NODATA_value} & \code{-9999}
2014\end{tabular}
2015
2016The remainder of the file contains the elevation data for each grid point
2017in the grid defined by the above information.
2018
2019A PRJ file is an ArcView file used in conjunction with an ASC file
2020to represent metadata for a DEM.
2021
2022
2023\subsection{DEM Format}
2024
2025A DEM file is a NetCDF representation of regular DEM data.
2026
2027
2028\subsection{Other Formats}
2029
2030
2031
2032
2033\subsection{Basic File Conversions}
2034
2035  \begin{funcdesc}{sww2dem}{basename_in, basename_out = None,
2036            quantity = None,
2037            timestep = None,
2038            reduction = None,
2039            cellsize = 10,
2040            NODATA_value = -9999,
2041            easting_min = None,
2042            easting_max = None,
2043            northing_min = None,
2044            northing_max = None,
2045            expand_search = False,
2046            verbose = False,
2047            origin = None,
2048            datum = 'WGS84',
2049            format = 'ers'}
2050  Module: \module{pyvolution.data\_manager}
2051
2052  Takes data from an SWW file \code{basename_in} and converts it to DEM format (ASC or
2053  ERS) of a desired grid size \code{cellsize} in metres.
2054  The easting and northing values are used if the user wished to clip the output
2055  file to a specified rectangular area. The \code{reduction} input refers to a function
2056  to reduce the quantities over all time step of the SWW file, example, maximum.
2057  \end{funcdesc}
2058
2059
2060  \begin{funcdesc}{dem2pts}{basename_in, basename_out=None,
2061            easting_min=None, easting_max=None,
2062            northing_min=None, northing_max=None,
2063            use_cache=False, verbose=False}
2064  Module: \module{pyvolution.data\_manager}
2065
2066  Takes DEM data (a NetCDF file representation of data from a regular Digital
2067  Elevation Model) and converts it to PTS format.
2068  \end{funcdesc}
2069
2070
2071
2072%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2073%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2074
2075\chapter{Basic \anuga Assumptions}
2076
2077(From pyvolution/documentation)
2078
2079
2080Physical model time cannot be earlier than 1 Jan 1970 00:00:00.
2081If one wished to recreate scenarios prior to that date it must be done
2082using some relative time (e.g. 0).
2083
2084
2085All spatial data relates to the WGS84 datum (or GDA94) and has been
2086projected into UTM with false easting of 500000 and false northing of
20871000000 on the southern hemisphere (0 on the northern).
2088
2089It is assumed that all computations take place within one UTM zone.
2090
2091DEMs, meshes and boundary conditions can have different origins within
2092one UTM zone. However, the computation will use that of the mesh for
2093numerical stability.
2094
2095
2096%OLD
2097%The dataflow is: (See data_manager.py and from scenarios)
2098%
2099%
2100%Simulation scenarios
2101%--------------------%
2102%%
2103%
2104%Sub directories contain scrips and derived files for each simulation.
2105%The directory ../source_data contains large source files such as
2106%DEMs provided externally as well as MOST tsunami simulations to be used
2107%as boundary conditions.
2108%
2109%Manual steps are:
2110%  Creation of DEMs from argcview (.asc + .prj)
2111%  Creation of mesh from pmesh (.tsh)
2112%  Creation of tsunami simulations from MOST (.nc)
2113%%
2114%
2115%Typical scripted steps are%
2116%
2117%  prepare_dem.py:  Convert asc and prj files supplied by arcview to
2118%                   native dem and pts formats%
2119%
2120%  prepare_pts.py: Convert netcdf output from MOST to an sww file suitable
2121%                  as boundary condition%
2122%
2123%  prepare_mesh.py: Merge DEM (pts) and mesh (tsh) using least squares
2124%                   smoothing. The outputs are tsh files with elevation data.%
2125%
2126%  run_simulation.py: Use the above together with various parameters to
2127%                     run inundation simulation.
2128
2129
2130%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2131%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2132
2133\appendix
2134
2135\chapter{Supporting Tools}
2136\label{ch:supportingtools}
2137
2138This section describes a number of supporting tools, supplied with \anuga, that offer a
2139variety of types of functionality and enhance the basic capabilities of \anuga.
2140
2141\section{caching}
2142\label{sec:caching}
2143
2144The \code{cache} function is used to provide supervised caching of function
2145results. A Python function call of the form
2146
2147      {\small \begin{verbatim}
2148      result = func(arg1,...,argn)
2149      \end{verbatim}}
2150
2151  can be replaced by
2152
2153      {\small \begin{verbatim}
2154      from caching import cache
2155      result = cache(func,(arg1,...,argn))
2156      \end{verbatim}}
2157
2158  which returns the same output but reuses cached
2159  results if the function has been computed previously in the same context.
2160  \code{result} and the arguments can be simple types, tuples, list, dictionaries or
2161  objects, but not unhashable types such as functions or open file objects.
2162  The function \code{func} may be a member function of an object or a module.
2163
2164  This type of caching is particularly useful for computationally intensive
2165  functions with few frequently used combinations of input arguments. Note that
2166  if the inputs or output are very large caching may not save time because
2167  disc access may dominate the execution time.
2168
2169  If the function definition changes after a result has been cached, this will be
2170  detected by examining the functions \code{bytecode (co\_code, co\_consts,
2171  func\_defaults, co\_argcount)} and the function will be recomputed.
2172  However, caching will not detect changes in modules used by \code{func}.
2173  In this case cache must be cleared manually.
2174
2175  Options are set by means of the function \code{set\_option(key, value)},
2176  where \code{key} is a key associated with a
2177  Python dictionary \code{options}. This dictionary stores settings such as the name of
2178  the directory used, the maximum
2179  number of cached files allowed, and so on.
2180
2181  The \code{cache} function allows the user also to specify a list of dependent files. If any of these
2182  have been changed, the function is recomputed and the results stored again.
2183
2184  %Other features include support for compression and a capability to \ldots
2185
2186
2187   \textbf{USAGE:} \nopagebreak
2188
2189    {\small \begin{verbatim}
2190    result = cache(func, args, kwargs, dependencies, cachedir, verbose,
2191                   compression, evaluate, test, return_filename)
2192    \end{verbatim}}
2193
2194
2195\section{swollen}
2196\label{sec:swollen}
2197 The output generated by \anuga may be viewed by
2198means of the visualisation tool \code{swollen}, which takes the
2199\code{SWW} file output by \anuga and creates a visual representation
2200of the data. Examples may be seen in Figures \ref{fig:runupstart}
2201and \ref{fig:bedslope2}. To view an \code{SWW} file with
2202\code{swollen} in the Windows environment, you can simply drag the
2203icon representing the file over an icon on the desktop for the
2204\code{swollen} executable file (or a shortcut to it), or set up a
2205file association to make files with the extension \code{.sww} open
2206with \code{swollen}. Alternatively, you can operate \code{swollen}
2207from the command line, in both Windows and Linux environments.
2208
2209On successful operation, you will see an interactive moving-picture
2210display. You can use keys and the mouse to slow down, speed up or
2211stop the display, change the viewing position or carry out a number
2212of other simple operations. Help is also displayed when you press
2213the \code{h} key.
2214
2215The main keys operating the interactive screen are:\\
2216
2217\begin{center}
2218\begin{tabular}{|ll|}   \hline
2219
2220\code{w} & toggle wireframe \\
2221
2222space bar & start/stop\\
2223
2224up/down arrows & increase/decrease speed\\
2225
2226left/right arrows & direction in time \emph{(when running)}\\
2227& step through simulation \emph{(when stopped)}\\
2228
2229left mouse button & rotate\\
2230
2231middle mouse button & pan\\
2232
2233right mouse button & zoom\\  \hline
2234
2235\end{tabular}
2236\end{center}
2237
2238\vfill
2239
2240The following table describes how to operate swollen from the command line:
2241
2242Usage: \code{swollen [options] swwfile \ldots}\\  \nopagebreak
2243Options:\\  \nopagebreak
2244\begin{tabular}{ll}
2245  \code{--display <type>} & \code{MONITOR | POWERWALL | REALITY\_CENTER |}\\
2246                                    & \code{HEAD\_MOUNTED\_DISPLAY}\\
2247  \code{--rgba} & Request a RGBA colour buffer visual\\
2248  \code{--stencil} & Request a stencil buffer visual\\
2249  \code{--stereo} & Use default stereo mode which is \code{ANAGLYPHIC} if not \\
2250                                    & overridden by environmental variable\\
2251  \code{--stereo <mode>} & \code{ANAGLYPHIC | QUAD\_BUFFER | HORIZONTAL\_SPLIT |}\\
2252                                    & \code{VERTICAL\_SPLIT | LEFT\_EYE | RIGHT\_EYE |}\\
2253                                     & \code{ON | OFF} \\
2254  \code{-alphamax <float 0-1>} & Maximum transparency clamp value\\
2255  \code{-alphamin <float 0-1>} & Transparency value at \code{hmin}\\
2256\end{tabular}
2257
2258\begin{tabular}{ll}
2259  \code{-cullangle <float angle 0-90>} & Cull triangles steeper than this value\\
2260  \code{-help} & Display this information\\
2261  \code{-hmax <float>} & Height above which transparency is set to
2262                                     \code{alphamax}\\
2263\end{tabular}
2264
2265\begin{tabular}{ll}
2266
2267  \code{-hmin <float>} & Height below which transparency is set to
2268                                     zero\\
2269\end{tabular}
2270
2271\begin{tabular}{ll}
2272  \code{-lightpos <float>,<float>,<float>} & $x,y,z$ of bedslope directional light ($z$ is
2273                                     up, default is overhead)\\
2274\end{tabular}
2275
2276\begin{tabular}{ll}
2277  \code{-loop}  & Repeated (looped) playback of \code{.swm} files\\
2278
2279\end{tabular}
2280
2281\begin{tabular}{ll}
2282  \code{-movie <dirname>} & Save numbered images to named directory and
2283                                     quit\\
2284
2285  \code{-nosky} & Omit background sky\\
2286
2287
2288  \code{-scale <float>} & Vertical scale factor\\
2289  \code{-texture <file>} & Image to use for bedslope topography\\
2290  \code{-tps <rate>} & Timesteps per second\\
2291  \code{-version} & Revision number and creation (not compile)
2292                                     date\\
2293\end{tabular}
2294
2295\section{utilities/polygons}
2296
2297  \declaremodule{standard}{utilities.polygon}
2298  \refmodindex{utilities.polygon}
2299
2300  \begin{classdesc}{Polygon\_function}{regions, default = 0.0, geo_reference = None}
2301  Module: \code{utilities.polygon}
2302
2303  Creates a callable object that returns one of a specified list of values when
2304  evaluated at a point \code{x, y}, depending on which polygon, from a specified list of polygons, the
2305  point belongs to. The parameter \code{regions} is a list of pairs
2306  \code{(P, v)}, where each \code{P} is a polygon and each \code{v}
2307  is either a constant value or a function of coordinates \code{x}
2308  and \code{y}, specifying the return value for a point inside \code{P}. The
2309  optional parameter \code{default} may be used to specify a value
2310  for a point not lying inside any of the specified polygons. When a
2311  point lies in more than one polygon, the return value is taken to
2312  be the value for whichever of these polygon appears later in the
2313  list.
2314  %FIXME (Howard): CAN x, y BE VECTORS?
2315
2316  \end{classdesc}
2317
2318  \begin{funcdesc}{read\_polygon}{filename}
2319  Module: \code{utilities.polygon}
2320
2321  Reads the specified file and returns a polygon. Each
2322  line of the file must contain exactly two numbers, separated by a comma, which are interpreted
2323  as coordinates of one vertex of the polygon.
2324  \end{funcdesc}
2325
2326  \begin{funcdesc}{populate\_polygon}{polygon, number_of_points, seed = None, exclude = None}
2327  Module: \code{utilities.polygon}
2328
2329  Populates the interior of the specified polygon with the specified number of points,
2330  selected by means of a uniform distribution function.
2331  \end{funcdesc}
2332
2333  \begin{funcdesc}{point\_in\_polygon}{polygon, delta=1e-8}
2334  Module: \code{utilities.polygon}
2335
2336  Returns a point inside the specified polygon and close to the edge. The distance between
2337  the returned point and the nearest point of the polygon is less than $\sqrt{2}$ times the
2338  second argument \code{delta}, which is taken as $10^{-8}$ by default.
2339  \end{funcdesc}
2340
2341  \begin{funcdesc}{inside\_polygon}{points, polygon, closed = True, verbose = False}
2342  Module: \code{utilities.polygon}
2343
2344  Used to test whether the members of a list of points
2345  are inside the specified polygon. Returns a Numeric
2346  array comprising the indices of the points in the list that lie inside the polygon.
2347  (If none of the points are inside, returns \code{zeros((0,), 'l')}.)
2348  Points on the edges of the polygon are regarded as inside if
2349  \code{closed} is set to \code{True} or omitted; otherwise they are regarded as outside.
2350  \end{funcdesc}
2351
2352  \begin{funcdesc}{outside\_polygon}{points, polygon, closed = True, verbose = False}
2353  Module: \code{utilities.polygon}
2354
2355  Exactly like \code{inside\_polygon}, but with the words `inside' and `outside' interchanged.
2356  \end{funcdesc}
2357
2358  \begin{funcdesc}{is\_inside\_polygon}{point, polygon, closed=True, verbose=False}
2359  Module: \code{utilities.polygon}
2360
2361  Returns \code{True} if \code{point} is inside \code{polygon} or
2362  \code{False} otherwise. Points on the edges of the polygon are regarded as inside if
2363  \code{closed} is set to \code{True} or omitted; otherwise they are regarded as outside.
2364  \end{funcdesc}
2365
2366  \begin{funcdesc}{is\_outside\_polygon}{point, polygon, closed=True, verbose=False}
2367  Module: \code{utilities.polygon}
2368
2369  Exactly like \code{is\_outside\_polygon}, but with the words `inside' and `outside' interchanged.
2370  \end{funcdesc}
2371
2372  \begin{funcdesc}{point\_on\_line}{x, y, x0, y0, x1, y1}
2373  Module: \code{utilities.polygon}
2374
2375  Returns \code{True} or \code{False}, depending on whether the point with coordinates
2376  \code{x, y} is on the line passing through the points with coordinates \code{x0, y0}
2377  and \code{x1, y1} (extended if necessary at either end).
2378  \end{funcdesc}
2379
2380  \begin{funcdesc}{separate\_points\_by\_polygon}{points, polygon, closed = True, verbose = False}
2381    \indexedcode{separate\_points\_by\_polygon}
2382  Module: \code{utilities.polygon}
2383
2384  \end{funcdesc}
2385
2386  \begin{funcdesc}{polygon\_area}{polygon}
2387  Module: \code{utilities.polygon}
2388
2389  Returns area of arbitrary polygon (reference http://mathworld.wolfram.com/PolygonArea.html)
2390  \end{funcdesc}
2391
2392  \begin{funcdesc}{plot\_polygons}{polygons, figname, verbose = False}
2393  Module: \code{utilities.polygon}
2394
2395  Plots each polygon contained in input polygon list, e.g.
2396 \code{polygons = [poly1, poly2, poly3]} where \code{poly1 = [[x11,y11],[x12,y12],[x13,y13]]}
2397 etc.  Each polygon is closed for plotting purposes and subsequent plot saved to \code{figname}.
2398  Returns list containing the minimum and maximum of \code{x} and \code{y},
2399  i.e. \code{[x_{min}, x_{max}, y_{min}, y_{max}]}.
2400  \end{funcdesc}
2401
2402\section{coordinate\_transforms}
2403
2404\section{geospatial\_data}
2405
2406This describes a class that represents arbitrary point data in UTM
2407coordinates along with named attribute values.
2408
2409%FIXME (Ole): This gives a LaTeX error
2410%\declaremodule{standard}{geospatial_data}
2411%\refmodindex{geospatial_data}
2412
2413
2414
2415\begin{classdesc}{Geospatial\_data}
2416  {data_points = None,
2417    attributes = None,
2418    geo_reference = None,
2419    default_attribute_name = None,
2420    file_name = None}
2421Module: \code{geospatial\_data}
2422
2423This class is used to store a set of data points and associated
2424attributes, allowing these to be manipulated by methods defined for
2425the class.
2426
2427The data points are specified either by reading them from a NetCDF
2428or XYA file, identified through the parameter \code{file\_name}, or
2429by providing their \code{x}- and \code{y}-coordinates in metres,
2430either as a sequence of 2-tuples of floats or as an $M \times 2$
2431Numeric array of floats, where $M$ is the number of points.
2432Coordinates are interpreted relative to the origin specified by the
2433object \code{geo\_reference}, which contains data indicating the UTM
2434zone, easting and northing. If \code{geo\_reference} is not
2435specified, a default is used.
2436
2437Attributes are specified through the parameter \code{attributes},
2438set either to a list or array of length $M$ or to a dictionary whose
2439keys are the attribute names and whose values are lists or arrays of
2440length $M$. One of the attributes may be specified as the default
2441attribute, by assigning its name to \code{default\_attribute\_name}.
2442If no value is specified, the default attribute is taken to be the
2443first one.
2444\end{classdesc}
2445
2446
2447\begin{methoddesc}{import\_points\_file}{delimiter = None, verbose = False}
2448
2449\end{methoddesc}
2450
2451
2452\begin{methoddesc}{export\_points\_file}{ofile, absolute=False}
2453
2454\end{methoddesc}
2455
2456
2457\begin{methoddesc}{get\_data\_points}{absolute = True}
2458
2459\end{methoddesc}
2460
2461
2462\begin{methoddesc}{set\_attributes}{attributes}
2463
2464\end{methoddesc}
2465
2466
2467\begin{methoddesc}{get\_attributes}{attribute_name = None}
2468
2469\end{methoddesc}
2470
2471
2472\begin{methoddesc}{get\_all\_attributes}{}
2473
2474\end{methoddesc}
2475
2476
2477\begin{methoddesc}{set\_default\_attribute\_name}{default_attribute_name}
2478
2479\end{methoddesc}
2480
2481
2482\begin{methoddesc}{set\_geo\_reference}{geo_reference}
2483
2484\end{methoddesc}
2485
2486
2487\begin{methoddesc}{add}{}
2488
2489\end{methoddesc}
2490
2491
2492\section{pmesh GUI}
2493
2494\section{alpha\_shape}
2495\emph{Alpha shapes} are used to generate close-fitting boundaries
2496around sets of points. The alpha shape algorithm produces a shape
2497that approximates to the `shape formed by the points'---or the shape
2498that would be seen by viewing the points from a coarse enough
2499resolution. For the simplest types of point sets, the alpha shape
2500reduces to the more precise notion of the convex hull. However, for
2501many sets of points the convex hull does not provide a close fit and
2502the alpha shape usually fits more closely to the original point set,
2503offering a better approximation to the shape being sought.
2504
2505In \anuga, an alpha shape is used to generate a polygonal boundary
2506around a set of points before mesh generation. The algorithm uses a
2507parameter $\alpha$ that can be adjusted to make the resultant shape
2508resemble the shape suggested by intuition more closely. An alpha
2509shape can serve as an initial boundary approximation that the user
2510can adjust as needed.
2511
2512The following paragraphs describe the class used to model an alpha
2513shape and some of the important methods and attributes associated
2514with instances of this class.
2515
2516\begin{classdesc}{Alpha\_Shape}{points, alpha = None}
2517Module: \code{alpha\_shape}
2518
2519To instantiate this class the user supplies the points from which
2520the alpha shape is to be created (in the form of a list of 2-tuples
2521\code{[[x1, y1],[x2, y2]}\ldots\code{]}, assigned to the parameter
2522\code{points}) and, optionally, a value for the parameter
2523\code{alpha}. The alpha shape is then computed and the user can then
2524retrieve details of the boundary through the attributes defined for
2525the class.
2526\end{classdesc}
2527
2528
2529\begin{funcdesc}{alpha\_shape\_via\_files}{point_file, boundary_file, alpha= None}
2530Module: \code{alpha\_shape}
2531
2532This function reads points from the specified point file
2533\code{point\_file}, computes the associated alpha shape (either
2534using the specified value for \code{alpha} or, if no value is
2535specified, automatically setting it to an optimal value) and outputs
2536the boundary to a file named \code{boundary\_file}. This output file
2537lists the coordinates \code{x, y} of each point in the boundary,
2538using one line per point.
2539\end{funcdesc}
2540
2541
2542\begin{methoddesc}{set\_boundary\_type}{self,raw_boundary=True,
2543                          remove_holes=False,
2544                          smooth_indents=False,
2545                          expand_pinch=False,
2546                          boundary_points_fraction=0.2}
2547Module: \code{alpha\_shape},  Class: \class{Alpha\_Shape}
2548
2549This function sets flags that govern the operation of the algorithm
2550that computes the boundary, as follows:
2551
2552\code{raw\_boundary = True} returns raw boundary, i.e. the regular edges of the
2553                alpha shape.\\
2554\code{remove\_holes = True} removes small holes (`small' is defined by
2555\code{boundary\_points\_fraction})\\
2556\code{smooth\_indents = True} removes sharp triangular indents in
2557boundary\\
2558\code{expand\_pinch = True} tests for pinch-off and
2559corrects---preventing a boundary vertex from having more than two edges.
2560\end{methoddesc}
2561
2562
2563\begin{methoddesc}{get\_boundary}{}
2564Module: \code{alpha\_shape},  Class: \class{Alpha\_Shape}
2565
2566Returns a list of tuples representing the boundary of the alpha
2567shape. Each tuple represents a segment in the boundary by providing
2568the indices of its two endpoints.
2569\end{methoddesc}
2570
2571
2572\begin{methoddesc}{write\_boundary}{file_name}
2573Module: \code{alpha\_shape},  Class: \class{Alpha\_Shape}
2574
2575Writes the list of 2-tuples returned by \code{get\_boundary} to the
2576file \code{file\_name}, using one line per tuple.
2577\end{methoddesc}
2578
2579\section{Numerical Tools}
2580
2581The following table describes some useful numerical functions that
2582may be found in the module \module{utilities.numerical\_tools}:
2583
2584\begin{tabular}{|p{8cm} p{8cm}|}  \hline
2585\code{angle(v1, v2=None)} & Angle between two-dimensional vectors
2586\code{v1} and \code{v2}, or between \code{v1} and the $x$-axis if
2587\code{v2} is \code{None}. Value is in range $0$ to $2\pi$. \\
2588
2589\code{normal\_vector(v)} & Normal vector to \code{v}.\\
2590
2591\code{mean(x)} & Mean value of a vector \code{x}.\\
2592
2593\code{cov(x, y=None)} & Covariance of vectors \code{x} and \code{y}.
2594If \code{y} is \code{None}, returns \code{cov(x, x)}.\\
2595
2596\code{err(x, y=0, n=2, relative=True)} & Relative error of
2597$\parallel$\code{x}$-$\code{y}$\parallel$ to
2598$\parallel$\code{y}$\parallel$ (2-norm if \code{n} = 2 or Max norm
2599if \code{n} = \code{None}). If denominator evaluates to zero or if
2600\code{y}
2601is omitted or if \code{relative = False}, absolute error is returned.\\
2602
2603\code{norm(x)} & 2-norm of \code{x}.\\
2604
2605\code{corr(x, y=None)} & Correlation of \code{x} and \code{y}. If
2606\code{y} is \code{None} returns autocorrelation of \code{x}.\\
2607
2608\code{ensure\_numeric(A, typecode = None)} & Returns a Numeric array
2609for any sequence \code{A}. If \code{A} is already a Numeric array it
2610will be returned unaltered. Otherwise, an attempt is made to convert
2611it to a Numeric array. (Needed because \code{array(A)} can
2612cause memory overflow.)\\
2613
2614\code{histogram(a, bins, relative=False)} & Standard histogram. If
2615\code{relative} is \code{True}, values will be normalised against
2616the total and thus represent frequencies rather than counts.\\
2617
2618\code{create\_bins(data, number\_of\_bins = None)} & Safely create
2619bins for use with histogram. If \code{data} contains only one point
2620or is constant, one bin will be created. If \code{number\_of\_bins}
2621is omitted, 10 bins will be created.\\  \hline
2622
2623\end{tabular}
2624
2625
2626\chapter{Modules available in \anuga}
2627
2628
2629\section{\module{pyvolution.general\_mesh} }
2630\declaremodule[pyvolution.generalmesh]{}{pyvolution.general\_mesh}
2631\label{mod:pyvolution.generalmesh}
2632
2633\section{\module{pyvolution.neighbour\_mesh} }
2634\declaremodule[pyvolution.neighbourmesh]{}{pyvolution.neighbour\_mesh}
2635\label{mod:pyvolution.neighbourmesh}
2636
2637\section{\module{pyvolution.domain} --- Generic module for 2D triangular domains for finite-volume computations of conservation laws}
2638\declaremodule{}{pyvolution.domain}
2639\label{mod:pyvolution.domain}
2640
2641
2642\section{\module{pyvolution.quantity}}
2643\declaremodule{}{pyvolution.quantity}
2644\label{mod:pyvolution.quantity}
2645
2646
2647Class Quantity - Implements values at each triangular element
2648
2649To create:
2650
2651   Quantity(domain, vertex_values)
2652
2653   domain: Associated domain structure. Required.
2654
2655   vertex_values: N x 3 array of values at each vertex for each element.
2656                  Default None
2657
2658   If vertex_values are None Create array of zeros compatible with domain.
2659   Otherwise check that it is compatible with dimenions of domain.
2660   Otherwise raise an exception
2661
2662
2663
2664\section{\module{pyvolution.shallow\_water} --- 2D triangular domains for finite-volume
2665computations of the shallow water wave equation. This module contains a specialisation
2666of class Domain from module domain.py consisting of methods specific to the Shallow Water
2667Wave Equation
2668}
2669\declaremodule[pyvolution.shallowwater]{}{pyvolution.shallow\_water}
2670\label{mod:pyvolution.shallowwater}
2671
2672
2673
2674
2675%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2676
2677\chapter{Frequently Asked Questions}
2678
2679
2680\section{General Questions}
2681
2682\subsubsection{What is \anuga?}
2683
2684\subsubsection{Why is it called \anuga?}
2685
2686\subsubsection{How do I obtain a copy of \anuga?}
2687
2688\subsubsection{What developments are expected for \anuga in the future?}
2689
2690\subsubsection{Are there any published articles about \anuga that I can reference?}
2691
2692\section{Modelling Questions}
2693
2694\subsubsection{Which type of problems are \anuga good for?}
2695
2696\subsubsection{Which type of problems are beyond the scope of \anuga?}
2697
2698\subsubsection{Can I start the simulation at an arbitrary time?}
2699Yes, using \code{domain.set\_time()} you can specify an arbitrary
2700starting time. This is for example useful in conjunction with a
2701file\_boundary, which may start hours before anything hits the model
2702boundary. By assigning a later time for the model to start,
2703computational resources aren't wasted.
2704
2705\subsubsection{Can I change values for any quantity during the simulation?}
2706Yes, using \code{domain.set\_quantity()} inside the domain.evolve
2707loop you can change values of any quantity. This is for example
2708useful if you wish to let the system settle for a while before
2709assigning an initial condition. Another example would be changing
2710the values for elevation to model e.g. erosion.
2711
2712\subsubsection{Can I change boundary conditions during the simulation?}
2713Not sure, but it would be nice :-)
2714
2715\subsubsection{Why does a file\_function return a list of numbers when evaluated?}
2716Currently, file\_function works by returning values for the conserved
2717quantities \code{stage}, \code{xmomentum} and \code{ymomentum} at a given point in time
2718and space as a triplet. To access e.g.\ \code{stage} one must specify element 0 of the
2719triplet returned by file\_function.
2720
2721\subsubsection{Which diagnostics are available to troubleshoot a simulation?}
2722
2723\subsubsection{How do I use a DEM in my simulation?}
2724You use \code{dem2pts} to convert your DEM to the required .pts format. This .pts file is then called
2725when setting the elevation data to the mesh in \code{domain.set_quantity}
2726
2727\subsubsection{What sort of DEM resolution should I use?}
2728Try and work with the \emph{best} you have available. Onshore DEMs
2729are typically available in 25m, 100m and 250m grids. Note, offshore
2730data is often sparse, or non-existent.
2731
2732\subsubsection{What sort of mesh resolution should I use?}
2733The mesh resolution should be commensurate with your DEM - it does not make sense to put in place a mesh which is finer than your DEM. As an example,
2734if your DEM is on a 25m grid, then the cell resolution should be of the order of 315$m^2$ (this represents half the area of the square grid). Ideally,
2735you need a fine mesh over regions where the DEM changes rapidly, and other areas of significant interest, such as the coast.
2736
2737
2738\subsubsection{How do I tag interior polygons?}
2739At the moment create_mesh_from_regions does not allow interior
2740polygons with symbolic tags. If tags are needed, the interior
2741polygons must be created subsequently. For example, given a filename
2742of polygons representing solid walls (in Arc Ungenerate format) can
2743be tagged as such using the code snippet:
2744\begin{verbatim}
2745  # Create mesh outline with tags
2746  mesh = create_mesh_from_regions(bounding_polygon,
2747                                  boundary_tags=boundary_tags)
2748  # Add buildings outlines with tags set to 'wall'. This would typically
2749  # bind to a Reflective boundary
2750  mesh.import_ungenerate_file(buildings_filename, tag='wall')
2751
2752  # Generate and write mesh to file
2753  mesh.generate_mesh(maximum_triangle_area=max_area)
2754  mesh.export_mesh_file(mesh_filename)
2755\end{verbatim}
2756
2757Note that a mesh object is returned from \code{create_mesh_from_regions}
2758when file name is omitted.
2759
2760\subsubsection{How often should I store the output?}
2761This will depend on what you are trying to answer with your model and how much memory you have available on your machine. If you need
2762to look in detail at the evolution, then you will need to balance your storage requirements and the duration of the simulation.
2763If the SWW file exceeds 1Gb, another SWW file will be created until the end of the simulation. As an example, to store all the conserved
2764quantities on a mesh with approximately 300000 triangles on a 2 min interval for 5 hours will result in approximately 350Mb SWW file
2765(as for the \file{run\_sydney\_smf.py} example).
2766
2767\subsection{Boundary Conditions}
2768
2769\subsubsection{How do I create a Dirichlet boundary condition?}
2770
2771A Dirichlet boundary condition sets a constant value for the
2772conserved quantities at the boundaries. A list containing
2773the constant values for stage, xmomentum and ymomentum is constructed
2774and used in the function call, e.g. \code{Dirichlet_boundary([0.2,0.,0.])}
2775
2776\subsubsection{How do I know which boundary tags are available?}
2777The method \code{domain.get\_boundary\_tags()} will return a list of
2778available tags for use with
2779\code{domain.set\_boundary\_condition()}.
2780
2781
2782
2783
2784
2785\chapter{Glossary}
2786
2787\begin{tabular}{|lp{10cm}|c|}  \hline
2788%\begin{tabular}{|llll|}  \hline
2789    \emph{Term} & \emph{Definition} & \emph{Page}\\  \hline
2790
2791    \indexedbold{\anuga} & Name of software (joint development between ANU and
2792    GA) & \pageref{def:anuga}\\
2793
2794    \indexedbold{bathymetry} & offshore elevation &\\
2795
2796    \indexedbold{conserved quantity} & conserved (stage, x and y
2797    momentum) & \\
2798
2799%    \indexedbold{domain} & The domain of a function is the set of all input values to the
2800%    function.&\\
2801
2802    \indexedbold{Digital Elevation Model (DEM)} & DEMs are digital files consisting of points of elevations,
2803sampled systematically at equally spaced intervals.& \\
2804
2805    \indexedbold{Dirichlet boundary} & A boundary condition imposed on a differential equation
2806 that specifies the values the solution is to take on the boundary of the
2807 domain. & \pageref{def:dirichlet boundary}\\
2808
2809    \indexedbold{edge} & A triangular cell within the computational mesh can be depicted
2810    as a set of vertices joined by lines (the edges). & \\
2811
2812    \indexedbold{elevation} & refers to bathymetry and topography &\\
2813
2814    \indexedbold{evolution} & integration of the shallow water wave equations
2815    over time &\\
2816
2817    \indexedbold{finite volume method} & The method evaluates the terms in the shallow water
2818    wave equation as fluxes at the surfaces of each finite volume. Because the
2819    flux entering a given volume is identical to that leaving the adjacent volume,
2820    these methods are conservative. Another advantage of the finite volume method is
2821    that it is easily formulated to allow for unstructured meshes. The method is used
2822    in many computational fluid dynamics packages. & \\
2823
2824    \indexedbold{forcing term} & &\\
2825
2826    \indexedbold{flux} & the amount of flow through the volume per unit
2827    time & \\
2828
2829    \indexedbold{grid} & Evenly spaced mesh & \\
2830
2831    \indexedbold{latitude} & The angular distance on a mericlear north and south of the
2832    equator, expressed in degrees and minutes. & \\
2833
2834    \indexedbold{longitude} & The angular distance east or west, between the meridian
2835    of a particular place on Earth and that of the Prime Meridian (located in Greenwich,
2836    England) expressed in degrees or time.& \\
2837
2838    \indexedbold{Manning friction coefficient} & &\\
2839
2840    \indexedbold{mesh} & Triangulation of domain &\\
2841
2842    \indexedbold{mesh file} & A TSH or MSH file & \pageref{def:mesh file}\\
2843
2844    \indexedbold{NetCDF} & &\\
2845
2846    \indexedbold{northing} & A rectangular (x,y) coordinate measurement of distance
2847    north from a north-south reference line, usually a meridian used as the axis of
2848    origin within a map zone or projection. Northing is a UTM (Universal Transverse
2849    Mercator) coordinate. & \\
2850
2851    \indexedbold{points file} & A PTS or XYA file & \\  \hline
2852
2853    \end{tabular}
2854
2855    \begin{tabular}{|lp{10cm}|c|}  \hline
2856
2857    \indexedbold{polygon} & A sequence of points in the plane. \anuga represents a polygon
2858    either as a list consisting of Python tuples or lists of length 2 or as an $N \times 2$
2859    Numeric array, where $N$ is the number of points.
2860
2861    The unit square, for example, would be represented either as
2862    \code{[ [0,0], [1,0], [1,1], [0,1] ]} or as \code{array( [0,0], [1,0], [1,1],
2863    [0,1] )}.
2864
2865    NOTE: For details refer to the module \module{utilities/polygon.py}. &
2866    \\     \indexedbold{resolution} &  The maximal area of a triangular cell in a
2867    mesh & \\
2868
2869
2870    \indexedbold{reflective boundary} & Models a solid wall. Returns same conserved
2871    quantities as those present in the neighbouring volume but reflected. Specific to the
2872    shallow water equation as it works with the momentum quantities assumed to be the
2873    second and third conserved quantities. & \pageref{def:reflective boundary}\\
2874
2875    \indexedbold{stage} & &\\
2876
2877%    \indexedbold{try this}
2878
2879    \indexedbold{swollen} & visualisation tool used with \anuga & \pageref{sec:swollen}\\
2880
2881    \indexedbold{time boundary} & Returns values for the conserved
2882quantities as a function of time. The user must specify
2883the domain to get access to the model time. & \pageref{def:time boundary}\\
2884
2885    \indexedbold{topography} & onshore elevation &\\
2886
2887    \indexedbold{transmissive boundary} & & \pageref{def:transmissive boundary}\\
2888
2889    \indexedbold{vertex} & A point at which edges meet. & \\
2890
2891    \indexedbold{xmomentum} & conserved quantity (note, two-dimensional SWW equations say
2892    only \code{x} and \code{y} and NOT \code{z}) &\\
2893
2894    \indexedbold{ymomentum}  & conserved quantity & \\  \hline
2895
2896    \end{tabular}
2897
2898
2899%The \code{\e appendix} markup need not be repeated for additional
2900%appendices.
2901
2902
2903%
2904%  The ugly "%begin{latexonly}" pseudo-environments are really just to
2905%  keep LaTeX2HTML quiet during the \renewcommand{} macros; they're
2906%  not really valuable.
2907%
2908%  If you don't want the Module Index, you can remove all of this up
2909%  until the second \input line.
2910%
2911
2912%begin{latexonly}
2913%\renewcommand{\indexname}{Module Index}
2914%end{latexonly}
2915\input{mod\jobname.ind}        % Module Index
2916%
2917%begin{latexonly}
2918%\renewcommand{\indexname}{Index}
2919%end{latexonly}
2920\input{\jobname.ind}            % Index
2921
2922
2923
2924\begin{thebibliography}{99}
2925\bibitem[nielsen2005]{nielsen2005} 
2926{\it Hydrodynamic modelling of coastal inundation}.
2927Nielsen, O., S. Roberts, D. Gray, A. McPherson and A. Hitchman.
2928In Zerger, A. and Argent, R.M. (eds) MODSIM 2005 International Congress on
2929Modelling and Simulation. Modelling and Simulation Society of Australia and
2930New Zealand, December 2005, pp. 518-523. ISBN: 0-9758400-2-9.\\
2931http://www.mssanz.org.au/modsim05/papers/nielsen.pdf
2932
2933
2934
2935
2936\end{thebibliography}{99}
2937
2938\end{document}
Note: See TracBrowser for help on using the repository browser.