1 | % Use the standard \LaTeXe\ article style in 12pt Computer Modern font |
---|
2 | % on A4 paper by |
---|
3 | \documentclass[12pt,a4paper]{article} |
---|
4 | |
---|
5 | % Do \emph{not} change the width nor the height of the text from the |
---|
6 | % defaults set by this document class. |
---|
7 | % |
---|
8 | % The alternative which is closer to what we actually use is |
---|
9 | %\documentclass[11pt,a5paper]{article} |
---|
10 | %\usepackage[a5paper]{geometry} |
---|
11 | % Because it is a great size for on screen reading |
---|
12 | % and prints nicely on a4paper either 2up or booklet. |
---|
13 | |
---|
14 | % The preamble is to contain your own \LaTeX\ commands and to say |
---|
15 | % what packages to use. Three useful packages are the following: |
---|
16 | \usepackage{graphicx} % avoid epsfig or earlier such packages |
---|
17 | |
---|
18 | |
---|
19 | \usepackage{url} % for URLs and DOIs |
---|
20 | \usepackage{amsmath} % many want amsmath extensions |
---|
21 | \usepackage{amsfonts} |
---|
22 | \usepackage{underscore} |
---|
23 | % Avoid loading unused packages (as done by some \LaTeX\ editors). |
---|
24 | |
---|
25 | % Create title and authors using \verb|\maketitle|. Separate authors by |
---|
26 | % \verb|\and| and put addresses in \verb|\thanks| command with |
---|
27 | % \verb|\url| command \verb|\protect|ed. |
---|
28 | \title{On The Validation of A Hydrodynamic Model} |
---|
29 | |
---|
30 | \author{ |
---|
31 | D.~S.~Gray\thanks{Natural Hazard Impacts Project, Geospatial and |
---|
32 | Earth Monitoring Division, Geoscience Australia, Symonston, |
---|
33 | \textsc{Australia}. \protect\url{mailto:Duncan.Gray@ga.gov.au}}\footnotemark[1] |
---|
34 | \and |
---|
35 | T.~Baldock\thanks{University of Queensland, Brisbane, \textsc{Australia}. |
---|
36 | \protect\url{mailto:tom.baldock@uq.edu.au}}\footnotemark[2] |
---|
37 | \and |
---|
38 | O.~M.~Nielsen\footnotemark[1] |
---|
39 | \and |
---|
40 | M.~J.~Sexton\footnotemark[1] |
---|
41 | \and |
---|
42 | N.~Bartzis\footnotemark[1] |
---|
43 | \and |
---|
44 | S.~G.~Roberts\thanks{Department of Mathematics, Australian National University, |
---|
45 | Canberra, \textsc{Australia}. \protect\url{mailto:stephen.roberts@anu.edu.au}}} |
---|
46 | |
---|
47 | \date{30 May 2008} |
---|
48 | |
---|
49 | \newcommand{\ANUGA}{\textsc{ANUGA}} |
---|
50 | \newcommand{\Python}{\textsc{Python}} |
---|
51 | \newcommand{\VPython}{\textsc{VPython}} |
---|
52 | \newcommand{\pypar}{\textsc{mpi}} |
---|
53 | \newcommand{\Metis}{\textsc{Metis}} |
---|
54 | \newcommand{\mpi}{\textsc{mpi}} |
---|
55 | |
---|
56 | \newcommand{\UU}{\mathbf{U}} |
---|
57 | \newcommand{\VV}{\mathbf{V}} |
---|
58 | \newcommand{\EE}{\mathbf{E}} |
---|
59 | \newcommand{\GG}{\mathbf{G}} |
---|
60 | \newcommand{\FF}{\mathbf{F}} |
---|
61 | \newcommand{\HH}{\mathbf{H}} |
---|
62 | \newcommand{\SSS}{\mathbf{S}} |
---|
63 | \newcommand{\nn}{\mathbf{n}} |
---|
64 | |
---|
65 | \newcommand{\code}[1]{\texttt{#1}} |
---|
66 | |
---|
67 | |
---|
68 | |
---|
69 | |
---|
70 | \begin{document} |
---|
71 | |
---|
72 | % Use default \verb|\maketitle|. |
---|
73 | \maketitle |
---|
74 | |
---|
75 | % Use the \verb|abstract| environment. |
---|
76 | \begin{abstract} |
---|
77 | Modelling the effects on the built environment of natural hazards such |
---|
78 | as riverine flooding, storm surges and tsunami is critical for |
---|
79 | understanding their economic and social impact on our urban |
---|
80 | communities. Geoscience Australia and the Australian National |
---|
81 | University have developed a hydrodynamic inundation modelling tool |
---|
82 | called \ANUGA{} to help simulate the impact of these hazards. |
---|
83 | The core of \ANUGA{} is a \Python{} implementation of a finite-volume method |
---|
84 | for solving the conservative form of the Shallow Water Wave equation. |
---|
85 | |
---|
86 | In this paper, a number of tests are performed to validate \ANUGA{}. These tests |
---|
87 | range from benchmark problems to wave and flume tank examples. |
---|
88 | \ANUGA{} is available as Open Source to enable |
---|
89 | free access to the software and allow the scientific community to |
---|
90 | use, validate and contribute to the software in the future. |
---|
91 | |
---|
92 | %This method allows the study area to be represented by an unstructured |
---|
93 | %mesh with variable resolution to suit the particular problem. The |
---|
94 | %conserved quantities are water level (stage) and horizontal momentum. |
---|
95 | %An important capability of ANUGA is that it can robustly model the |
---|
96 | %process of wetting and drying as water enters and leaves an area. This |
---|
97 | %means that it is suitable for simulating water flow onto a beach or |
---|
98 | %dry land and around structures such as buildings. |
---|
99 | |
---|
100 | |
---|
101 | \end{abstract} |
---|
102 | |
---|
103 | % By default we include a table of contents in each paper. |
---|
104 | %\tableofcontents |
---|
105 | |
---|
106 | % Use \verb|\section|, \verb|\subsection|, \verb|\subsubsection| and |
---|
107 | % possibly \verb|\paragraph| to structure your document. Make sure |
---|
108 | % you \verb|\label| them for cross-referencing with \verb|\ref|\,. |
---|
109 | |
---|
110 | %\clearpage |
---|
111 | \section{Introduction} |
---|
112 | \label{sec:intro} |
---|
113 | |
---|
114 | Hydrodynamic modelling allows impacts from flooding, storm-surge and |
---|
115 | tsunami to be better understood, their impacts to be anticipated and, |
---|
116 | with appropriate planning, their effects to be mitigated. A significant |
---|
117 | proportion of the Australian population reside in the coastal |
---|
118 | corridors, thus the potential of significant disruption and loss |
---|
119 | is real. The extent of |
---|
120 | inundation is critically linked to the event, tidal conditions, |
---|
121 | bathymetry and topography and it not feasible to make impact |
---|
122 | predictions using heuristics alone. |
---|
123 | Geoscience |
---|
124 | Australia in collaboration with the Mathematical Sciences Institute, |
---|
125 | Australian National University, is developing a software application |
---|
126 | called \ANUGA{} to model the hydrodynamics of floods, storm surges and |
---|
127 | tsunami. These hazards are modelled using the conservative shallow |
---|
128 | water equations which are described in section~\ref{sec:model}. In |
---|
129 | \ANUGA{} these equations are solved using a finite volume method as |
---|
130 | described in section~\ref{sec:model}. A more complete discussion of the |
---|
131 | method can be found in \cite{modsim2005} where the model and solution |
---|
132 | technique is validated on a standard tsunami benchmark data set |
---|
133 | or in \cite{Roberts2007} where parallelisation of ANUGA is discussed. |
---|
134 | This modelling capability is part of |
---|
135 | Geoscience Australia's ongoing research effort to model and |
---|
136 | understand the potential impact from natural hazards in order to |
---|
137 | reduce their impact on Australian communities (see \cite{Nielsen2006}). |
---|
138 | \ANUGA{} is currently being trialled for flood |
---|
139 | modelling (see \cite{Rigby2008}). |
---|
140 | |
---|
141 | The validity of other hydrodynamic models have been reported elsewhere, |
---|
142 | with Hubbard and Dodd (2002) \cite{Hubbard02} providing |
---|
143 | an excellent review of 1D and 2D models and associated validation tests. They |
---|
144 | described the evolution of these models from fixed, nested to adaptive grids |
---|
145 | and the ability of the solvers to cope with the moving shoreline. They highlighted the |
---|
146 | difficulty in verify the nonlinear shallow water equations themselves as the only |
---|
147 | standard analytical solution is that of Carrier and Greenspan (1958) |
---|
148 | \cite{Carrier58} that is strictly |
---|
149 | for non-breaking waves. Further, whilst there is a 2D analytic solution from Thacker (1981), it |
---|
150 | appears that the circular island wave tank example of Briggs et al will become |
---|
151 | the standard data set to verify the equations. |
---|
152 | |
---|
153 | This paper will describe the validation outputs in a similar way to Hubbard and Dodd |
---|
154 | \cite{Hubbard02} to |
---|
155 | present an exhaustive validation of the numerical model. Further to these tests, we will |
---|
156 | incorporate a test to verify friction values. The six tests are: |
---|
157 | (1) verification against the 1D analytical solution of Carrier and Greenspan; |
---|
158 | (2) testing against 1D (flume) data sets to verify wave height and velocity |
---|
159 | (3) determining friction values from 1D flume data sets; |
---|
160 | (4) validation against a genuinely 2D analytical solution of the model equations; |
---|
161 | (5) testing against the 2D Okushiri benchmark problem; and |
---|
162 | (6) testing against the 2D data sets modelling wave run-up around a circular island by Briggs et al. |
---|
163 | Throughout the paper, qualitative comparisons will be drawn against other models. |
---|
164 | |
---|
165 | %Hubbard and Dodd's model, OTT-2D, has some similarities to \ANUGA{}, and |
---|
166 | %whilst the mesh can be refined, it is based on rectangular mesh. |
---|
167 | |
---|
168 | The \ANUGA{} model and numerical scheme is briefly described in section~\ref{sec:model}. |
---|
169 | A detailed description of the numerical scheme and software implementation can be found in |
---|
170 | the MODSIM, CTAC etc papers. The six case studies to validation and verify \ANUGA{} will be |
---|
171 | presented in section~\ref{sec:validation}, with the conclusions |
---|
172 | outlined in section~\ref{sec:conclusions}. |
---|
173 | |
---|
174 | {\bf question - if the Okushiri result has already been presented in the |
---|
175 | MODSIM paper, how should it be presented in this paper? - simply refer to it I think} |
---|
176 | |
---|
177 | \section{Model} |
---|
178 | \label{sec:model} |
---|
179 | |
---|
180 | The shallow water wave equations are a system of differential |
---|
181 | conservation equations which describe the flow of a thin layer of |
---|
182 | fluid over terrain. The form of the equations are: |
---|
183 | \[ |
---|
184 | \frac{\partial \UU}{\partial t}+\frac{\partial \EE}{\partial |
---|
185 | x}+\frac{\partial \GG}{\partial y}=\SSS |
---|
186 | \] |
---|
187 | where $\UU=\left[ {{\begin{array}{*{20}c} |
---|
188 | h & {uh} & {vh} \\ |
---|
189 | \end{array} }} \right]^T$ is the vector of conserved quantities; water depth |
---|
190 | $h$, $x$-momentum $uh$ and $y$-momentum $vh$. Other quantities |
---|
191 | entering the system are bed elevation $z$ and stage (absolute water |
---|
192 | level) $w$, where the relation $w = z + h$ holds true at all times. |
---|
193 | The fluxes in the $x$ and $y$ directions, $\EE$ and $\GG$ are given |
---|
194 | by |
---|
195 | \[ |
---|
196 | \EE=\left[ {{\begin{array}{*{20}c} |
---|
197 | {uh} \hfill \\ |
---|
198 | {u^2h+gh^2/2} \hfill \\ |
---|
199 | {uvh} \hfill \\ |
---|
200 | \end{array} }} \right]\mbox{ and }\GG=\left[ {{\begin{array}{*{20}c} |
---|
201 | {vh} \hfill \\ |
---|
202 | {vuh} \hfill \\ |
---|
203 | {v^2h+gh^2/2} \hfill \\ |
---|
204 | \end{array} }} \right] |
---|
205 | \] |
---|
206 | and the source term (which includes gravity and friction) is given |
---|
207 | by |
---|
208 | \[ |
---|
209 | \SSS=\left[ {{\begin{array}{*{20}c} |
---|
210 | 0 \hfill \\ |
---|
211 | -{gh(z_{x} + S_{fx} )} \hfill \\ |
---|
212 | -{gh(z_{y} + S_{fy} )} \hfill \\ |
---|
213 | \end{array} }} \right] |
---|
214 | \] |
---|
215 | where $S_f$ is the bed friction. The friction term is modelled using |
---|
216 | Manning's resistance law |
---|
217 | \[ |
---|
218 | S_{fx} =\frac{u\eta ^2\sqrt {u^2+v^2} }{h^{4/3}}\mbox{ and }S_{fy} |
---|
219 | =\frac{v\eta ^2\sqrt {u^2+v^2} }{h^{4/3}} |
---|
220 | \] |
---|
221 | in which $\eta$ is the Manning resistance coefficient. |
---|
222 | |
---|
223 | As demonstrated in our papers, \cite{modsim2005,Rob99l} these |
---|
224 | equations provide an excellent model of flows associated with |
---|
225 | inundation such as dam breaks and tsunamis. Question - how do we |
---|
226 | know it is excellent? |
---|
227 | |
---|
228 | \ANUGA{} uses a finite-volume method as |
---|
229 | described in \cite{Rob991} where the study area is represented by an |
---|
230 | unstructured triangular mesh. The flexibility afforded by this approach |
---|
231 | is the ability for the user to refine the mesh in areas of interest. |
---|
232 | \ANUGA{} is mostly written in the object-oriented programming |
---|
233 | language \Python{} with computationally intensive parts implemented |
---|
234 | as highly optimised shared objects written in C. The API is a |
---|
235 | \Python{} script where the user sets up the scenario. This script |
---|
236 | defines the study area, mesh refinement as well as initial and boundary conditions. |
---|
237 | The user is free to update quantity values or boundary conditions through |
---|
238 | the simulation. Reference to user manual |
---|
239 | |
---|
240 | Could include here a brief overview of the numerical |
---|
241 | technique and reference the CTAC 2006 paper and has Steve written something |
---|
242 | that could be included? |
---|
243 | |
---|
244 | \section{Finite Volume Method} |
---|
245 | \label{sec:fvm} |
---|
246 | |
---|
247 | {\bf Jane: I don't think this section is needed here, but the |
---|
248 | content is referred to at the end of section 1} |
---|
249 | |
---|
250 | We use a finite-volume method for solving the shallow water wave |
---|
251 | equations \cite{Rob99l}. The study area is represented by a mesh of |
---|
252 | triangular cells as in Figure~\ref{fig:mesh} in which the conserved |
---|
253 | quantities of water depth $h$, and horizontal momentum $(uh, vh)$, |
---|
254 | in each volume are to be determined. The size of the triangles may |
---|
255 | be varied within the mesh to allow greater resolution in regions of |
---|
256 | particular interest. |
---|
257 | |
---|
258 | \begin{figure} |
---|
259 | \begin{center} |
---|
260 | %\includegraphics[width=5.0cm,keepaspectratio=true]{step-five} |
---|
261 | \caption{Triangular mesh used in our finite volume method. Conserved |
---|
262 | quantities $h$, $uh$ and $vh$ are associated with the centroid of |
---|
263 | each triangular cell.} \label{fig:mesh} |
---|
264 | \end{center} |
---|
265 | \end{figure} |
---|
266 | |
---|
267 | The equations constituting the finite-volume method are obtained by |
---|
268 | integrating the differential conservation equations over each |
---|
269 | triangular cell of the mesh. Introducing some notation we use $i$ to |
---|
270 | refer to the $i$th triangular cell $T_i$, and ${\cal N}(i)$ to the |
---|
271 | set of indices referring to the cells neighbouring the $i$th cell. |
---|
272 | Then $A_i$ is the area of the $i$th triangular cell and $l_{ij}$ is |
---|
273 | the length of the edge between the $i$th and $j$th cells. |
---|
274 | |
---|
275 | By applying the divergence theorem we obtain for each volume an |
---|
276 | equation which describes the rate of change of the average of the |
---|
277 | conserved quantities within each cell, in terms of the fluxes across |
---|
278 | the edges of the cells and the effect of the source terms. In |
---|
279 | particular, rate equations associated with each cell have the form |
---|
280 | $$ |
---|
281 | \frac{d\UU_i }{dt}+ \frac1{A_i}\sum\limits_{j\in{\cal N}(i)} \HH_{ij} l_{ij} = \SSS_i |
---|
282 | $$ |
---|
283 | where |
---|
284 | \begin{itemize} |
---|
285 | \item $\UU_i$ the vector of conserved quantities averaged over the $i$th cell, |
---|
286 | \item $\SSS_i$ is the source term associated with the $i$th cell, |
---|
287 | and |
---|
288 | \item $\HH_{ij}$ is the outward normal flux of |
---|
289 | material across the \textit{ij}th edge. |
---|
290 | \end{itemize} |
---|
291 | |
---|
292 | |
---|
293 | %\item $l_{ij}$ is the length of the edge between the $i$th and $j$th |
---|
294 | %cells |
---|
295 | %\item $m_{ij}$ is the midpoint of |
---|
296 | %the \textit{ij}th edge, |
---|
297 | %\item |
---|
298 | %$\mathbf{n}_{ij} = (n_{ij,1} , n_{ij,2})$is the outward pointing |
---|
299 | %normal along the \textit{ij}th edge, and The |
---|
300 | |
---|
301 | The flux $\HH_{ij}$ is evaluated using a numerical flux function |
---|
302 | $\HH(\cdot, \cdot ; \ \cdot)$ which is consistent with the shallow |
---|
303 | water flux in the sense that for all conservation vectors $\UU$ and normal vectors $\nn$ |
---|
304 | $$ |
---|
305 | H(\UU,\UU;\ \nn) = \EE(\UU) n_1 + \GG(\UU) n_2 . |
---|
306 | $$ |
---|
307 | |
---|
308 | Then |
---|
309 | $$ |
---|
310 | \HH_{ij} = \HH(\UU_i(m_{ij}), |
---|
311 | \UU_j(m_{ij}); \mathbf{n}_{ij}) |
---|
312 | $$ |
---|
313 | where $m_{ij}$ is the midpoint of the \textit{ij}th edge and |
---|
314 | $\mathbf{n}_{ij}$ is the outward pointing normal, with respect to the $i$th cell, on the |
---|
315 | \textit{ij}th edge. The function $\UU_i(x)$ for $x \in |
---|
316 | T_i$ is obtained from the vector $\UU_k$ of conserved average values for the $i$th and |
---|
317 | neighbouring cells. |
---|
318 | |
---|
319 | We use a second order reconstruction to produce a piece-wise linear |
---|
320 | function construction of the conserved quantities for all $x \in |
---|
321 | T_i$ for each cell (see Figure~\ref{fig:mesh:reconstruct}. This |
---|
322 | function is allowed to be discontinuous across the edges of the |
---|
323 | cells, but the slope of this function is limited to avoid |
---|
324 | artificially introduced oscillations. |
---|
325 | |
---|
326 | Godunov's method (see \cite{Toro-92}) involves calculating the |
---|
327 | numerical flux function $\HH(\cdot, \cdot ; \ \cdot)$ by exactly |
---|
328 | solving the corresponding one dimensional Riemann problem normal to |
---|
329 | the edge. We use the central-upwind scheme of \cite{KurNP2001} to |
---|
330 | calculate an approximation of the flux across each edge. |
---|
331 | |
---|
332 | \begin{figure} |
---|
333 | \begin{center} |
---|
334 | %\includegraphics[width=5.0cm,keepaspectratio=true]{step-reconstruct} |
---|
335 | \caption{From the values of the conserved quantities at the centroid |
---|
336 | of the cell and its neighbouring cells, a discontinuous piecewise |
---|
337 | linear reconstruction of the conserved quantities is obtained.} |
---|
338 | \label{fig:mesh:reconstruct} |
---|
339 | \end{center} |
---|
340 | \end{figure} |
---|
341 | |
---|
342 | In the computations presented in this paper we use an explicit Euler |
---|
343 | time stepping method with variable timestepping adapted to the |
---|
344 | observed CFL condition. |
---|
345 | |
---|
346 | |
---|
347 | \section{Software Implementation} |
---|
348 | \label{sec:software} |
---|
349 | |
---|
350 | {\bf Jane: I don't think this section is needed here, but the |
---|
351 | content is referred to at the end of section 1} |
---|
352 | |
---|
353 | \ANUGA{} is mostly written in the object-oriented programming |
---|
354 | language \Python{} with computationally intensive parts implemented |
---|
355 | as highly optimised shared objects written in C. |
---|
356 | |
---|
357 | \Python{} is known for its clarity, elegance, efficiency and |
---|
358 | reliability. Complex software can be built in \Python{} without undue |
---|
359 | distractions arising from idiosyncrasies of the underlying software |
---|
360 | language syntax. In addition, \Python{}'s automatic memory management, |
---|
361 | dynamic typing, object model and vast number of libraries means that |
---|
362 | software can be produced quickly and can be readily adapted to |
---|
363 | changing requirements throughout its lifetime. |
---|
364 | |
---|
365 | The fundamental object in \ANUGA{} is the \code{Domain} which |
---|
366 | inherits functionality from a hierarchy of increasingly specialised |
---|
367 | classes starting with a basic structural Mesh to classes implementing |
---|
368 | the finite-volume scheme described in section \ref{sec:fvm}. Other classes |
---|
369 | are \code{Quantity} which represents values of one variable across the mesh |
---|
370 | along with their associated operations, \code{Geospatial_data} which |
---|
371 | represents georeferenced elevation data and a collection of \code{Boundary} |
---|
372 | classes which allows for a 'pluggable' way of driving the model. |
---|
373 | The conserved quantities updated automatically by the numerical scheme |
---|
374 | are stage (water level) $w$, $x$-momentum $uh$ and $y$-momentum |
---|
375 | $vh$. The quanitites elevation $z$ and friction $\eta$ are |
---|
376 | quantities that are not updated automatically but can be changed explicitly |
---|
377 | during run-time if the user wishes to do so. |
---|
378 | |
---|
379 | To set up a scenario the user specifies the study area along with any internal |
---|
380 | regions where increased mesh resolution is required. External edges may |
---|
381 | be labelled using symbolic tags which are subsequently used to bind |
---|
382 | boundary condition objects to tagged segments of the mesh boundary. |
---|
383 | The mesh is then generated using \ANUGA{}'s built-in mesh generator and |
---|
384 | converted into the \code{Domain} object which provides all methods used to |
---|
385 | setup and run the flow simulation. Figure \ref{fig:anuga mesh} shows an example of a mesh generated by \ANUGA{}. |
---|
386 | \begin{figure} |
---|
387 | \begin{center} |
---|
388 | \includegraphics[width=4in,keepaspectratio=true]{triangular-mesh} |
---|
389 | \caption{Triangular mesh used in our finite volume method. Conserved |
---|
390 | quantities $h$, $uh$ and $vh$ are associated with each triangular cell.} |
---|
391 | \label{fig:anuga mesh} |
---|
392 | \end{center} |
---|
393 | \end{figure} |
---|
394 | |
---|
395 | |
---|
396 | |
---|
397 | Next step is to setup initial conditions for each \code{Quantity} object. For |
---|
398 | the elevation $z$ this is typically obtained from bathymetric and |
---|
399 | topographic data sets. Setting initial values for quantities is done |
---|
400 | through the method \code{domain.set_quantity(name, X, location, |
---|
401 | region)} where name is the name of the quantity (e.g.\ 'stage', |
---|
402 | 'xmomentum', 'ymomentum', 'elevation' or 'friction'). The variable X |
---|
403 | represents the source data for populating the quantity and may take |
---|
404 | one of the following forms: |
---|
405 | |
---|
406 | \begin{itemize} |
---|
407 | \item A constant value as in \code{domain.set_quantity('stage', 1)} which |
---|
408 | will set the initial water level to 1 m everywhere. |
---|
409 | \item Another quantity or a linear combination of quantities. If \code{q1} |
---|
410 | and \code{q2} are two arbitrary quantities defined within the same domain, |
---|
411 | the expression \code{domain.set_quantity('stage', q1*(3*q2 + 5))} will set the stage |
---|
412 | quantity accordingly. One common application of this would be to |
---|
413 | assign the stage as a constant depth above the bed elevation. |
---|
414 | \item An arbitrary function (or a callable object), \code{f(x, y)}, where \code{x} |
---|
415 | and \code{y} are assumed to be vectors. The quantity will be assigned values by |
---|
416 | evaluating \code{f} at each location within the mesh. |
---|
417 | \item An arbitrary set of points and associated values (wrapped into a |
---|
418 | Geospatial_data object). The points need not coincide with triangle |
---|
419 | vertices or centroids and a penalised least squares technique is |
---|
420 | employed to populate the quantity in a smooth and stable way. |
---|
421 | Since the least squares technique can be time consuming for large |
---|
422 | problems, \code{set_quantity} employs a caching technique which automatically |
---|
423 | decides whether to perform the computations or retrieve them from a |
---|
424 | cache. This will typically speed up the build by several orders of |
---|
425 | magnitude after each computation has been performed once. |
---|
426 | \item A filename containing points and attributes. |
---|
427 | \item A Numerical Python array (or a list of numbers) ordered |
---|
428 | according to the internal data structure. |
---|
429 | \end{itemize} |
---|
430 | The parameter \code{location} determines whether the values should be |
---|
431 | assigned to triangle edge, midpoints or vertices and \code{region} allows the |
---|
432 | operation to be restricted to a region specified by a symbolic tag or |
---|
433 | a set of indices. |
---|
434 | |
---|
435 | Boundary conditions are bound to symbolic tags through the method |
---|
436 | \code{domain.set_boundary} which takes as input a lookup table (implemented |
---|
437 | as a Python dictionary) of the form \code{\{tag:~boundary_object\}}. |
---|
438 | The boundary objects are all assumed to be callable functions of vectors x |
---|
439 | and y. Several predefined standard boundary objects are available and |
---|
440 | it is relatively straightforward to define problem-specific custom |
---|
441 | boundaries if needed. The predefined boundary conditions include |
---|
442 | Dirichlet, Reflective, Transmissive, Temporal, and Spatio-Temporal |
---|
443 | boundaries. |
---|
444 | |
---|
445 | Forcing terms can be written according to a fixed protocol and added |
---|
446 | to the model using the idiom \code{domain.forcing_terms.append(F)} where \code{F} is |
---|
447 | assumed to be a user-defined callable object. |
---|
448 | |
---|
449 | When the simulation is running, the length of each time step is |
---|
450 | determined from the maximal speeds encountered and the sizes of |
---|
451 | triangles in order not to violate the CFL condition which specifies |
---|
452 | that no information should skip any triangles in one time step. With |
---|
453 | large speeds and small triangles, time steps can become very small. |
---|
454 | In order to access the state of the simulation at regular time |
---|
455 | intervals, \ANUGA{} uses the method evolve: |
---|
456 | \begin{verbatim} |
---|
457 | For t in domain.evolve(yieldstep, duration): |
---|
458 | <model interrogation and modification> |
---|
459 | \end{verbatim} |
---|
460 | The parameter \code{duration} specifies the time period over which |
---|
461 | evolve operates, and control is passed to the body of the for-loop at |
---|
462 | each fixed time step called \code{yieldstep}. The internal workings |
---|
463 | of the numerical scheme and its variable time stepping are thus |
---|
464 | decoupled from the fixed time stepping of the evolve loop. This means |
---|
465 | that the user of the API may access the model at fixed timesteps to |
---|
466 | e.g.\ store model outputs, interrogate quantities or change the model |
---|
467 | itself at runtime. The evolve method has been implemented using a |
---|
468 | Python generator hence the reference to 'yield' in the parameter name. |
---|
469 | |
---|
470 | %Figure \ref{fig:beach runup} shows a simulation of water flowing onto a |
---|
471 | %hypothetical beach with obstacles. |
---|
472 | %A number of complex patterns are captured in this example including a shock where water reflected off the wall far (at the right hand side) meets the main flow. Other physical features are the standing waves and interference patterns. |
---|
473 | %See the \ANUGA{} User Manual at \url{http://sourceforge.net/projects/anuga} for more details and examples. |
---|
474 | %\begin{figure} |
---|
475 | %\begin{center} |
---|
476 | %\includegraphics[width=4in,keepaspectratio=true]{runup} |
---|
477 | %\caption{A hypothetical runup scenario.} |
---|
478 | %\label{fig:beach runup} |
---|
479 | %\end{center} |
---|
480 | %\end{figure} |
---|
481 | |
---|
482 | |
---|
483 | |
---|
484 | \section{Validation} |
---|
485 | \label{sec:validation} Validation is an ongoing process and the purpose of this paper |
---|
486 | is to describe a range of tests that validate \ANUGA{} as a hydrodynamic model. |
---|
487 | This section will describe the six tests outlined in section~\ref{sec:intro}. |
---|
488 | |
---|
489 | \subsection{1D analytical validation} |
---|
490 | |
---|
491 | Tom Baldock has done something here for that NSW report |
---|
492 | |
---|
493 | \subsection{Stage and Velocity Validation in a Flume} |
---|
494 | This section will describe flume tank experiments that were |
---|
495 | conducted at the Gordon McKay Hydraulics Laboratory at the University of |
---|
496 | Queensland that confirm \ANUGA{}'s ability to estimate wave height |
---|
497 | and velocity. The same flume tank simulations were also used |
---|
498 | to explore Manning's friction and this will be described in the next section. |
---|
499 | |
---|
500 | The flume was set up for dam-break experiments, having a |
---|
501 | water reservior at one end. The flume was glass-sided, 3m long, 0.4m |
---|
502 | in wide, and 0.4m deep, with a PVC bottom. The reservoir in the flume |
---|
503 | was 0.75m long. For this experiment the reservoir water was 0.2m |
---|
504 | deep. At time zero the reservoir gate is opened and the water flows |
---|
505 | into the other side of the flume. The water ran up a flume slope of |
---|
506 | 0.03 m/m. To accurately model the bed surface a Manning's friction |
---|
507 | value of 0.01, representing PVC was used. |
---|
508 | |
---|
509 | % Neale, L.C. and R.E. Price. Flow characteristics of PVC sewer pipe. |
---|
510 | % Journal of the Sanitary Engineering Division, Div. Proc 90SA3, ASCE. |
---|
511 | % pp. 109-129. 1964. |
---|
512 | |
---|
513 | Acoustic displacement sensors that produced a voltage that changed |
---|
514 | with the water depth was positioned 0.4m from the reservoir gate. The |
---|
515 | water velocity was measured with an Acoustic Doppler Velocimeter 0.45m |
---|
516 | from the reservoir gate. This sensor only produced reliable results 4 |
---|
517 | seconds after the reservoir gate opened, due to limitations of the sensor. |
---|
518 | |
---|
519 | |
---|
520 | % Validation UQ flume |
---|
521 | % at X:\anuga_validation\uq_sloped_flume_2008 |
---|
522 | % run run_dam.py to create sww file and .csv files |
---|
523 | % run plot.py to create graphs heere automatically |
---|
524 | % The Coasts and Ports '2007 paper is in TRIM d2007-17186 |
---|
525 | \begin{figure}[htbp] |
---|
526 | \centerline{\includegraphics[width=4in]{uq-flume-depth}} |
---|
527 | \caption{Comparison of wave tank and \ANUGA{} water height at .4 m |
---|
528 | from the gate}\label{fig:uq-flume-depth} |
---|
529 | \end{figure} |
---|
530 | |
---|
531 | \begin{figure}[htbp] |
---|
532 | \centerline{\includegraphics[width=4in]{uq-flume-velocity}} |
---|
533 | \caption{Comparison of wave tank and \ANUGA{} water velocity at .45 m |
---|
534 | from the gate}\label{fig:uq-flume-velocity} |
---|
535 | \end{figure} |
---|
536 | |
---|
537 | Figure~\ref{fig:uq-flume-depth} shows that ANUGA predicts the actual |
---|
538 | water depth very well, with the exception of the fluid tip-region {\bf |
---|
539 | Duncan - what does that mean? About where on the graph is that). |
---|
540 | Water depth and velocity are coupled as described by the nonlinear shallow water equations, thus |
---|
541 | if one of these quantities accurately estimates the measured values, we would expect |
---|
542 | the same for the other quantity. This is demonstrated in figure~\ref{fig:uq-flume-velocity) |
---|
543 | where the water velocity is also predicted accurately. Sediment transport studies |
---|
544 | rely on water velocity estimates in the region where the sensors cannot provide this data. |
---|
545 | With water velocity being accurately predicted, studies such as sediment transport can now use |
---|
546 | reliable estimates. |
---|
547 | |
---|
548 | |
---|
549 | \subsection{1D flume tank to verify friction} |
---|
550 | |
---|
551 | The same flume tank experimental setup was used to obtain friction values for |
---|
552 | use in hydrodynamic models. A number of bed friction scenarios were simulated in |
---|
553 | the flume tank. The PVC bottom of the tank is equivalent to a friction value of 0 (i.e |
---|
554 | completely smooth) and small pebbles were used to cover the base of the tank and the |
---|
555 | aim of the experiment was to determine what the Manning's friction value is for |
---|
556 | this case. |
---|
557 | |
---|
558 | As described in the model equations in \section~\ref{sec:model}, the bed |
---|
559 | friction is modelled using the Manning's model. |
---|
560 | Validation of this model was carried out by comparing results |
---|
561 | from ANUGA against experimental results from flume wave tanks. |
---|
562 | |
---|
563 | % Validation UQ friction |
---|
564 | % at X:\anuga_validation\uq_friction_2007 |
---|
565 | % run run_dam.py to create sww file and .csv files |
---|
566 | % run plot.py to create graphs, and move them here |
---|
567 | \begin{figure}[htbp] |
---|
568 | \centerline{\includegraphics[width=4in]{uq-friction-depth}} |
---|
569 | \caption{Comparison of wave tank and \ANUGA{} water height at .4 m |
---|
570 | from the gate, simulated using a Mannings friction of 0.0 and 0.1.}\label{fig:uq-friction-depth} |
---|
571 | \end{figure} |
---|
572 | |
---|
573 | \subsection{Okushiri Wavetank Validation} |
---|
574 | |
---|
575 | As part of the Third International Workshop on Long-wave Runup |
---|
576 | Models in 2004 (\url{http://www.cee.cornell.edu/longwave}), four |
---|
577 | benchmark problems were specified to allow the comparison of |
---|
578 | numerical, analytical and physical models with laboratory and field |
---|
579 | data. One of these problems describes a wave tank simulation of the |
---|
580 | 1993 Okushiri Island tsunami off Hokkaido, Japan \cite{MatH2001}. A |
---|
581 | significant feature of this tsunami was a maximum run-up of 32~m |
---|
582 | observed at the head of the Monai Valley. This run-up was not |
---|
583 | uniform along the coast and is thought to have resulted from a |
---|
584 | particular topographic effect. Among other features, simulations of |
---|
585 | the Hokkaido tsunami should capture this run-up phenomenon. |
---|
586 | |
---|
587 | This dataset has been used by to validate tsunami models by |
---|
588 | a number of tsunami scientists. Examples include Titov ... lit review |
---|
589 | here on who has used this example for verification |
---|
590 | |
---|
591 | \begin{figure}[htbp] |
---|
592 | %\centerline{\includegraphics[width=4in]{okushiri-gauge-5.eps}} |
---|
593 | \centerline{\includegraphics[width=4in]{ch5.png}} |
---|
594 | \centerline{\includegraphics[width=4in]{ch7.png}} |
---|
595 | \centerline{\includegraphics[width=4in]{ch9.png}} |
---|
596 | \caption{Comparison of wave tank and \ANUGA{} water stages at gauge |
---|
597 | 5,7 and 9.}\label{fig:val} |
---|
598 | \end{figure} |
---|
599 | |
---|
600 | |
---|
601 | \begin{figure}[htbp] |
---|
602 | \centerline{\includegraphics[width=4in]{okushiri-model.jpg}} |
---|
603 | \caption{Complex reflection patterns and run-up into Monai Valley |
---|
604 | simulated by \ANUGA{} and visualised using our netcdf OSG |
---|
605 | viewer.}\label{fig:run} |
---|
606 | \end{figure} |
---|
607 | |
---|
608 | The wave tank simulation of the Hokkaido tsunami was used as the |
---|
609 | first scenario for validating \ANUGA{}. The dataset provided |
---|
610 | bathymetry and topography along with initial water depth and the |
---|
611 | wave specifications. The dataset also contained water depth time |
---|
612 | series from three wave gauges situated offshore from the simulated |
---|
613 | inundation area. The \ANUGA{} model comprised $41404$ triangles |
---|
614 | and took about $2000$ s to run on a standard PC or $1500$ s |
---|
615 | on a 64-bit Opteron 2000 series Linux server. |
---|
616 | |
---|
617 | Figure~\ref{fig:val} compares the observed wave tank and modelled |
---|
618 | \ANUGA{} water depth (stage height) at one of the gauges. The plots |
---|
619 | show good agreement between the two time series, with \ANUGA{ |
---|
620 | closely modelling the initial draw down, the wave shoulder and the |
---|
621 | subsequent reflections. The discrepancy between modelled and |
---|
622 | simulated data in the first 10 seconds is due to the initial |
---|
623 | condition in the physical tank not being uniformly zero. Similarly |
---|
624 | good comparisons are evident with data from the other two gauges. |
---|
625 | Additionally, \ANUGA{} replicates exceptionally well the 32~m Monai |
---|
626 | Valley run-up, and demonstrates its occurrence to be due to the |
---|
627 | interaction of the tsunami wave with two juxtaposed valleys above |
---|
628 | the coastline. The run-up is depicted in Figure~\ref{fig:run}. |
---|
629 | |
---|
630 | This successful replication of the tsunami wave tank simulation on a |
---|
631 | complex 3D beach is a positive first step in validating the \ANUGA{} |
---|
632 | modelling capability. |
---|
633 | |
---|
634 | \subsection{Runup of solitary wave on circular island wavetank validation} |
---|
635 | |
---|
636 | This section will describe the ANUGA results for the experiments conducted |
---|
637 | by Briggs et al (1995). Here, a 30x25m basin with a conical island is situated near |
---|
638 | the centre and a directional wavemaker is used to produce planar solitary waves of |
---|
639 | specified crest lenghts and heights. A series of gauges were distributed within the |
---|
640 | experimental setup. As described by Hubbard and Dodd \cite{Hubbard02}, a number of researchers |
---|
641 | have used this benchmark problem to test their numerical models. {\bf Jane: check |
---|
642 | whether these results are now avilable as they were not in 2002}. Hubbard and Dodd |
---|
643 | \cite{Hubbard02} note that |
---|
644 | a particular 3D model appears to obtain slightly better results than the 2D ones reported |
---|
645 | but that 3D models are unlikely to be competitive in terms of computing power for |
---|
646 | applications in coastal engineering at least. Choi et al \cite{Choi07) use a 3D RANS model |
---|
647 | (based on the Navier-Stokes equations) |
---|
648 | for the same problem and find a very good comparison with laboratory and 2D numerical |
---|
649 | results. An obvious advantage of the 3D model is its ability to investigate the |
---|
650 | velocity field and Choi et al also report on the limitation of depth-averaged |
---|
651 | 2D models for run-up simulations of this type. |
---|
652 | |
---|
653 | Once results are availble, need to compare to Hubbard and Dodd and draw any conclusions |
---|
654 | from nested rectangular grid vs unstructured gird. |
---|
655 | Figure \ref{fig:circular screenshots} shows a sequence of screenshots depicting the evolution of the solitary wave as it hits the circular island. |
---|
656 | |
---|
657 | \begin{figure}[htbp] |
---|
658 | \centerline{ |
---|
659 | \includegraphics[width=5cm]{circular1.png} |
---|
660 | \includegraphics[width=5cm]{circular2.png}} |
---|
661 | \centerline{ |
---|
662 | \includegraphics[width=5cm]{circular3.png} |
---|
663 | \includegraphics[width=5cm]{circular4.png}} |
---|
664 | \centerline{ |
---|
665 | \includegraphics[width=5cm]{circular5.png} |
---|
666 | \includegraphics[width=5cm]{circular6.png}} |
---|
667 | \centerline{ |
---|
668 | \includegraphics[width=5cm]{circular7.png} |
---|
669 | \includegraphics[width=5cm]{circular8.png}} |
---|
670 | \centerline{ |
---|
671 | \includegraphics[width=5cm]{circular9.png} |
---|
672 | \includegraphics[width=5cm]{circular10.png}} |
---|
673 | \caption{Screenshots of the evolution of solitary wave around circular island.} |
---|
674 | \label{fig:circular screenshots} |
---|
675 | \end{figure} |
---|
676 | |
---|
677 | |
---|
678 | \clearpage |
---|
679 | |
---|
680 | \section{Conclusions} |
---|
681 | \label{sec:6} |
---|
682 | \ANUGA{} is a flexible and robust modelling system |
---|
683 | that simulates hydrodynamics by solving the shallow water wave |
---|
684 | equation in a triangular mesh. It can model the process of wetting |
---|
685 | and drying as water enters and leaves an area and is capable of |
---|
686 | capturing hydraulic shocks due to the ability of the finite-volume |
---|
687 | method to accommodate discontinuities in the solution. |
---|
688 | \ANUGA{} can take as input bathymetric and topographic datasets and |
---|
689 | simulate the behaviour of riverine flooding, storm surge, |
---|
690 | tsunami or even dam breaks. |
---|
691 | Initial validation using wave tank data supports \ANUGA{}'s |
---|
692 | ability to model complex scenarios. Further validation will be |
---|
693 | pursued as additional datasets become available. |
---|
694 | The \ANUGA{} source code is available |
---|
695 | at \url{http://sourceforge.net/projects/anuga}. |
---|
696 | |
---|
697 | something about use on flood modelling community and their validation initiatives |
---|
698 | |
---|
699 | \bibliographystyle{plain} |
---|
700 | \bibliography{anuga-bibliography} |
---|
701 | |
---|
702 | \end{document} |
---|