Changeset 2463


Ignore:
Timestamp:
Mar 1, 2006, 4:58:28 PM (19 years ago)
Author:
howard
Message:

Further updates to Chapter 2 examples. Other minor changes.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • documentation/user_manual/anuga_user_manual.tex

    r2457 r2463  
    213213that describes each step of the program and explains it significance.
    214214
    215 
    216 {\scriptsize \begin{verbatim}
    217 from pyvolution.mesh_factory import rectangular
    218 from pyvolution.shallow_water import Domain, Reflective_boundary,
    219      Dirichlet_boundary, Time_boundary, Transmissive_boundary
    220 
    221 #Create basic mesh
    222 points, vertices, boundary = rectangular(10,10)
    223 
    224 #Create shallow water domain
    225 domain = Domain(points, vertices,boundary)
    226 domain.set_name('bedslope')
    227 
    228 
    229 #######################
    230 # Initial conditions
    231 def f(x,y):
    232     return -x/2
    233 
    234 domain.set_quantity('elevation', f)
    235 domain.set_quantity('friction', 0.1)
    236 
    237 h = 0.05  # Constant depth
    238 domain.set_quantity('stage', expression = 'elevation + %f' %h)
    239 
    240 
    241 # Boundary conditions
    242 from math import sin, pi
    243 Br = Reflective_boundary(domain)
    244 Bt = Transmissive_boundary(domain)
    245 Bd = Dirichlet_boundary([0.2,0.,0.])
    246 
    247 Bw = Time_boundary(domain=domain,
    248                    f=lambda t: [(0.1*sin(t*2*pi)), 0.0, 0.0])
    249 
    250 
    251 domain.set_boundary({'left': Bd, 'right': Br, 'top': Br, 'bottom': Br})
    252 
    253 
    254 ######################
    255 #Evolution
    256 
    257 domain.check_integrity()
    258 
    259 for t in domain.evolve(yieldstep = 0.1, finaltime = 4.0):
    260     domain.write_time()
    261 \end{verbatim}}
    262 
     215\verbatiminput{../../inundation/examples/bedslope.py}
    263216
    264217\subsection{Establishing the Mesh}
    … …  
    272225
    273226The function \code{rectangular} is imported from a module
    274 \code{mesh\_factory} defined elsewhere. \anuga also
    275 contains several other schemes that can be used for setting up
    276 meshes, but we shall not discuss these now.) The above assignment
    277 sets up a $10 \times 10$ rectangular mesh, triangulated in a
    278 specific way. In general, the assignment
     227\code{mesh\_factory} defined elsewhere. (\anuga also contains
     228several other schemes that can be used for setting up meshes, but we
     229shall not discuss these now.) The above assignment sets up a $10
     230\times 10$ rectangular mesh, triangulated in a specific way. In
     231general, the assignment
    279232
    280233{\small \begin{verbatim}
    … …  
    286239\begin{itemize}
    287240
    288    \item a list \code{points} % of length $N$, where $N = (m + 1)(n + 1)$,
    289 %comprising the coordinates $(x, y)$ of each of the $N$ mesh
    290 %points,
    291 
    292    \item a list \code{vertices} %of length $2mn$ (each entry specifies the three
    293 %vertices of one of the triangles used in the triangulation) , and
     241   \item a list \code{points} of length $N$, where $N = (m + 1)(n + 1)$,
     242comprising the coordinates $(x, y)$ of each of the $N$ mesh points,
     243
     244   \item a list \code{vertices} of length $2mn$ (each entry specifies the three
     245vertices of one of the triangles used in the triangulation) , and
    294246
    295247   \item a dictionary \code{boundary}, used to tag the triangle edges on
    … …  
    318270are set at this point. One of them is to set the basename for the output file:
    319271
    320 {\scriptsize \begin{verbatim}
     272{\small \begin{verbatim}
    321273    domain.set_name('bedslope')
    322274\end{verbatim}}
    … …  
    342294\subsubsection{Elevation}
    343295
    344 The elevation is set using a function, defined through the
    345 statements below, which is specific to this example and specifies
    346 a particularly simple initial configuration for demonstration
    347 purposes:
     296The elevation, or height of the bed, is set using a function,
     297defined through the statements below, which is specific to this
     298example and specifies a particularly simple initial configuration
     299for demonstration purposes:
    348300
    349301{\small \begin{verbatim}
    … …  
    378330to 0.1 at every mesh point.
    379331
    380 \subsubsection{Depth}
    381 
    382 Assigning depth illustrates a more complex way to use
    383 \code{set\_quantity}, introducing an expression involving other
    384 quantities:
     332\subsubsection{Stage}
     333
     334The stage (the height of the water surface) is related to the
     335elevation and the depth at any time by the equation \[\code{stage} =
     336\code{elevation} + \code{depth}\] To specify a constant depth of a
     337particular value, therefore, we need to specify that \code{stage} is
     338everywhere obtained by adding that constant value to the already
     339specified \code{elevation}. Doing this allows us to illustrate
     340another way to use \code{set\_quantity}, introducing an expression
     341involving other quantities:
    385342
    386343{\small \begin{verbatim}
    … …  
    389346\end{verbatim}}
    390347
    391 Here the quantity \code{stage} is defined by taking the quantity
    392 elevation already defined and adding a constant value $h = 0.05$
    393 to it everywhere. This expresses the fact that the water depth is
    394 everywhere constant, so the surface is a constant height above the
    395 elevation of the bed.
     348Here we are stipulating that the quantity \code{stage} is defined by
     349taking the quantity elevation already defined and adding a constant
     350value $h = 0.05$ to it everywhere.
    396351
    397352\subsubsection{Boundary Conditions}
    … …  
    477432
    478433\begin{itemize}
    479 \item{from a Windows command line} as in \code{python bedslope.py}
    480 
    481 \item{within the Python IDLE environment}
    482 
    483 \item{within emacs}
    484 
    485 \item{from a Linux command line} as in \code{python bedslope.py}
     434  \item{from a Windows command line} as in \code{python bedslope.py}
     435  \item{within the Python IDLE environment}
     436  \item{within emacs}
     437  \item{from a Linux command line} as in \code{python bedslope.py}
    486438\end{itemize}
    487439
     440%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    488441
    489442\section{An Example with Real Data}
    … …  
    491444The following discussion builds on the concepts introduced through
    492445the \code{bedslope.py} example and introduces a second example,
    493 \code{run_sydney_smf.py}, that follows the same basic outline, but
     446\code{run\_sydney\_smf.py}, that follows the same basic outline, but
    494447incorporates more complex features and refers to a real-life
    495448scenario, rather than the artificial illustrative one used in
    496 \code{bedslope.py}.
     449\code{bedslope.py}. [Include details of the scenario to which it
     450refers??]
    497451
    498452\subsection{Overview}
    … …  
    523477\subsection{The Code}
    524478
    525 Here is the code for \code{run_sydney_smf.py}:
    526 
    527 {\scriptsize \begin{verbatim}
    528 """Script for running a tsunami
    529 inundation scenario for Sydney, NSW, Australia.
    530 
    531 Source data such as elevation and boundary data is assumed to be
    532 available in directories specified by project.py The output sww file
    533 is stored in project.outputdir
    534 
    535 The scenario is defined by a triangular mesh created from
    536 project.polygon, the elevation data and a simulated submarine
    537 landslide.
    538 
    539 Ole Nielsen and Duncan Gray, GA - 2005 and Adrian Hitchman and Jane
    540 Sexton, GA - 2006 """
    541 
    542 
    543 #-------------------------------------------------------------------------------#
    544 Import necessary modules
    545 #-------------------------------------------------------------------------------
    546 
    547 # Standard modules
    548 import os
    549 import time
    550 
    551 # Related major packages
    552 from pyvolution.shallow_water import
    553 Domain, Reflective_boundary
    554 from pyvolution.data_manager import
    555 convert_dem_from_ascii2netcdf, dem2pts
    556 from pyvolution.combine_pts
    557 import combine_rectangular_points_files
    558 from pyvolution.pmesh2domain
    559 import pmesh_to_domain_instance
    560 
    561 # Application specific imports
    562 import project                 #
    563 Definition of file names and polygons
    564 from smf import slump_tsunami
    565 # Function for submarine mudslide
    566 \end{verbatim}}
    567 
    568 
    569 {\scriptsize \begin{verbatim}
    570 #-------------------------------------------------------------------------------
    571 # Preparation of topographic data # # Convert ASC 2 DEM 2 PTS using
    572 source data and store result in source data # Do for coarse and fine
    573 data
    574 
    575 # Fine pts file to be clipped to area of interest
    576 #-------------------------------------------------------------------------------
    577 
    578 # filenames
    579 coarsedemname = project.coarsedemname
    580 finedemname =
    581 project.finedemname meshname = project.meshname+'.msh'
    582 
    583 # coarse data convert_dem_
    584 from_ascii2netcdf(coarsedemname,
    585 use_cache=True, verbose=True) dem2pts(coarsedemname, use_cache=True,
    586 verbose=True)
    587 
    588 # fine data (clipping the points file to smaller area)
    589 convert_dem_from_ascii2netcdf(finedemname, use_cache=True,
    590 verbose=True) dem2pts(finedemname,
    591         easting_min=project.eastingmin,
    592         easting_max=project.eastingmax,
    593         northing_min=project.northingmin,
    594         northing_max= project.northingmax,
    595         use_cache=True,
    596         verbose=True)
    597 
    598 
    599 # combining the coarse and fine data
    600 combine_rectangular_points_files(project.finedemname + '.pts',
    601                                  project.coarsedemname + '.pts',
    602                                  project.combineddemname + '.pts')
    603 
    604 
    605 #-------------------------------------------------------------------------------
    606 # Create the triangular mesh based on overall clipping polygon with
    607 a tagged # boundary and interior regions defined in project.py along
    608 with # resolutions (maximal area of per triangle) for each polygon
    609 #-------------------------------------------------------------------------------
    610 
    611 interior_regions = [[project.harbour_polygon_2, interior_res],
    612 [project.botanybay_polygon_2, interior_res]]
    613 
    614 #FIXME: Fix caching of this one once the mesh_interface is ready
    615 from caching import cache
    616 _ = cache(create_mesh_from_regions,
    617           project.diffpolygonall,
    618           {'boundary_tags': {'bottom': [0],
    619                              'right1': [1], 'right0': [2],
    620                              'right2': [3], 'top': [4], 'left1': [5],
    621                              'left2': [6], 'left3': [7]},
    622            'resolution': 100000,
    623            'filename': meshname,
    624            'interior_regions': interior_regions,
    625            'UTM': True,
    626            'refzone': project.refzone},
    627           verbose = True)
    628 
    629 
    630 #-------------------------------------------------------------------------------
    631 # Set up computational domain
    632 #-------------------------------------------------------------------------------
    633 
    634 domain = pmesh_to_domain_instance(meshname, Domain,
    635                                   use_cache = True,
    636                                   verbose = True)
    637 
    638 print 'Number of triangles = ', len(domain) print 'The extent is ',
    639 domain.get_extent()
    640 
    641 domain.set_name(project.basename)
     479Here is the code for \code{run\_sydney\_smf.py}:
     480
     481\verbatiminput{"runsydneysmf.py"}
     482
     483In discussing the details of this example, we follow the outline
     484given above, discussing each major step of the code in turn.
     485
     486\subsection{Establishing the Mesh}
     487
     488One obvious way that the present example differs from
     489\code{bedslope.py} is in the use of a more complex method to create
     490the mesh. Instead of imposing a mesh structure on a rectangular
     491grid, the technique used for this example involves building mesh
     492structures inside polygons specified by the user, using a
     493mesh-generator referred to as \code{pmesh}.
     494
     495The following remarks may help the reader understand how
     496\code{pmesh} is used.
     497
     498In its simplest form, \code{pmesh} creates the mesh within a single
     499polygon whose vertices are at geographical locations specified by
     500the user. The user specifies the \emph{resolution}---that is, the
     501maximal area of a triangle used for triangulation---and mesh points
     502are created inside the polygon through a random process. Figure XXX
     503shows a simple example of this, in which the triangulation is
     504carried out within a pentagon.
     505
     506Boundary tags are not restricted to \code{`left'}, \code{`right'},
     507\code{`bottom'} and \code{`top'}, as in the case of
     508\code{bedslope.py}. Instead the user specifies a list of tags
     509appropriate to the configuration being modelled.
     510
     511While a mesh created inside a polygon offers more flexibility than
     512one based on a rectangular grid, using \code{pmesh} in the limited
     513form we have described so far still doesn't provide a way to adapt
     514to geographic or other features in the landscape, whose presence may
     515require us to vary the resolution in the neighbourhood of the
     516features. However, this capability is in fact also provided in
     517\code{pmesh}. It provides for these more general situations by
     518allowing the user to specify a number of \emph{interior polygons},
     519which are triangulated separately, each according to a separately
     520specified resolution. See Figure XXX.
     521
     522In its general form, \code{pmesh} takes for its input a bounding
     523polygon and (optionally) a list of interior polygons. The user
     524specifies resolutions, both for the bounding polygon and for each of
     525the interior polygons. Given this data, \code{pmesh} first creates a
     526triangular mesh inside the bounding polygon, using the specified
     527resolution, and then creates a separate triangulation inside each of
     528the interior polygons, using the resolution specified for that
     529polygon.
     530
     531The function used to implement this process is
     532\code{create\_mesh\_from\_regions}. Its arguments include the
     533bounding polygon and its resolution, a list of boundary tags, and a
     534list of pairs \code{[polygon, resolution]}, specifying the interior
     535polygons and their resolutions.
     536
     537In practice, the details of the polygons used are read from a
     538separate file \code{project.py}. The resulting mesh is output to a
     539\emph{meshfile}\index{meshfile}. This term is used to describe a
     540file of a specific format used to store the data specifying a mesh.
     541(There are in fact two possible formats for such a file: it can
     542either be a binary file, with extension \code{.msh}, or an ASCII
     543file, with extension \code{.tsh}. In the present case, the binary
     544file format \code{.msh} is used. See Section \ref{sec:file formats}
     545for more on file formats.) \code{pmesh} assigns a name to the file
     546by appending the extension \code{.msh} to the name specified in the
     547input file \code{project.py}. This name is stored in the variable
     548\code{meshname}.
     549
     550
     551\subsection{Initialising the Domain}
     552
     553As with \code{bedslope.py}, once we have created the mesh, the next
     554step is to create the data structure \code{domain}. We did this for
     555\code{bedslope.py} by inputting lists of points and triangles and
     556specifying the boundary tags directly. However, in the present case,
     557we use a method that works directly with the meshfile
     558\code{meshname}, as follows:
     559
     560{\small \begin{verbatim}
     561domain = pmesh_to_domain_instance(meshname,
     562Domain, use_cache = True, verbose = True)
     563\end{verbatim}}
     564
     565The function \code{pmesh\_to\_domain\_instance} converts a meshfile
     566\code{meshname} into an instance of the data structure
     567\code{domain}, allowing us to use methods like \code{set\_quantity}
     568to set quantities and to apply other operations. (In principle, the
     569second argument of \code{pmesh\_to\_domain\_instance} can be any
     570subclass of \code{Domain}, but for applications involving the
     571shallow-water wave equation, the second argument of
     572\code{pmesh\_to\_domain\_instance} can always be set simply to
     573\code{Domain}.)
     574
     575\subsection{Specifying the Quantities}
     576Setting quantities for \code{run\_sydney\_smf.py} is accomplished
     577using similar methods to those used in \code{bedslope.py}. However,
     578in this case, many of the values are read from the auxiliary file
     579\code{project.py} or, in the case of \code{elevation}, from an
     580ancillary points file.
     581
     582
     583\subsubsection{Fundamental Quantities}
     584
     585The following statements specify a basename and data directory, and
     586identify quantities to be stored. For the first two, values are
     587taken from \code{project.py}.
     588
     589{\small \begin{verbatim} domain.set_name(project.basename)
    642590domain.set_datadir(project.outputdir)
    643591domain.set_quantities_to_be_stored(['stage', 'xmomentum',
    644592'ymomentum'])
    645 
    646 
     593\end{verbatim}}
     594
     595\subsubsection{Stage}
     596
     597For the scenario we are modelling in this case, we use a callable
     598object \code{tsunami\_source}, assigned by means of a function
     599\code{slump\_tsunami}. This is similar to how we set elevation in
     600\code{bedslope.py} using a function---however, in this case the
     601function is both more complex and more interesting.
     602
     603{\small \begin{verbatim}
    647604#-------------------------------------------------------------------------------
    648605# Set up scenario (tsunami_source is a callable object used with
    … …  
    663620
    664621#-------------------------------------------------------------------------------
    665 # Setup initial conditions
     622# Set up initial conditions
    666623#-------------------------------------------------------------------------------
    667624
    668625domain.set_quantity('stage', tsunami_source)
    669 domain.set_quantity('friction', 0.03) # supplied by Benfield
     626\end{verbatim}}
     627
     628\subsubsection{Friction}
     629
     630Assigning the friction is exactly parallel to what we did for
     631\code{bedslope.py}:
     632
     633{\small \begin{verbatim}
     634domain.set_quantity('friction', 0.03)
     635#supplied by Benfield
     636\end{verbatim}}
     637
     638
     639\subsubsection{Elevation}
     640
     641In this example, the elevation is specified by reading data from a
     642file.
     643
     644{\small \begin{verbatim}
    670645domain.set_quantity('elevation',
    671646                    filename = project.combineddemname + '.pts',
    672647                    use_cache = True,
    673648                    verbose = True)
    674 
    675 
    676 #-------------------------------------------------------------------------------
    677 # Setup boundary conditions (all reflective)
    678 #-------------------------------------------------------------------------------
    679 
    680 print 'Available boundary tags', domain.get_boundary_tags()
    681 
    682 Br = Reflective_boundary(domain) domain.set_boundary( {'bottom': Br,
    683 'right1': Br, 'right0': Br,
    684                       'right2': Br, 'top': Br, 'left1': Br,
    685                       'left2': Br, 'left3': Br} )
    686 
    687 
    688 #-------------------------------------------------------------------------------
    689 # Evolve system through time
    690 #-------------------------------------------------------------------------------
    691 
    692 import time t0 = time.time()
    693 
    694 for t in domain.evolve(yieldstep = 120, finaltime = 18000):
    695     domain.write_time()
    696     domain.write_boundary_statistics(tags = 'bottom')
    697 
    698 print 'That took %.2f seconds' %(time.time()-t0)
    699 
    700 \end{verbatim}}
    701 
    702 Before discussing the details of this example, let us take a more
    703 high-level perspective of the various tasks it undertakes.
    704 
    705 \subsection{Establishing the Mesh}
    706 
    707 Comparing the present example with \code{bedslope.py}, we note that
    708 one important difference is the use of a more complex method to
    709 create the mesh. Instead of imposing a mesh structure on a
    710 rectangular grid, the technique used for this example involves
    711 building mesh structures inside polygons specified by the user.
    712 
    713 The following remarks may help the reader understand the basic ideas
    714 behind the mesh-creation technique used in this program.
    715 
    716 In its simplest form, the technique creates the mesh within a single
    717 polygon whose vertices are at geographical locations specified by
    718 the user. The user specifies the \emph{resolution}---that is, the
    719 maximal area of a triangle used for triangulation---and mesh points
    720 are created inside the polygon through a random process. Figure XXX
    721 shows a simple example of this, in which the triangulation is
    722 carried out within a pentagon.
    723 
    724 Boundary tags are not restricted to \code{`left'}, \code{`right'},
    725 \code{`bottom'} and \code{`top'}, as in the case of
    726 \code{bedslope.py}. Instead the user specifies a list of tags
    727 appropriate to the configuration being modelled.
    728 
    729 While a mesh created inside a polygon offers more flexibility than
    730 one based on a rectangular grid, the technique of mesh creation that
    731 we have so far described doesn't provide a way to adapt to
    732 geographic or other features in the landscape, whose presence may
    733 require us to vary the resolution in the neighbourhood of the
    734 features. To achieve more flexibility, an extension of the technique
    735 is adopted, which allows the user to specify a number of interior
    736 polygons to be triangulated separately, possibly using different
    737 resolutions. See Figure XXX.
    738 
    739 The upshot is a general method that takes for its input a bounding
    740 polygon and (optionally) a list of interior polygons. It creates a
    741 triangular mesh inside the bounding polygon, using a user-specified
    742 resolution, but then creates a separate triangulation inside each of
    743 the interior polygons, according to the resolution specified for
    744 that polygon.
    745 
    746 The function used to implement this process is
    747 \code{create_mesh_from_regions}. Its arguments include the bounding
    748 polygon and its resolution, a list of boundary tags, and a list of
    749 pairs \code{[polygon, resolution]}, specifying the interior polygons
    750 and their resolutions.
    751 
    752 In practice, the details of the polygons used are read from a
    753 separate file \code{project.py}. The resulting mesh is output to a
    754 \emph{meshfile}, which is a file of a specific format used to store
    755 the data specifying a mesh. (There are two formats for such a file:
    756 it can either be a binary file, with extension \code{.msh}, or an
    757 ASCII file, with extension \code{.tsh}. In the present case, the
    758 binary file format \code{.msh} is used. See Section \ref{sec:file
    759 formats} for more on file formats.) and assigned a name obtained by
    760 appending the extension \code{.msh} to the name specified in the
    761 input file \code{project.py}. This name is stored in the variable
    762 \code{meshname}.
    763 
    764 
    765 \subsection{Initialising the Domain}
    766 
    767 As with \code{bedslope.py}, once we have created the mesh, the next
    768 step is to create the data structure \code{domain}. This is
    769 accomplished through the following statements:
    770 
    771 {\scriptsize \begin{verbatim}
    772 domain =
    773 pmesh_to_domain_instance(meshname, Domain,
    774                                   use_cache = True,
    775                                   verbose = True)
    776 \end{verbatim}}
    777 
    778 The function \code{pmesh_to_domain_instance} converts a meshfile
    779 \code{meshname} into an instance of the data structure
    780 \code{domain}, allowing us to use methods like \code{set\_quantity}
    781 to set quantities and to apply other operations. In the case of
    782 \code{bedslope.py},
    783 
    784 \subsection{Specifying the Quantities}
    785 
    786 
    787 
    788 {\scriptsize \begin{verbatim}
    789 domain.set_name(project.basename)
    790 domain.set_datadir(project.outputdir)
    791 domain.set_quantities_to_be_stored(['stage', 'xmomentum',
    792 'ymomentum'])
    793 
    794 
    795 #-------------------------------------------------------------------------------
    796 # Set up scenario (tsunami_source is a callable object used with
    797 set_quantity)
    798 #-------------------------------------------------------------------------------
    799 
    800 tsunami_source = slump_tsunami(length=30000.0,
    801                                depth=400.0,
    802                                slope=6.0,
    803                                thickness=176.0,
    804                                radius=3330,
    805                                dphi=0.23,
    806                                x0=project.slump_origin[0],
    807                                y0=project.slump_origin[1],
    808                                alpha=0.0,
    809                                domain=domain)
    810 
    811 
    812 #-------------------------------------------------------------------------------
    813 # Set up initial conditions
    814 #-------------------------------------------------------------------------------
    815 
    816 domain.set_quantity('stage', tsunami_source)
    817 domain.set_quantity('friction', 0.03) # supplied by Benfield
    818 domain.set_quantity('elevation',
    819                     filename = project.combineddemname + '.pts',
    820                     use_cache = True,
    821                     verbose = True)
    822 \end{verbatim}}
     649\end{verbatim}}
     650
     651However, before this step can be executed, some preliminary steps
     652are needed to prepare the file from which the data is taken. Two
     653source files are used for this data---their names are specified in
     654the file \code{project.py}, in the variables \code{coarsedemname}
     655and \code{finedemname}. They contain `coarse' and `fine' data,
     656respectively---that is, data sampled at widely spaced points over a
     657large region and data sampled at closely spaced points over a
     658smaller subregion. The data in these files is combined through the
     659statement
     660
     661{\small \begin{verbatim}
     662combine_rectangular_points_files(project.finedemname + '.pts',
     663                                 project.coarsedemname + '.pts',
     664                                 project.combineddemname + '.pts')
     665\end{verbatim}}
     666
     667The effect of this is simply to combine the datasets by eliminating
     668any coarse data associated with points inside the smaller region
     669common to both datasets. The name to be assigned to the resulting
     670dataset is also derived from the name stored in the variable
     671\code{combinedname} in the file \code{project.py}.
    823672
    824673\subsection{Boundary Conditions}
    825674
    826 {\scriptsize \begin{verbatim}
     675Setting boundaries in \code{run\_sydney\_smf.py} is very similar to
     676what we did in \code{bedslope.py}, except that we now have a larger
     677number of boundary tags:
     678
     679{\small \begin{verbatim}
    827680Br = Reflective_boundary(domain)
    828681domain.set_boundary( {'bottom': Br, 'right1': Br, 'right0': Br,
    … …  
    831684\end{verbatim}}
    832685
     686\subsection{Evolution}
     687
     688With the basics established, the running of the `evolve' step is
     689very similar to the corresponding step in \code{bedslope.py}:
     690
     691{\small \begin{verbatim}
     692import time t0 = time.time()
     693
     694for t in domain.evolve(yieldstep = 120, finaltime = 18000):
     695    domain.write_time()
     696    domain.write_boundary_statistics(tags = 'bottom')
     697
     698print 'That took %.2f seconds' %(time.time()-t0)
     699\end{verbatim}}
     700
    833701\chapter{\anuga Public Interface}
    834702
    … …  
    837705\section{Functions and Classes}
    838706
    839 \indexedcodeheader{create_mesh_from_region}
     707\indexedcodeheader{create\_mesh\_from\_region}\label{codehdr:create\_mesh\_from\_region}
    840708
    841709  Creates a triangular mesh based on a bounding polygon and
    … …  
    856724
    857725
    858 \indexedcodeheader{pmesh_to_domain_instance}
     726\indexedcodeheader{pmesh\_to\_domain\_instance}
    859727
    860728 Converts a generated mesh file to a domain object.
    … …  
    863731
    864732  \begin{itemize}
    865     \item \code{file_name} is the name of the mesh file to convert, including the extension
     733    \item \code{file\_name} is the name of the mesh file to convert, including the extension
    866734    \item \code{DomainClass} is the Class that will be returned.
    867735    It must be a subclass of \code{Domain}, with the same interface as domain.
    868     \item \code{use_cache}: \code{True} means that caching is attempted for the computed domain.
     736    \item \code{use\_cache}: \code{True} means that caching is attempted for the computed domain.
    869737  \end{itemize}
    870738
    … …  
    875743  \end{itemize}
    876744
    877 \indexedcodeheader{file_function} %in util.py "High priority"
     745\indexedcodeheader{file\_function} %in util.py "High priority"
    878746
    879747  Reads the time history of spatial data from NetCDF file and returns a callable object.
    880748
    881   \textbf{Input variables:}
    882 
    883     \code{filename} - Name of \code{sww} or \code{tms} file (see
     749  \textbf{Arguments:}
     750
     751    \code{filename} -- Name of \code{sww} or \code{tms} file (see
    884752    Section \ref{sec:file formats} on page \pageref{sec:file formats}
    885753    for details about file formats).
    … …  
    894762
    895763       Either form will return interpolated values based on the input file
    896        using the underlying \code{interpolation_function}.
     764       using the underlying \code{interpolation\_function}.
    897765       \end{quote}
    898766
    899     \code{domain} - Associated domain object
     767    \code{domain} -- Associated domain object
    900768       If domain is specified, model time (\code{domain.starttime})
    901769       will be checked and possibly modified.
    … …  
    907775       \end{quote}
    908776
    909     \code{quantities} - the name of the quantity to be interpolated or a
     777    \code{quantities} -- the name of the quantity to be interpolated or a
    910778                 list of quantity names. The resulting function will return
    911779                 a tuple of values -- one for each quantity.
    912780
    913     \code{interpolation_points} - list of absolute UTM coordinates for points at
     781    \code{interpolation\_points} -- list of absolute UTM coordinates for points at
    914782    which values are sought
    915783
    916     \code{use_cache}: \code{True} means that caching of intermediate result of
    917                \code{Interpolation_function} is attempted
     784    \code{use\_cache}: \code{True} means that caching of intermediate result of
     785               \code{Interpolation\_function} is attempted
    918786
    919787
    920788  %  See Interpolation function for further documentation
    921 \indexedcodeheader{Interpolation_function}
     789\indexedcodeheader{Interpolation\_function}
    922790 Creates a callable object \code{f(t, id)} or \code{f(t,x,y)}
    923791    which is interpolated from time series defined at vertices of
    … …  
    927795    and $p$ the number of time steps.
    928796
    929     \textbf{Mandatory input:}
     797    \textbf{Mandatory Arguments:}
    930798
    931799        \begin{tabular}{ll}
    … …  
    939807
    940808
    941     \textbf{Optional input:}
     809    \textbf{Optional Arguments:}
    942810
    943811        \begin{tabular}{ll}
    944         \code{quantity_names}: & List of keys into the quantities dictionary\\
    945 
    946         \code{vertex_coordinates}: & $m \times 2$ array of coordinates (Float)\\
    947 
    948         \code{triangles}: & $n \times 3$ array of indices into \code{vertex_coordinates} (Int)\\
    949 
    950         \code{interpolation_points}: & $N \times 2$ array of coordinates to be interpolated to \\
     812        \code{quantity\_names}: & List of keys into the quantities dictionary\\
     813
     814        \code{vertex\_coordinates}: & $m \times 2$ array of coordinates (Float)\\
     815
     816        \code{triangles}: & $n \times 3$ array of indices into \code{vertex\_coordinates} (Int)\\
     817
     818        \code{interpolation\_points}: & $N \times 2$ array of coordinates to be interpolated to \\
    951819
    952820        \code{verbose}: & Level of reporting\\
    … …  
    964832
    965833
    966 \indexedcodeheader{set_region} ``Low priority. Will be merged into set\_quantity''
    967 
    968 \indexedcodeheader{set_quantity} ``Pretty mature''
    969 
    970 \indexedcodeheader{set_boundary} ``Pretty mature''
     834\indexedcodeheader{set\_region} ``Low priority. Will be merged into
     835set\_quantity''
     836
     837\indexedcodeheader{set\_quantity} ``Pretty mature''
     838
     839\indexedcodeheader{set\_boundary} ``Pretty mature''
    971840
    972841
    … …  
    975844\section{Diagnostics}
    976845\begin{itemize}
    977   \item \indexedcode{write_time}
    978   \item \indexedcode{write_boundary_statistics}
     846  \item \indexedcode{write\_time}
     847  \item \indexedcode{write\_boundary\_statistics}
    979848
    980849
    … …  
    984853\section{Boundary Conditions}
    985854
    986 \anuga provides a large number of predefined boundary conditions to be used with
    987 \code{set_boundary}
     855\anuga provides a large number of predefined boundary conditions to
     856be used with \code{set\_boundary}
    988857
    989858What do they do
    … …  
    1008877  \item \indexedcode{}
    1009878
    1010 
    1011879  \item \indexedcode{User defined boundary conditions.}
    1012880  How to roll your own
    1013 
    1014 
    1015 
    1016881\end{itemize}
    1017882
    … …  
    1058923\label{sec:file formats}
    1059924
    1060 
    1061 \[
    1062   \left[
    1063   \begin{array}{ccr}
    1064   2 & 4 & 4\\
    1065   1 & 1 & 1
    1066   \end{array}
    1067   \right]
    1068 \]
    1069 
     925This module takes care of reading and writing datafiles such as
     926topograhies, model output, etc
     927
     928
     929Formats used within AnuGA:
     930
     931.sww: Netcdf format for storing model output f(t,x,y) .tms: Netcdf
     932format for storing time series f(t)
     933
     934.xya: ASCII format for storing arbitrary points and associated
     935attributes .pts: NetCDF format for storing arbitrary points and
     936associated attributes
     937
     938.asc: ASCII format of regular DEMs as output from ArcView .prj:
     939Associated ArcView file giving more meta data for asc format .ers:
     940ERMapper header format of regular DEMs for ArcView
     941
     942.dem: NetCDF representation of regular DEM data
     943
     944.tsh: ASCII format for storing meshes and associated boundary and
     945region info .msh: NetCDF format for storing meshes and associated
     946boundary and region info
     947
     948.nc: Native ferret NetCDF format .geo: Houdinis ascii geometry
     949format (?)
     950
     951
     952A typical dataflow can be described as follows
     953
     954Manually created files: ASC, PRJ:     Digital elevation models
     955(gridded) TSH:          Triangular meshes (e.g. created from pmesh)
     956NC            Model outputs for use as boundary conditions (e.g from
     957MOST)
     958
     959
     960AUTOMATICALLY CREATED FILES:
     961
     962ASC, PRJ  ->  DEM  ->  PTS: Conversion of DEM's to native pts file
     963
     964NC -> SWW: Conversion of MOST bundary files to boundary sww
     965
     966PTS + TSH -> TSH with elevation: Least squares fit
     967
     968TSH -> SWW:  Conversion of TSH to sww viewable using Swollen
     969
     970TSH + Boundary SWW -> SWW: Simluation using pyvolution
     971
     972
     973%\[
     974%  \left[
     975%  \begin{array}{ccr}
     976%  2 & 4 & 4\\
     977%  1 & 1 & 1
     978%  \end{array}
     979%  \right]
     980%\]
     981
     982%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     983%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    1070984
    1071985\chapter{Basic \anuga Assumptions}
    … …  
    11241038
    11251039
    1126 
     1040%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     1041%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    11271042
    11281043\appendix
    … …  
    11361051  function call of the form
    11371052
    1138       {\scriptsize \begin{verbatim}
     1053      {\small \begin{verbatim}
    11391054      result = func(arg1,...,argn)
    11401055      \end{verbatim}}
    … …  
    11421057  can be replaced by
    11431058
    1144       {\scriptsize \begin{verbatim}
     1059      {\small \begin{verbatim}
    11451060      from caching import cache
    11461061      result = cache(func,(arg1,...,argn))
    … …  
    11751090   \textbf{USAGE:}
    11761091
    1177     {\scriptsize \begin{verbatim}
     1092    {\small \begin{verbatim}
    11781093    result = cache(func, args, kwargs, dependencies, cachedir, verbose,
    11791094                   compression, evaluate, test, return_filename)}
    … …  
    12851200\end{itemize}
    12861201
    1287 
    1288 
    1289 
    1290 
    12911202\section{coordinate_transforms}
    12921203
    … …  
    13061217\end{itemize}
    13071218
     1219%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     1220
    13081221\chapter{Glossary}
    13091222
    13101223\begin{itemize}
    1311     \item \indexedbold{\anuga} name of software (joint development between ANU and GA)
     1224    \item \indexedbold{\anuga} Name of software (joint development between ANU and GA)
    13121225
    13131226    \item \indexedbold{Conserved quantity}
    13141227
    1315     \item \indexedbold{Default order} is this really needed?
    1316 
    13171228    \item \indexedbold{Domain}
    13181229
    … …  
    13211232    \item \indexedbold{Elevation} - refers to bathymetry and topography
    13221233
    1323     \item \indexedbold{bathymetry} offshore
    1324 
    1325     \item \indexedbold{topography} onshore
     1234    \item \indexedbold{Bathymetry} offshore
     1235
     1236    \item \indexedbold{Topography} onshore
    13261237
    13271238    \item \indexedbold{Evolution} integration of the shallow water wave equations over time
    … …  
    13291240    \item \indexedbold{Forcing term}
    13301241
    1331     \item \indexedbold{IDLE} development environment shipped with Python
     1242    \item \indexedbold{IDLE} Development environment shipped with Python
    13321243
    13331244    \item \indexedbold{Manning friction coefficient}
    13341245
    1335     \item \indexedbold{Mesh}    triangulation of domain
    1336 
    1337     \item \indexedbold{meshfile}  [generic word for either .tsh or
     1246    \item \indexedbold{Mesh}    Triangulation of domain
     1247
     1248    \item \indexedbold{Meshfile}  [generic word for either .tsh or
    13381249    .msh file]
    13391250
    1340     \item \indexedbold{points file}  [generic word for either .pts or
     1251    \item \indexedbold{Points file}  [generic word for either .pts or
    13411252    .xya file]
    13421253
    … …  
    13491260    \item \indexedbold{pyvolution} does this really need to be here? it's a class/module?
    13501261
    1351     \item \indexedbold{Quantity} conserved (state, x and y momentum)
     1262    \item \indexedbold{Conserved Quantity} conserved (state, x and y momentum)
    13521263
    13531264    \item \indexedbold{Reflective boundary}
    … …  
    13671278    \item \indexedbold{ymomentum}  conserved quantity
    13681279
    1369     \item \indexedbold{resolution}   refers to the maximal area of each triangular cell in the mesh
    1370 
    1371     \item \indexedbold{polygon} A sequence of points in the plane. (Arbitrary polygons can be created
    1372     in this way )
    1373     ANUGA represents polygons as either a list of 2-tuples, where the latter are either Python tuples
    1374     or Python lists of length 2. The unit square, for example, would be represented by the polygon
    1375     [ [0,0], [1,0], [1,1], [0,1] ]. Alternatively, polygons can be represented as $N \times 2$ Numeric
    1376     arrays, where $N$ is the number of points.
     1280    \item \indexedbold{Resolution}   The maximal area of a triangular cell in a mesh
     1281
     1282    \item \indexedbold{Polygon} A sequence of points in the plane. (Arbitrary polygons can be created
     1283    in this way.)
     1284    \anuga represents a polygon in one of two ways. One way is to represent it as a
     1285    list whose members are either Python tuples
     1286    or Python lists of length 2. The unit square, for example, would be represented by the
     1287    list
     1288    [ [0,0], [1,0], [1,1], [0,1] ]. The alternative is to represent it as an
     1289    $N \times 2$ Numeric array, where $N$ is the number of points.
    13771290
    13781291    NOTE: More can be read in the module utilities/polygon.py ....
    13791292
    1380     \item \indexedbold{easting}
    1381 
    1382     \item \indexedbold{northing}
    1383 
    1384     \item \indexedbold{latitude}
    1385 
    1386     \item \indexedbold{longitude}
    1387 
    1388     \item \indexedbold{edge}
    1389 
    1390     \item \indexedbold{vertex}
    1391 
    1392     \item \indexedbold{finite volume}
    1393 
    1394     \item \indexedbold{flux}
     1293    \item \indexedbold{Easting}
     1294
     1295    \item \indexedbold{Northing}
     1296
     1297    \item \indexedbold{Latitude}
     1298
     1299    \item \indexedbold{Longitude}
     1300
     1301    \item \indexedbold{Edge}
     1302
     1303    \item \indexedbold{Vertex}
     1304
     1305    \item \indexedbold{Finite volume}
     1306
     1307    \item \indexedbold{Flux}
    13951308
    13961309    \item \indexedbold{Digital Elevation Model (DEM)}
Note: See TracChangeset for help on using the changeset viewer.