Changeset 7778


Ignore:
Timestamp:
Jun 4, 2010, 10:03:48 PM (12 years ago)
Author:
hudson
Message:

Cleaned up unit tests to use API.

Location:
trunk/anuga_core/source/anuga
Files:
1 added
30 edited
1 copied

Legend:

Unmodified
Added
Removed
  • trunk/anuga_core/source/anuga/__init__.py

    r7776 r7778  
    2222#-----------------------------------------------------
    2323
    24 from anuga.shallow_water import Domain
     24from anuga.shallow_water.shallow_water_domain import Domain
    2525
    2626from anuga.abstract_2d_finite_volumes.util import file_function
     
    2828from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
    2929
    30 from anuga.file.csv_file import csv2building_polygons
     30from anuga.file.csv_file import load_csv_as_building_polygons,  \
     31                                load_csv_as_polygons
    3132
    3233from anuga.file.sts import create_sts_boundary
     
    3536from anuga.geometry.polygon import Polygon_function
    3637
    37 from anuga.abstract_2d_finite_volumes.pmesh2domain import pmesh_to_domain_instance
     38from anuga.abstract_2d_finite_volumes.pmesh2domain import \
     39                                            pmesh_to_domain_instance
    3840
    3941
     
    219221    """
    220222
    221     from anuga.shallow_water import Domain
     223    from anuga.shallow_water.shallow_water_domain import Domain
    222224    from anuga.pmesh.mesh_interface import create_mesh_from_regions
    223225   
  • trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/test_domain.py

    r7737 r7778  
    33import unittest
    44from math import sqrt
     5
     6import anuga
    57
    68from anuga.abstract_2d_finite_volumes.generic_domain import *
    79from anuga.pmesh.mesh_interface import create_mesh_from_regions
    810from anuga.config import epsilon
    9 from anuga.shallow_water import Reflective_boundary
    10 from anuga.shallow_water import Dirichlet_boundary
    1111import numpy as num
    1212from anuga.pmesh.mesh import Segment, Vertex, Mesh
     
    484484        domain = Generic_Domain(points, vertices)
    485485
    486         domain.set_boundary( {default_boundary_tag: Dirichlet_boundary([5,2,1])} )
     486        domain.set_boundary( \
     487                {default_boundary_tag: anuga.Dirichlet_boundary([5,2,1])} )
    487488
    488489
     
    524525
    525526
    526         domain.set_boundary( {'First': Dirichlet_boundary([5,2,1]),
    527                               'Second': Transmissive_boundary(domain)} )
     527        domain.set_boundary( {'First': anuga.Dirichlet_boundary([5,2,1]),
     528                              'Second': anuga.Transmissive_boundary(domain)} )
    528529
    529530        domain.update_boundary()
     
    590591
    591592
    592         domain.set_boundary( {'First': Dirichlet_boundary([5,2,1,4,6]),
    593                               'Second': Transmissive_boundary(domain)} )
     593        domain.set_boundary( {'First': anuga.Dirichlet_boundary([5,2,1,4,6]),
     594                              'Second': anuga.Transmissive_boundary(domain)} )
    594595
    595596        try:
  • trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/test_gauge.py

    r7776 r7778  
    77from os import access, F_OK,sep, removedirs,remove,mkdir,getcwd
    88
     9import anuga
    910#from anuga.abstract_2d_finite_volumes.util import *
    1011from anuga.abstract_2d_finite_volumes.gauge import *
     
    1617
    1718from anuga.pmesh.mesh import Mesh
    18 from anuga.shallow_water import Domain, Transmissive_boundary
     19
    1920from anuga.file.sww import SWW_file
    2021from anuga.file_conversion.file_conversion import timefile2netcdf
    2122from anuga.utilities.file_utils import del_dir
    22 from csv import reader,writer
     23from csv import reader, writer
    2324import time
    2425import string
     
    4748       
    4849        # Create shallow water domain
    49         domain = Domain(mesh_file)
     50        domain = anuga.Domain(mesh_file)
    5051        os.remove(mesh_file)
    5152 
     
    6465        ######################
    6566        # Boundary conditions
    66         B = Transmissive_boundary(domain)
     67        B = anuga.Transmissive_boundary(domain)
    6768        domain.set_boundary( {'exterior': B})
    6869
     
    161162    def test_sww2csv_gauges1(self):
    162163        from anuga.pmesh.mesh import Mesh
    163         from anuga.shallow_water import Domain, Transmissive_boundary
    164164        from csv import reader,writer
    165165        import time
     
    304304    def test_sww2csv_gauge_point_off_mesh(self):
    305305        from anuga.pmesh.mesh import Mesh
    306         from anuga.shallow_water import Domain, Transmissive_boundary
    307306        from csv import reader,writer
    308307        import time
  • trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/test_util.py

    r7776 r7778  
    1818
    1919from anuga.pmesh.mesh import Mesh
    20 from anuga.shallow_water import Domain, Transmissive_boundary
     20from anuga.shallow_water.shallow_water_domain import Domain
     21from generic_boundary_conditions import Transmissive_boundary
    2122from anuga.file.sww import SWW_file
    2223from csv import reader,writer
  • trunk/anuga_core/source/anuga/culvert_flows/test_culvert_class.py

    r7711 r7778  
    1212from anuga.abstract_2d_finite_volumes.quantity import Quantity
    1313
    14 from anuga.shallow_water import Domain, Reflective_boundary,\
    15     Dirichlet_boundary,\
    16     Transmissive_boundary, Time_boundary
    17 
    18 from anuga.culvert_flows.culvert_class import Culvert_flow, Culvert_flow_rating, Culvert_flow_energy
     14import anuga
     15
     16from anuga.culvert_flows.culvert_class import Culvert_flow, \
     17                            Culvert_flow_rating, Culvert_flow_energy
    1918from anuga.culvert_flows.culvert_routines import boyd_generalised_culvert_model
    2019     
    21 from math import pi,pow,sqrt
     20from math import pi, pow, sqrt
    2221
    2322import numpy as num
     
    3837                                                   len1=length,
    3938                                                   len2=width)
    40     domain = Domain(points, vertices, boundary)   
     39    domain = anuga.Domain(points, vertices, boundary)   
    4140    domain.set_name('Test_culvert_shallow') # Output name
    4241    domain.set_default_order(2)
     
    105104    # Inflow based on Flow Depth and Approaching Momentum
    106105
    107     Br = Reflective_boundary(domain)              # Solid reflective wall
     106    Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    108107    domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    109108
     
    152151                                                       len1=length,
    153152                                                       len2=width)
    154         domain = Domain(points, vertices, boundary)   
     153        domain = anuga.Domain(points, vertices, boundary)   
    155154        domain.set_name('Test_culvert')                 # Output name
    156155        domain.set_default_order(2)
     
    214213
    215214        # Inflow based on Flow Depth and Approaching Momentum
    216         Bi = Dirichlet_boundary([0.0, 0.0, 0.0])
    217         Br = Reflective_boundary(domain)              # Solid reflective wall
    218         Bo = Dirichlet_boundary([-5, 0, 0])           # Outflow
     215        Bi = anuga.Dirichlet_boundary([0.0, 0.0, 0.0])
     216        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
     217        Bo = anuga.Dirichlet_boundary([-5, 0, 0])           # Outflow
    219218       
    220219        # Upstream and downstream conditions that will exceed the rating curve
    221220        # I.e produce delta_h outside the range [0, 10] specified in the the
    222221        # file example_rating_curve.csv
    223         Btus = Time_boundary(domain, lambda t: [100*num.sin(2*pi*(t-4)/10), 0.0, 0.0])
    224         Btds = Time_boundary(domain, lambda t: [-5*(num.cos(2*pi*(t-4)/20)), 0.0, 0.0])
     222        Btus = anuga.Time_boundary(domain, \
     223                    lambda t: [100*num.sin(2*pi*(t-4)/10), 0.0, 0.0])
     224        Btds = anuga.Time_boundary(domain, \
     225                    lambda t: [-5*(num.cos(2*pi*(t-4)/20)), 0.0, 0.0])
    225226        domain.set_boundary({'left': Btus, 'right': Btds, 'top': Br, 'bottom': Br})
    226227
     
    268269                                                       len1=length,
    269270                                                       len2=width)
    270         domain = Domain(points, vertices, boundary)   
     271        domain = anuga.Domain(points, vertices, boundary)   
    271272        domain.set_name('Test_culvert_dry') # Output name
    272273        domain.set_default_order(2)
     
    331332        # Inflow based on Flow Depth and Approaching Momentum
    332333
    333         Br = Reflective_boundary(domain)              # Solid reflective wall
     334        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    334335        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    335336
     
    388389                                                       len1=length,
    389390                                                       len2=width)
    390         domain = Domain(points, vertices, boundary)   
     391        domain = anuga.Domain(points, vertices, boundary)   
    391392        domain.set_name('Test_culvert_shallow') # Output name
    392393        domain.set_default_order(2)
     
    463464        # Inflow based on Flow Depth and Approaching Momentum
    464465
    465         Br = Reflective_boundary(domain)              # Solid reflective wall
     466        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    466467        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    467468
     
    525526                                                       len1=length,
    526527                                                       len2=width)
    527         domain = Domain(points, vertices, boundary)   
     528        domain = anuga.Domain(points, vertices, boundary)   
    528529        domain.set_name('Test_culvert_dry') # Output name
    529530        domain.set_default_order(2)
     
    594595        # Inflow based on Flow Depth and Approaching Momentum
    595596
    596         Br = Reflective_boundary(domain)              # Solid reflective wall
     597        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    597598        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    598599
     
    635636                                                       len1=length,
    636637                                                       len2=width)
    637         domain = Domain(points, vertices, boundary)   
     638        domain = anuga.Domain(points, vertices, boundary)   
    638639       
    639640        domain.set_name('test_culvert')                 # Output name
     
    711712                                                       len1=length,
    712713                                                       len2=width)
    713         domain = Domain(points, vertices, boundary)   
     714        domain = anuga.Domain(points, vertices, boundary)   
    714715        domain.set_name('Test_culvert_shallow') # Output name
    715716        domain.set_default_order(2)
     
    782783
    783784
    784         Br = Reflective_boundary(domain)              # Solid reflective wall
     785        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    785786        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    786787
  • trunk/anuga_core/source/anuga/culvert_flows/test_new_culvert_class.py

    r7727 r7778  
    66import sys
    77
     8import anuga
     9
    810from anuga.utilities.system_tools import get_pathname_from_package
    911from anuga.geometry.polygon import Polygon_function
     
    1113from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
    1214from anuga.abstract_2d_finite_volumes.quantity import Quantity
    13 
    14 from anuga.shallow_water import Domain, Reflective_boundary,\
    15     Dirichlet_boundary,\
    16     Transmissive_boundary, Time_boundary
    1715
    1816from anuga.culvert_flows.new_culvert_class import Culvert_flow
     
    3836                                                   len1=length,
    3937                                                   len2=width)
    40     domain = Domain(points, vertices, boundary)   
     38    domain = anuga.Domain(points, vertices, boundary)   
    4139    domain.set_name('Test_culvert_shallow') # Output name
    4240    domain.set_default_order(2)
     
    105103    # Inflow based on Flow Depth and Approaching Momentum
    106104
    107     Br = Reflective_boundary(domain)              # Solid reflective wall
     105    Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    108106    domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    109107
     
    152150                                                       len1=length,
    153151                                                       len2=width)
    154         domain = Domain(points, vertices, boundary)   
     152        domain = anuga.Domain(points, vertices, boundary)   
    155153        domain.set_name('Test_culvert')                 # Output name
    156154        domain.set_default_order(2)
     
    214212
    215213        # Inflow based on Flow Depth and Approaching Momentum
    216         Bi = Dirichlet_boundary([0.0, 0.0, 0.0])
    217         Br = Reflective_boundary(domain)              # Solid reflective wall
    218         Bo = Dirichlet_boundary([-5, 0, 0])           # Outflow
     214        Bi = anuga.Dirichlet_boundary([0.0, 0.0, 0.0])
     215        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
     216        Bo = anuga.Dirichlet_boundary([-5, 0, 0])           # Outflow
    219217       
    220218        # Upstream and downstream conditions that will exceed the rating curve
    221219        # I.e produce delta_h outside the range [0, 10] specified in the the
    222220        # file example_rating_curve.csv
    223         Btus = Time_boundary(domain, lambda t: [100*num.sin(2*pi*(t-4)/10), 0.0, 0.0])
    224         Btds = Time_boundary(domain, lambda t: [-5*(num.cos(2*pi*(t-4)/20)), 0.0, 0.0])
     221        Btus = anuga.Time_boundary(domain, lambda t: [100*num.sin(2*pi*(t-4)/10), 0.0, 0.0])
     222        Btds = anuga.Time_boundary(domain, lambda t: [-5*(num.cos(2*pi*(t-4)/20)), 0.0, 0.0])
    225223        domain.set_boundary({'left': Btus, 'right': Btds, 'top': Br, 'bottom': Br})
    226224
     
    268266                                                       len1=length,
    269267                                                       len2=width)
    270         domain = Domain(points, vertices, boundary)   
     268        domain = anuga.Domain(points, vertices, boundary)   
    271269        domain.set_name('Test_culvert_dry') # Output name
    272270        domain.set_default_order(2)
     
    331329        # Inflow based on Flow Depth and Approaching Momentum
    332330
    333         Br = Reflective_boundary(domain)              # Solid reflective wall
     331        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    334332        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    335333
     
    388386                                                       len1=length,
    389387                                                       len2=width)
    390         domain = Domain(points, vertices, boundary)   
     388        domain = anuga.Domain(points, vertices, boundary)   
    391389        domain.set_name('Test_culvert_shallow') # Output name
    392390        domain.set_default_order(2)
     
    463461        # Inflow based on Flow Depth and Approaching Momentum
    464462
    465         Br = Reflective_boundary(domain)              # Solid reflective wall
     463        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    466464        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    467465
     
    525523                                                       len1=length,
    526524                                                       len2=width)
    527         domain = Domain(points, vertices, boundary)   
     525        domain = anuga.Domain(points, vertices, boundary)   
    528526        domain.set_name('Test_culvert_dry') # Output name
    529527        domain.set_default_order(2)
     
    594592        # Inflow based on Flow Depth and Approaching Momentum
    595593
    596         Br = Reflective_boundary(domain)              # Solid reflective wall
     594        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    597595        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    598596
     
    635633                                                       len1=length,
    636634                                                       len2=width)
    637         domain = Domain(points, vertices, boundary)   
     635        domain = anuga.Domain(points, vertices, boundary)   
    638636       
    639637        domain.set_name('test_culvert')                 # Output name
     
    711709                                                       len1=length,
    712710                                                       len2=width)
    713         domain = Domain(points, vertices, boundary)   
     711        domain = anuga.Domain(points, vertices, boundary)   
    714712        domain.set_name('Test_culvert_shallow') # Output name
    715713        domain.set_default_order(2)
     
    782780
    783781
    784         Br = Reflective_boundary(domain)              # Solid reflective wall
     782        Br = anuga.Reflective_boundary(domain)              # Solid reflective wall
    785783        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    786784
  • trunk/anuga_core/source/anuga/damage_modelling/test_inundation_damage.py

    r7776 r7778  
    1313from anuga.pmesh.mesh import Mesh
    1414from anuga.coordinate_transforms.geo_reference import Geo_reference
    15 from anuga.shallow_water import Domain, Transmissive_boundary
    1615from anuga.utilities.numerical_tools import mean
    1716from anuga.file.sww import SWW_file
  • trunk/anuga_core/source/anuga/file/csv_file.py

    r7776 r7778  
    242242
    243243
    244 def csv2building_polygons(file_name,
     244def load_csv_as_building_polygons(file_name,
    245245                          floor_height=3,
    246246                          clipping_polygons=None):
     
    287287# @brief Convert CSV file into a dictionary of polygons and associated values.
    288288# @param filename The path to the file to read, value_name name for the 4th column
    289 def csv2polygons(file_name,
     289def load_csv_as_polygons(file_name,
    290290                 value_name='value',
    291291                 clipping_polygons=None):
  • trunk/anuga_core/source/anuga/file/pts.py

    r7776 r7778  
    1616
    1717    from anuga.pmesh.mesh import Mesh, importMeshFromFile
    18     from anuga.shallow_water import Domain
     18    from anuga.shallow_water.shallow_water_domain import Domain
    1919    from anuga.pmesh.mesh_interface import create_mesh_from_regions
    2020
  • trunk/anuga_core/source/anuga/file/sts.py

    r7766 r7778  
    11import numpy as num
     2import log
     3from Scientific.IO.NetCDF import NetCDFFile
     4
    25from anuga.config import max_float
    36from anuga.config import netcdf_float, netcdf_float32, netcdf_int
     
    58from anuga.coordinate_transforms.geo_reference import Geo_reference, \
    69     ensure_geo_reference
    7 
     10from anuga.config import netcdf_mode_r, netcdf_mode_w, netcdf_mode_a
    811
    912class Write_sts:
     
    268271# @param stsname Stem of path to the STS file to read.
    269272# @return A list of boundary points.
    270 def create_sts_boundary(stsname):
     273def create_sts_boundary(sts_filename):
    271274    """Create a list of points defining a boundary from an STS file.
    272275
     
    280283    """
    281284
     285    if sts_filename.endswith('.sts'):
     286        stsname_postfixed = sts_filename
     287    else:
     288        stsname_postfixed = sts_filename + '.sts'
     289
    282290    try:
    283         fid = NetCDFFile(stsname + '.sts', netcdf_mode_r)
    284     except:
    285         msg = 'Cannot open %s' % stsname + '.sts'
     291        fid = NetCDFFile(stsname_postfixed, netcdf_mode_r)
     292    except IOError:
     293        msg = 'Cannot open %s' % stsname_postfixed
    286294        raise IOError(msg)
    287295
  • trunk/anuga_core/source/anuga/file_conversion/__init__.py

    r7776 r7778  
    11""" Modules for performing conversions between file types, or for
    2     resampling a given file.
     2    resampling a given file.
     3   
     4    In general, the naming convention follows this rule:
     5        <file_ext_in>2<file_ext_out>('filename_in', 'filename_out')
     6       
     7    for example:
     8        sww2dem('northbeach.sww', 'outfile.dem')
    39"""
    410
  • trunk/anuga_core/source/anuga/file_conversion/dem2dem.py

    r7776 r7778  
    22"""
    33
     4import numpy as num
     5
     6from anuga.config import netcdf_mode_r, netcdf_mode_w, netcdf_float
    47
    58def dem2dem(basename_in, stencil, cellsize_new, basename_out=None,
  • trunk/anuga_core/source/anuga/file_conversion/sww2dem.py

    r7776 r7778  
     1"""
     2    Module to convert SWW to DEM files.
     3"""
     4
    15# external modules
    26import numpy as num
     
    417421
    418422
    419 
    420 ##
    421 # @brief
    422 # @param basename_in
    423 # @param extra_name_out
    424 # @param quantities
    425 # @param timestep
    426 # @param reduction
    427 # @param cellsize
    428 # @param number_of_decimal_places
    429 # @param NODATA_value
    430 # @param easting_min
    431 # @param easting_max
    432 # @param northing_min
    433 # @param northing_max
    434 # @param verbose
    435 # @param origin
    436 # @param datum
    437 # @param format
    438 # @return
    439423def sww2dem_batch(basename_in, extra_name_out=None,
    440424                quantities=None, # defaults to elevation
  • trunk/anuga_core/source/anuga/file_conversion/test_urs2sts.py

    r7770 r7778  
    1313from anuga.file.csv_file import load_csv_as_dict, load_csv_as_array
    1414
     15from anuga.shallow_water.shallow_water_domain import Domain
     16
    1517# boundary functions
    1618from anuga.shallow_water.boundaries import Reflective_boundary, \
     
    2022     import Transmissive_boundary, Dirichlet_boundary, \
    2123            Time_boundary, File_boundary, AWI_boundary
     24
     25from anuga.pmesh.mesh_interface import create_mesh_from_regions
    2226
    2327from urs2sts import urs2sts
     
    13161320        """
    13171321       
    1318         # FIXME (Ole): These tests should really move to
    1319         # test_generic_boundaries.py
    1320        
    1321         from anuga.shallow_water import Domain
    1322         from anuga.shallow_water import Reflective_boundary
    1323         from anuga.shallow_water import Dirichlet_boundary
    1324         from anuga.shallow_water import File_boundary
    1325         from anuga.pmesh.mesh_interface import create_mesh_from_regions
    1326 
    13271322        bounding_polygon=[[6.0,97.0],[6.01,97.0],[6.02,97.0],[6.02,97.02],[6.00,97.02]]
    13281323        tide = 0.37
     
    14581453        import warnings
    14591454        warnings.simplefilter('ignore')
    1460        
    1461         from anuga.shallow_water import Domain
    1462         from anuga.shallow_water import Reflective_boundary
    1463         from anuga.shallow_water import Dirichlet_boundary
    1464         from anuga.shallow_water import File_boundary
    1465         from anuga.pmesh.mesh_interface import create_mesh_from_regions
    14661455
    14671456        bounding_polygon=[[6.0,97.0],[6.01,97.0],[6.02,97.0],
     
    16051594        import warnings
    16061595        warnings.simplefilter('ignore')
    1607        
    1608         from anuga.shallow_water import Domain
    1609         from anuga.shallow_water import Reflective_boundary
    1610         from anuga.shallow_water import Dirichlet_boundary
    1611         from anuga.shallow_water import File_boundary
    1612         from anuga.pmesh.mesh_interface import create_mesh_from_regions
    16131596
    16141597        bounding_polygon=[[6.0,97.0],[6.01,97.0],[6.02,97.0],
     
    17231706         the correct order).
    17241707         """
    1725          
    1726         from anuga.shallow_water import Domain
    1727         from anuga.shallow_water import Reflective_boundary
    1728         from anuga.shallow_water import Dirichlet_boundary
    1729         from anuga.shallow_water import File_boundary
    1730         from anuga.pmesh.mesh_interface import create_mesh_from_regions
    17311708
    17321709        bounding_polygon=[[6.01,97.0],[6.02,97.0],[6.02,97.02],[6.00,97.02],[6.0,97.0]]
     
    18311808        Read correct points from ordering file and apply sts to boundary
    18321809        """
    1833         from anuga.shallow_water import Domain
    1834         from anuga.shallow_water import Reflective_boundary
    1835         from anuga.shallow_water import Dirichlet_boundary
    1836         from anuga.shallow_water import File_boundary
    1837         from anuga.pmesh.mesh_interface import create_mesh_from_regions
    18381810
    18391811        lat_long_points=[[6.01,97.0],[6.02,97.0],[6.05,96.9],[6.0,97.0]]
     
    18591831                                           ua=ua,
    18601832                                           va=va)
     1833        # base name will not exist, but 3 other files are created
    18611834
    18621835        # Write order file
     
    18861859        self.delete_mux(files)
    18871860
    1888         assert(os.access(base_name, os.R_OK))
     1861        assert(os.access(sts_file+'.sts', os.F_OK))
    18891862
    18901863        boundary_polygon = create_sts_boundary(base_name)
     
    19791952        #assert allclose(domain_drchlt.quantities['stage'].vertex_values[6], 2)       
    19801953        #assert allclose(domain_fbound.quantities['stage'].vertex_values[6], 2)
    1981        
    1982        
    19831954
    19841955        try:
    19851956            os.remove(sts_file+'.sts')
    1986         except:
     1957        except IOError:
    19871958            # Windoze can't remove this file for some reason
    19881959            pass
  • trunk/anuga_core/source/anuga/file_conversion/urs2sts.py

    r7765 r7778  
    11import numpy as num
     2import log
    23
    34# ANUGA modules
     
    315316
    316317    outfile.close()
    317 
     318   
     319    if verbose:
     320        log.critical('Wrote sts file ' + stsname)   
     321   
     322
  • trunk/anuga_core/source/anuga/fit_interpolate/benchmark_least_squares.py

    r7751 r7778  
    2828from anuga.pmesh.mesh import Mesh
    2929from anuga.geospatial_data.geospatial_data import Geospatial_data
    30 from anuga.shallow_water import Domain
     30from anuga.shallow_water.shallow_water_domain import Domain
    3131from anuga.fit_interpolate.fit import Fit, fit_to_mesh
    3232from anuga.fit_interpolate.interpolate import benchmark_interpolate
  • trunk/anuga_core/source/anuga/fit_interpolate/interpolate.py

    r7717 r7778  
    3333from anuga.geospatial_data.geospatial_data import Geospatial_data
    3434from anuga.geospatial_data.geospatial_data import ensure_absolute
    35 from anuga.geometry.mesh_quadtree import MeshQuadtree
     35from anuga.pmesh.mesh_quadtree import MeshQuadtree
    3636from anuga.fit_interpolate.general_fit_interpolate import FitInterpolate
    3737from anuga.abstract_2d_finite_volumes.file_function import file_function
  • trunk/anuga_core/source/anuga/fit_interpolate/test_fit.py

    r7276 r7778  
    1919from anuga.utilities.numerical_tools import ensure_numeric
    2020from anuga.geospatial_data.geospatial_data import Geospatial_data
    21 from anuga.shallow_water import Domain
     21from anuga.shallow_water.shallow_water_domain import Domain
    2222
    2323import numpy as num
  • trunk/anuga_core/source/anuga/fit_interpolate/test_interpolate.py

    r7776 r7778  
    1818
    1919
    20 
    2120# ANUGA code imports
     21import anuga
    2222from interpolate import *
    2323from anuga.coordinate_transforms.geo_reference import Geo_reference
    24 from anuga.shallow_water import Domain, Transmissive_boundary
    2524from anuga.utilities.numerical_tools import mean, NAN
    2625from anuga.file.sww import SWW_file
  • trunk/anuga_core/source/anuga/geometry/polygon.py

    r7711 r7778  
    77from math import sqrt
    88from anuga.utilities.numerical_tools import ensure_numeric
    9 from anuga.geospatial_data.geospatial_data import ensure_absolute, Geospatial_data
     9from anuga.geospatial_data.geospatial_data import ensure_absolute, \
     10                                                    Geospatial_data
    1011from anuga.config import netcdf_float
    1112import anuga.utilities.log as log
     
    6970                         % (str(p1), str(p2), str(p3), str(p4)))
    7071
    71 #                     0s1    0e1    1s0    1e0   # line 0 starts on 1, 0 ends 1, 1 starts 0, 1 ends 0
    72 collinear_result = { (False, False, False, False): lines_dont_coincide,
    73                      (False, False, False, True ): lines_error,
    74                      (False, False, True,  False): lines_error,
    75                      (False, False, True,  True ): lines_1_fully_included_in_0,
    76                      (False, True,  False, False): lines_error,
    77                      (False, True,  False, True ): lines_overlap_opposite_direction2,
    78                      (False, True,  True,  False): lines_overlap_same_direction2,
    79                      (False, True,  True,  True ): lines_1_fully_included_in_0,
    80                      (True,  False, False, False): lines_error,
    81                      (True,  False, False, True ): lines_overlap_same_direction,
    82                      (True,  False, True,  False): lines_overlap_opposite_direction,
    83                      (True,  False, True,  True ): lines_1_fully_included_in_0,
    84                      (True,  True,  False, False): lines_0_fully_included_in_1,
    85                      (True,  True,  False, True ): lines_0_fully_included_in_1,
    86                      (True,  True,  True,  False): lines_0_fully_included_in_1,
    87                      (True,  True,  True,  True ): lines_0_fully_included_in_1
    88                    }
     72collinear_result = {
     73# line 0 starts on 1, 0 ends 1, 1 starts 0, 1 ends 0
     74#       0s1    0e1    1s0    1e0   
     75       (False, False, False, False): lines_dont_coincide,
     76       (False, False, False, True ): lines_error,
     77       (False, False, True,  False): lines_error,
     78       (False, False, True,  True ): lines_1_fully_included_in_0,
     79       (False, True,  False, False): lines_error,
     80       (False, True,  False, True ): lines_overlap_opposite_direction2,
     81       (False, True,  True,  False): lines_overlap_same_direction2,
     82       (False, True,  True,  True ): lines_1_fully_included_in_0,
     83       (True,  False, False, False): lines_error,
     84       (True,  False, False, True ): lines_overlap_same_direction,
     85       (True,  False, True,  False): lines_overlap_opposite_direction,
     86       (True,  False, True,  True ): lines_1_fully_included_in_0,
     87       (True,  True,  False, False): lines_0_fully_included_in_1,
     88       (True,  True,  False, True ): lines_0_fully_included_in_1,
     89       (True,  True,  True,  False): lines_0_fully_included_in_1,
     90       (True,  True,  True,  True ): lines_0_fully_included_in_1
     91   }
    8992
    9093##
     
    275278   
    276279    # Quickly reject points that are clearly outside
    277     if point[0] < min(triangle[:,0]): return False
    278     if point[0] > max(triangle[:,0]): return False   
    279     if point[1] < min(triangle[:,1]): return False
    280     if point[1] > max(triangle[:,1]): return False       
     280    if point[0] < min(triangle[:, 0]): return False
     281    if point[0] > max(triangle[:, 0]): return False   
     282    if point[1] < min(triangle[:, 1]): return False
     283    if point[1] > max(triangle[:, 1]): return False       
    281284
    282285
     
    367370                (type, point) = intersection(leftmost, l_x[cmp])
    368371                comparisons += 1
    369                 if type != 0 and type != 4 or (type == 2 and list(point[0]) != list(point[1])):
     372                if type != 0 and type != 4 or (type == 2 and list(point[0]) !=\
     373                                                list(point[1])):
    370374                    if verbose:
    371                         print 'Self-intersecting polygon found, type ', type, ' point', point,
     375                        print 'Self-intersecting polygon found, type ', type
     376                        print 'point', point,
    372377                        print 'vertices: ', leftmost, ' - ', l_x[cmp]               
    373378                    return True           
     
    586591    if count == len(indices):
    587592        # No points are outside
    588         return indices[:count],[]
     593        return indices[:count], []
    589594    else:
    590595        return  indices[:count], indices[count:][::-1]  #return reversed
    591596
    592 ##
    593 # @brief Sort a list of points into contiguous points inside+outside a polygon.
    594 # @param points A set of points (tuple, list or array).
    595 # @param polygon A set of points defining a polygon (tuple, list or array).
    596 # @param closed True if points on boundary are considered 'inside' polygon.
    597 # @param verbose True if this function is to be verbose.
    598 # @return (indices, count) where indices are point indices and count is the
    599 #          delimiter index between point inside (on left) and others.
     597
     598
    600599def separate_points_by_polygon(points, polygon,
    601600                               closed=True,
     
    644643    if check_input:
    645644        #Input checks
    646         assert isinstance(closed, bool), 'Keyword argument "closed" must be boolean'
    647         assert isinstance(verbose, bool), 'Keyword argument "verbose" must be boolean'
     645        assert isinstance(closed, bool), \
     646                    'Keyword argument "closed" must be boolean'
     647        assert isinstance(verbose, bool), \
     648                    'Keyword argument "verbose" must be boolean'
    648649
    649650        try:
     
    737738    return abs(poly_area/2)
    738739
    739 ##
    740 # @brief Plot a set of polygons.
    741 # @param polygons_points List of polygons to plot.
    742 # @param style List of styles for each polygon.
    743 # @param figname Name to save figure to.
    744 # @param label Title for the plot.
    745 # @param verbose True if this function is to be verbose.
    746 # @return A list of min/max x and y values [minx, maxx, miny, maxy].
    747 # @note A style value is 'line' for polygons, 'outside' for points outside.
     740
    748741def plot_polygons(polygons_points,
    749742                  style=None,
     
    818811            if s == 'outside': colour.append('r.')
    819812            if s == 'point': colour.append('g.')
    820             if s <> 'line':
    821                 if s <> 'outside':
    822                     if s <> 'point':
     813            if s != 'line':
     814                if s != 'outside':
     815                    if s != 'point':
    823816                        colour.append(s)
    824817
     
    836829        title(label)
    837830
    838     #raw_input('wait 1')
    839     #FIXME(Ole): This makes for some strange scalings sometimes.
    840     #if minx <> 0:
    841     #    axis([minx*0.9,maxx*1.1,miny*0.9,maxy*1.1])
    842     #else:
    843     #    if miny == 0:
    844     #        axis([-maxx*.01,maxx*1.1,-maxy*0.01,maxy*1.1])
    845     #    else:
    846     #        axis([-maxx*.01,maxx*1.1,miny*0.9,maxy*1.1])
    847 
    848831    if figname is not None:
    849832        savefig(figname)
     
    856839    return vec
    857840
    858 ##
    859 # @brief
    860 # @param polygon A set of points defining a polygon.
    861 # @param verbose True if this function is to be verbose.
    862 # @return A tuple (x, y) of X and Y coordinates of the polygon.
    863 # @note We duplicate the first point so can have closed polygon in plot.
     841
    864842def poly_xy(polygon, verbose=False):
    865843    """ this is used within plot_polygons so need to duplicate
    866844        the first point so can have closed polygon in plot
     845        # @param polygon A set of points defining a polygon.
     846        # @param verbose True if this function is to be verbose.
     847        # @return A tuple (x, y) of X and Y coordinates of the polygon.
     848        # @note We duplicate the first point so can have closed polygon in plot.
    867849    """
    868850
     
    884866
    885867
    886 ##
    887 # @brief Define a class that defines a callable object for a polygon.
    888 # @note Object created is function: f: x,y -> z
    889 #       where x, y and z are vectors and z depends on whether x,y belongs
    890 #       to specified polygons.
    891868class Polygon_function:
    892869    """Create callable object f: x,y -> z, where a,y,z are vectors and
     
    923900    """
    924901
    925     ##
    926     # @brief Create instance of a polygon function.
    927     # @param regions A list of (x,y) tuples defining a polygon.
    928     # @param default Value or function returning value for points outside poly.
    929     # @param geo_reference ??
    930902    def __init__(self, regions, default=0.0, geo_reference=None):
     903        """
     904        # @brief Create instance of a polygon function.
     905        # @param regions A list of (x,y) tuples defining a polygon.
     906        # @param default Value or function returning value for points outside poly.
     907        # @param geo_reference ??
     908        """
     909
    931910        try:
    932911            len(regions)
     
    967946            self.regions.append((P, value))
    968947
    969     ##
    970     # @brief Implement the 'callable' property of Polygon_function.
    971     # @param x List of x coordinates of points ot interest.
    972     # @param y List of y coordinates of points ot interest.
    973948    def __call__(self, x, y):
     949        """
     950        # @brief Implement the 'callable' property of Polygon_function.
     951        # @param x List of x coordinates of points ot interest.
     952        # @param y List of y coordinates of points ot interest.
     953        """
    974954        x = num.array(x, num.float)
    975955        y = num.array(y, num.float)
     
    1014994################################################################################
    1015995
    1016 ##
    1017 # @brief Read polygon data from a file.
     996def read_polygon(filename, delimiter=','):
     997    """Read points assumed to form a polygon.
     998
    1018999# @param filename Path to file containing polygon data.
    10191000# @param delimiter Delimiter to split polygon data with.
    10201001# @return A list of point data from the polygon file.
    1021 def read_polygon(filename, delimiter=','):
    1022     """Read points assumed to form a polygon.
     1002
    10231003
    10241004    There must be exactly two numbers in each line separated by the delimiter.
     
    10361016    # check this is a valid polygon.
    10371017    if is_complex(polygon, verbose=True):   
    1038         msg = 'ERROR: Self-intersecting polygon detected in file ' + filename +'. '
    1039         msg += 'A complex polygon will not necessarily break the algorithms within ANUGA, '
    1040         msg += 'but it usually signifies pathological data. Please fix this file.'
     1018        msg = 'ERROR: Self-intersecting polygon detected in file '
     1019        msg += filename +'. A complex polygon will not '
     1020        msg += 'necessarily break the algorithms within ANUGA, but it'
     1021        msg += 'usually signifies pathological data. Please fix this file.'
    10411022        raise Exception, msg
    10421023   
     
    10691050    pass
    10701051
    1071 ##
    1072 # @brief Populate given polygon with uniformly distributed points.
    1073 # @param polygon Polygon to uniformly fill.
    1074 # @param number_of_points Number of points required in polygon.
    1075 # @param seed Seed for random number generator.
    1076 # @param exclude List of polygons inside main where points should be excluded.
    1077 # @return List of random points inside input polygon.
    1078 # @note Delimiter is assumed to be a comma.
     1052
    10791053def populate_polygon(polygon, number_of_points, seed=None, exclude=None):
    10801054    """Populate given polygon with uniformly distributed points.
     
    11311105    return points
    11321106
    1133 ##
    1134 # @brief Get a point inside a polygon that is close to an edge.
    1135 # @param polygon List  of vertices of polygon.
    1136 # @param delta Maximum distance from an edge is delta * sqrt(2).
    1137 # @return A point that is inside polgon and close to the polygon edge.
     1107
    11381108def point_in_polygon(polygon, delta=1e-8):
    11391109    """Return a point inside a given polygon which will be close to the
     
    11861156    return point
    11871157
    1188 ##
    1189 # @brief Calculate approximate number of triangles inside a bounding polygon.
    1190 # @param interior_regions
    1191 # @param bounding_poly
    1192 # @param remainder_res
    1193 # @return The number of triangles.
     1158
    11941159def number_mesh_triangles(interior_regions, bounding_poly, remainder_res):
    11951160    """Calculate the approximate number of triangles inside the
     
    12471212    return int(total_number_of_triangles)
    12481213
     1214
     1215def decimate_polygon(polygon, factor=10):
     1216    """Reduce number of points in polygon by the specified
     1217    factor (default=10, hence the name of the function) such that
     1218    the extrema in both axes are preserved.
     1219
    12491220##
    12501221# @brief Reduce number of points in polygon by the specified factor.
    12511222# @param polygon The polygon to reduce.
    12521223# @param factor The factor to reduce polygon points by (default 10).
    1253 # @return The reduced polygon points list.
    12541224# @note The extrema of both axes are preserved.
    1255 def decimate_polygon(polygon, factor=10):
    1256     """Reduce number of points in polygon by the specified
    1257     factor (default=10, hence the name of the function) such that
    1258     the extrema in both axes are preserved.
    12591225
    12601226    Return reduced polygon
     
    12851251    return reduced_polygon
    12861252
    1287 ##
    1288 # @brief Interpolate linearly from polyline nodes to midpoints of triangles.
    1289 # @param data The data on the polyline nodes.
    1290 # @param polyline_nodes ??
    1291 # @param gauge_neighbour_id ??  FIXME(Ole): I want to get rid of this
    1292 # @param point_coordinates ??
    1293 # @param verbose True if this function is to be verbose.
     1253
    12941254def interpolate_polyline(data,
    12951255                         polyline_nodes,
     
    13701330               
    13711331    return points, vertices
    1372 ##
    1373 # @brief
    1374 # @param data
    1375 # @param polyline_nodes
    1376 # @param gauge_neighbour_id
    1377 # @param interpolation_points
    1378 # @param interpolated_values
    1379 # @param rtol
    1380 # @param atol
    1381 # @return
    1382 # @note OBSOLETED BY C-EXTENSION
    1383 def _interpolate_polyline(data,
    1384                           polyline_nodes,
    1385                           gauge_neighbour_id,
    1386                           interpolation_points,
    1387                           interpolated_values,
    1388                           rtol=1.0e-6,
    1389                           atol=1.0e-8):
    1390     """Auxiliary function used by interpolate_polyline
    1391 
    1392     NOTE: OBSOLETED BY C-EXTENSION
    1393     """
    1394 
    1395     number_of_nodes = len(polyline_nodes)
    1396     number_of_points = len(interpolation_points)
    1397 
    1398     for j in range(number_of_nodes):
    1399         neighbour_id = gauge_neighbour_id[j]
    1400 
    1401         # FIXME(Ole): I am convinced that gauge_neighbour_id can be discarded,
    1402         #             but need to check with John J.
    1403         # Keep it for now (17 Jan 2009)
    1404         # When gone, we can simply interpolate between neighbouring nodes,
    1405         # i.e. neighbour_id = j+1.
    1406         # and the test below becomes something like: if j < number_of_nodes...
    1407 
    1408         if neighbour_id >= 0:
    1409             x0, y0 = polyline_nodes[j,:]
    1410             x1, y1 = polyline_nodes[neighbour_id,:]
    1411 
    1412             segment_len = sqrt((x1-x0)**2 + (y1-y0)**2)
    1413             segment_delta = data[neighbour_id] - data[j]
    1414             slope = segment_delta/segment_len
    1415 
    1416             for i in range(number_of_points):
    1417                 x, y = interpolation_points[i,:]
    1418                 if point_on_line([x, y], [[x0, y0], [x1, y1]],
    1419                                  rtol=rtol, atol=atol):
    1420                     dist = sqrt((x-x0)**2 + (y-y0)**2)
    1421                     interpolated_values[i] = slope*dist + data[j]
    14221332
    14231333
  • trunk/anuga_core/source/anuga/pmesh/mesh_quadtree.py

    r7751 r7778  
    1111from anuga.utilities.numerical_tools import ensure_numeric
    1212from anuga.config import max_float
    13 from quad import Cell
    14 from aabb import AABB
     13
     14from anuga.geometry.quad import Cell
     15from anuga.geometry.aabb import AABB
    1516
    1617from anuga.utilities import compile
  • trunk/anuga_core/source/anuga/shallow_water/__init__.py

    r7736 r7778  
    1 """Make directory available as a Python package
     1""" Shallow Water Wave Module.
     2
     3    This module performs the shallow water wave simulation, and is a
     4    specialisation of the ANUGA domain. It contains supporting modules and
     5    extra boundary classes to do with this specific domain.
    26"""
    37
     
    610sys.path += __path__
    711
    8 # Make selected classes available directly
    9 from boundaries import Reflective_boundary,\
    10      Transmissive_momentum_set_stage_boundary,\
    11      Dirichlet_discharge_boundary,\
    12      Field_boundary,\
    13      Transmissive_stage_zero_momentum_boundary,\
    14      Transmissive_n_momentum_zero_t_momentum_set_stage_boundary
    15 
    16 from anuga.abstract_2d_finite_volumes.generic_boundary_conditions\
    17      import Transmissive_boundary, Dirichlet_boundary, \
    18             Time_boundary, File_boundary, AWI_boundary
    19 
    20 from shallow_water_domain import Domain
    21 
    22 
    23 #from shallow_water_balanced_domain import Swb_domain
    24 
    25 
    2612 
  • trunk/anuga_core/source/anuga/shallow_water/shallow_water_domain.py

    r7776 r7778  
    511511
    512512
    513 
    514     ##
    515     # @brief Get the total flow through an arbitrary poly line.
    516     # @param polyline Representation of desired cross section.
    517     # @param verbose True if this method is to be verbose.
    518     # @note 'polyline' may contain multiple sections allowing complex shapes.
    519     # @note Assume absolute UTM coordinates.
    520513    def get_flow_through_cross_section(self, polyline, verbose=False):
    521514        """Get the total flow through an arbitrary poly line.
    522515
    523516        This is a run-time equivalent of the function with same name
    524         in data_manager.py
     517        in sww_interrogate.py
    525518
    526519        Input:
     
    540533
    541534
    542 
    543 
    544     ##
    545     # @brief
    546     # @param polyline Representation of desired cross section.
    547     # @param kind Select energy type to compute ('specific' or 'total').
    548     # @param verbose True if this method is to be verbose.
    549     # @note 'polyline' may contain multiple sections allowing complex shapes.
    550     # @note Assume absolute UTM coordinates.
    551535    def get_energy_through_cross_section(self, polyline,
    552536                                         kind='total',
  • trunk/anuga_core/source/anuga/shallow_water/sww_interrogate.py

    r7771 r7778  
    99from anuga.coordinate_transforms.geo_reference import Geo_reference
    1010from anuga.abstract_2d_finite_volumes.util import file_function
    11 from anuga.config import netcdf_mode_r, netcdf_mode_w, netcdf_mode_a
    1211from anuga.geometry.polygon import is_inside_polygon
    1312from anuga.file.sww import get_mesh_and_quantities_from_file
    1413from anuga.abstract_2d_finite_volumes.neighbour_mesh import segment_midpoints
    15 from anuga.config import g
    1614
    1715##
     
    374372    maximal_runup_location = None
    375373
    376     for file, swwfile in enumerate (iterate_over):
     374    for _, swwfile in enumerate (iterate_over):
    377375        # Read sww file
    378376        filename = os.path.join(dir, swwfile+'.sww')
     
    387385        try:
    388386            geo_reference = Geo_reference(NetCDFObject=fid)
    389         except AttributeError, e:
     387        except AttributeError:
    390388            geo_reference = Geo_reference() # Default georef object
    391389
    392390        xllcorner = geo_reference.get_xllcorner()
    393391        yllcorner = geo_reference.get_yllcorner()
    394         zone = geo_reference.get_zone()
    395392
    396393        # Get extent
     
    417414            # FIXME (Ole): Make a generic polygon input check in polygon.py
    418415            # and call it here
    419             points = num.ascontiguousarray(num.concatenate((x[:,num.newaxis],
    420                                                             y[:,num.newaxis]),
     416            points = num.ascontiguousarray(num.concatenate((x[:, num.newaxis],
     417                                                            y[:, num.newaxis]),
    421418                                                            axis=1))
    422419            point_indices = inside_polygon(points, polygon)
     
    497494                wet_x = num.take(x, wet_nodes, axis=0)
    498495                wet_y = num.take(y, wet_nodes, axis=0)
    499                 maximal_runup_location = [wet_x[runup_index],wet_y[runup_index]]
     496                maximal_runup_location =    [wet_x[runup_index], \
     497                                            wet_y[runup_index]]
    500498
    501499    return maximal_runup, maximal_runup_location
    502500
    503 
    504 
    505     def test_get_maximum_inundation_from_sww(self):
    506         """test_get_maximum_inundation_from_sww(self)
    507 
    508         Test of get_maximum_inundation_elevation()
    509         and get_maximum_inundation_location() from data_manager.py
    510 
    511         This is based on test_get_maximum_inundation_3(self) but works with the
    512         stored results instead of with the internal data structure.
    513 
    514         This test uses the underlying get_maximum_inundation_data for tests
    515         """
    516 
    517         from anuga.abstract_2d_finite_volumes.mesh_factory \
    518                 import rectangular_cross
    519         from data_manager import get_maximum_inundation_elevation
    520         from data_manager import get_maximum_inundation_location
    521         from data_manager import get_maximum_inundation_data
    522 
    523         initial_runup_height = -0.4
    524         final_runup_height = -0.3
    525 
    526         #--------------------------------------------------------------
    527         # Setup computational domain
    528         #--------------------------------------------------------------
    529         N = 10
    530         points, vertices, boundary = rectangular_cross(N, N)
    531         domain = Domain(points, vertices, boundary)
    532         domain.set_name('runup_test')
    533         domain.set_maximum_allowed_speed(1.0)
    534 
    535         # FIXME: This works better with old limiters so far
    536         domain.tight_slope_limiters = 0
    537 
    538         #--------------------------------------------------------------
    539         # Setup initial conditions
    540         #--------------------------------------------------------------
    541         def topography(x, y):
    542             return -x/2                             # linear bed slope
    543 
    544         # Use function for elevation
    545         domain.set_quantity('elevation', topography)
    546         domain.set_quantity('friction', 0.)                # Zero friction
    547         # Constant negative initial stage
    548         domain.set_quantity('stage', initial_runup_height)
    549 
    550         #--------------------------------------------------------------
    551         # Setup boundary conditions
    552         #--------------------------------------------------------------
    553         Br = Reflective_boundary(domain)                       # Reflective wall
    554         Bd = Dirichlet_boundary([final_runup_height, 0, 0])    # Constant inflow
    555 
    556         # All reflective to begin with (still water)
    557         domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
    558 
    559         #--------------------------------------------------------------
    560         # Test initial inundation height
    561         #--------------------------------------------------------------
    562         indices = domain.get_wet_elements()
    563         z = domain.get_quantity('elevation').\
    564                 get_values(location='centroids', indices=indices)
    565         assert num.alltrue(z < initial_runup_height)
    566 
    567         q_ref = domain.get_maximum_inundation_elevation()
    568         # First order accuracy
    569         assert num.allclose(q_ref, initial_runup_height, rtol=1.0/N)
    570 
    571         #--------------------------------------------------------------
    572         # Let triangles adjust
    573         #--------------------------------------------------------------
    574         for t in domain.evolve(yieldstep = 0.1, finaltime = 1.0):
    575             pass
    576 
    577         #--------------------------------------------------------------
    578         # Test inundation height again
    579         #--------------------------------------------------------------
    580         q_ref = domain.get_maximum_inundation_elevation()
    581         q = get_maximum_inundation_elevation('runup_test.sww')
    582         msg = 'We got %f, should have been %f' % (q, q_ref)
    583         assert num.allclose(q, q_ref, rtol=1.0/N), msg
    584 
    585         q = get_maximum_inundation_elevation('runup_test.sww')
    586         msg = 'We got %f, should have been %f' % (q, initial_runup_height)
    587         assert num.allclose(q, initial_runup_height, rtol = 1.0/N), msg
    588 
    589         # Test error condition if time interval is out
    590         try:
    591             q = get_maximum_inundation_elevation('runup_test.sww',
    592                                                  time_interval=[2.0, 3.0])
    593         except ValueError:
    594             pass
    595         else:
    596             msg = 'should have caught wrong time interval'
    597             raise Exception, msg
    598 
    599         # Check correct time interval
    600         q, loc = get_maximum_inundation_data('runup_test.sww',
    601                                              time_interval=[0.0, 3.0])
    602         msg = 'We got %f, should have been %f' % (q, initial_runup_height)
    603         assert num.allclose(q, initial_runup_height, rtol = 1.0/N), msg
    604         assert num.allclose(-loc[0]/2, q)    # From topography formula
    605 
    606         #--------------------------------------------------------------
    607         # Update boundary to allow inflow
    608         #--------------------------------------------------------------
    609         domain.set_boundary({'right': Bd})
    610 
    611         #--------------------------------------------------------------
    612         # Evolve system through time
    613         #--------------------------------------------------------------
    614         q_max = None
    615         for t in domain.evolve(yieldstep = 0.1, finaltime = 3.0,
    616                                skip_initial_step = True):
    617             q = domain.get_maximum_inundation_elevation()
    618             if q > q_max:
    619                 q_max = q
    620 
    621         #--------------------------------------------------------------
    622         # Test inundation height again
    623         #--------------------------------------------------------------
    624         indices = domain.get_wet_elements()
    625         z = domain.get_quantity('elevation').\
    626                 get_values(location='centroids', indices=indices)
    627 
    628         assert num.alltrue(z < final_runup_height)
    629 
    630         q = domain.get_maximum_inundation_elevation()
    631         # First order accuracy
    632         assert num.allclose(q, final_runup_height, rtol=1.0/N)
    633 
    634         q, loc = get_maximum_inundation_data('runup_test.sww',
    635                                              time_interval=[3.0, 3.0])
    636         msg = 'We got %f, should have been %f' % (q, final_runup_height)
    637         assert num.allclose(q, final_runup_height, rtol=1.0/N), msg
    638         assert num.allclose(-loc[0]/2, q)    # From topography formula
    639 
    640         q = get_maximum_inundation_elevation('runup_test.sww')
    641         loc = get_maximum_inundation_location('runup_test.sww')
    642         msg = 'We got %f, should have been %f' % (q, q_max)
    643         assert num.allclose(q, q_max, rtol=1.0/N), msg
    644         assert num.allclose(-loc[0]/2, q)    # From topography formula
    645 
    646         q = get_maximum_inundation_elevation('runup_test.sww',
    647                                              time_interval=[0, 3])
    648         msg = 'We got %f, should have been %f' % (q, q_max)
    649         assert num.allclose(q, q_max, rtol=1.0/N), msg
    650 
    651         # Check polygon mode
    652         # Runup region
    653         polygon = [[0.3, 0.0], [0.9, 0.0], [0.9, 1.0], [0.3, 1.0]]
    654         q = get_maximum_inundation_elevation('runup_test.sww',
    655                                              polygon = polygon,
    656                                              time_interval=[0, 3])
    657         msg = 'We got %f, should have been %f' % (q, q_max)
    658         assert num.allclose(q, q_max, rtol=1.0/N), msg
    659 
    660         # Offshore region
    661         polygon = [[0.9, 0.0], [1.0, 0.0], [1.0, 1.0], [0.9, 1.0]]
    662         q, loc = get_maximum_inundation_data('runup_test.sww',
    663                                              polygon = polygon,
    664                                              time_interval=[0, 3])
    665         msg = 'We got %f, should have been %f' % (q, -0.475)
    666         assert num.allclose(q, -0.475, rtol=1.0/N), msg
    667         assert is_inside_polygon(loc, polygon)
    668         assert num.allclose(-loc[0]/2, q)    # From topography formula
    669 
    670         # Dry region
    671         polygon = [[0.0, 0.0], [0.4, 0.0], [0.4, 1.0], [0.0, 1.0]]
    672         q, loc = get_maximum_inundation_data('runup_test.sww',
    673                                              polygon = polygon,
    674                                              time_interval=[0, 3])
    675         msg = 'We got %s, should have been None' % (q)
    676         assert q is None, msg
    677         msg = 'We got %s, should have been None' % (loc)
    678         assert loc is None, msg
    679 
    680         # Check what happens if no time point is within interval
    681         try:
    682             q = get_maximum_inundation_elevation('runup_test.sww',
    683                                                  time_interval=[2.75, 2.75])
    684         except AssertionError:
    685             pass
    686         else:
    687             msg = 'Time interval should have raised an exception'
    688             raise Exception, msg
    689 
    690         # Cleanup
    691         try:
    692             os.remove(domain.get_name() + '.sww')
    693         except:
    694             pass
    695             #FIXME(Ole): Windows won't allow removal of this
    696 
  • trunk/anuga_core/source/anuga/shallow_water/test_data_manager.py

    r7776 r7778  
    3333from anuga.config import netcdf_float, epsilon, g
    3434
    35 from anuga.file.csv_file import load_csv_as_dict, load_csv_as_array
     35from anuga.file.csv_file import load_csv_as_dict, load_csv_as_array, \
     36                                load_csv_as_building_polygons, \
     37                                load_csv_as_polygons
    3638from anuga.file.sts import create_sts_boundary
    3739
     
    12881290
    12891291
    1290 
    1291         ##################################
    1292         #Import the file as a new domain
    1293         ##################################
    1294         from data_manager import load_sww_as_domain
    1295         import os
    1296 
    12971292        filename = domain.datadir + os.sep + domain.get_name() + '.sww'
    12981293
     
    14861481            assert num.allclose(eval('domain.'+bit),eval('domain2.'+bit))
    14871482
    1488 
    1489     def test_decimate_dem(self):
    1490         """Test decimation of dem file
    1491         """
    1492 
    1493         import os
    1494         from Scientific.IO.NetCDF import NetCDFFile
    1495 
    1496         #Write test dem file
    1497         root = 'decdemtest'
    1498 
    1499         filename = root + '.dem'
    1500         fid = NetCDFFile(filename, netcdf_mode_w)
    1501 
    1502         fid.institution = 'Geoscience Australia'
    1503         fid.description = 'NetCDF DEM format for compact and portable ' +\
    1504                           'storage of spatial point data'
    1505 
    1506         nrows = 15
    1507         ncols = 18
    1508 
    1509         fid.ncols = ncols
    1510         fid.nrows = nrows
    1511         fid.xllcorner = 2000.5
    1512         fid.yllcorner = 3000.5
    1513         fid.cellsize = 25
    1514         fid.NODATA_value = -9999
    1515 
    1516         fid.zone = 56
    1517         fid.false_easting = 0.0
    1518         fid.false_northing = 0.0
    1519         fid.projection = 'UTM'
    1520         fid.datum = 'WGS84'
    1521         fid.units = 'METERS'
    1522 
    1523         fid.createDimension('number_of_points', nrows*ncols)
    1524 
    1525         fid.createVariable('elevation', netcdf_float, ('number_of_points',))
    1526 
    1527         elevation = fid.variables['elevation']
    1528 
    1529         elevation[:] = (num.arange(nrows*ncols))
    1530 
    1531         fid.close()
    1532 
    1533         #generate the elevation values expected in the decimated file
    1534         ref_elevation = [(  0+  1+  2+ 18+ 19+ 20+ 36+ 37+ 38) / 9.0,
    1535                          (  4+  5+  6+ 22+ 23+ 24+ 40+ 41+ 42) / 9.0,
    1536                          (  8+  9+ 10+ 26+ 27+ 28+ 44+ 45+ 46) / 9.0,
    1537                          ( 12+ 13+ 14+ 30+ 31+ 32+ 48+ 49+ 50) / 9.0,
    1538                          ( 72+ 73+ 74+ 90+ 91+ 92+108+109+110) / 9.0,
    1539                          ( 76+ 77+ 78+ 94+ 95+ 96+112+113+114) / 9.0,
    1540                          ( 80+ 81+ 82+ 98+ 99+100+116+117+118) / 9.0,
    1541                          ( 84+ 85+ 86+102+103+104+120+121+122) / 9.0,
    1542                          (144+145+146+162+163+164+180+181+182) / 9.0,
    1543                          (148+149+150+166+167+168+184+185+186) / 9.0,
    1544                          (152+153+154+170+171+172+188+189+190) / 9.0,
    1545                          (156+157+158+174+175+176+192+193+194) / 9.0,
    1546                          (216+217+218+234+235+236+252+253+254) / 9.0,
    1547                          (220+221+222+238+239+240+256+257+258) / 9.0,
    1548                          (224+225+226+242+243+244+260+261+262) / 9.0,
    1549                          (228+229+230+246+247+248+264+265+266) / 9.0]
    1550 
    1551         # generate a stencil for computing the decimated values
    1552         stencil = num.ones((3,3), num.float) / 9.0
    1553 
    1554         decimate_dem(root, stencil=stencil, cellsize_new=100)
    1555 
    1556         # Open decimated NetCDF file
    1557         fid = NetCDFFile(root + '_100.dem', netcdf_mode_r)
    1558 
    1559         # Get decimated elevation
    1560         elevation = fid.variables['elevation']
    1561 
    1562         # Check values
    1563         assert num.allclose(elevation, ref_elevation)
    1564 
    1565         # Cleanup
    1566         fid.close()
    1567 
    1568         os.remove(root + '.dem')
    1569         os.remove(root + '_100.dem')
    1570 
    1571     def test_decimate_dem_NODATA(self):
    1572         """Test decimation of dem file that includes NODATA values
    1573         """
    1574 
    1575         import os
    1576         from Scientific.IO.NetCDF import NetCDFFile
    1577 
    1578         # Write test dem file
    1579         root = 'decdemtest'
    1580 
    1581         filename = root + '.dem'
    1582         fid = NetCDFFile(filename, netcdf_mode_w)
    1583 
    1584         fid.institution = 'Geoscience Australia'
    1585         fid.description = 'NetCDF DEM format for compact and portable ' +\
    1586                           'storage of spatial point data'
    1587 
    1588         nrows = 15
    1589         ncols = 18
    1590         NODATA_value = -9999
    1591 
    1592         fid.ncols = ncols
    1593         fid.nrows = nrows
    1594         fid.xllcorner = 2000.5
    1595         fid.yllcorner = 3000.5
    1596         fid.cellsize = 25
    1597         fid.NODATA_value = NODATA_value
    1598 
    1599         fid.zone = 56
    1600         fid.false_easting = 0.0
    1601         fid.false_northing = 0.0
    1602         fid.projection = 'UTM'
    1603         fid.datum = 'WGS84'
    1604         fid.units = 'METERS'
    1605 
    1606         fid.createDimension('number_of_points', nrows*ncols)
    1607 
    1608         fid.createVariable('elevation', netcdf_float, ('number_of_points',))
    1609 
    1610         elevation = fid.variables['elevation']
    1611 
    1612         # Generate initial elevation values
    1613         elevation_tmp = (num.arange(nrows*ncols))
    1614 
    1615         # Add some NODATA values
    1616         elevation_tmp[0]   = NODATA_value
    1617         elevation_tmp[95]  = NODATA_value
    1618         elevation_tmp[188] = NODATA_value
    1619         elevation_tmp[189] = NODATA_value
    1620         elevation_tmp[190] = NODATA_value
    1621         elevation_tmp[209] = NODATA_value
    1622         elevation_tmp[252] = NODATA_value
    1623 
    1624         elevation[:] = elevation_tmp
    1625 
    1626         fid.close()
    1627 
    1628         # Generate the elevation values expected in the decimated file
    1629         ref_elevation = [NODATA_value,
    1630                          (  4+  5+  6+ 22+ 23+ 24+ 40+ 41+ 42) / 9.0,
    1631                          (  8+  9+ 10+ 26+ 27+ 28+ 44+ 45+ 46) / 9.0,
    1632                          ( 12+ 13+ 14+ 30+ 31+ 32+ 48+ 49+ 50) / 9.0,
    1633                          ( 72+ 73+ 74+ 90+ 91+ 92+108+109+110) / 9.0,
    1634                          NODATA_value,
    1635                          ( 80+ 81+ 82+ 98+ 99+100+116+117+118) / 9.0,
    1636                          ( 84+ 85+ 86+102+103+104+120+121+122) / 9.0,
    1637                          (144+145+146+162+163+164+180+181+182) / 9.0,
    1638                          (148+149+150+166+167+168+184+185+186) / 9.0,
    1639                          NODATA_value,
    1640                          (156+157+158+174+175+176+192+193+194) / 9.0,
    1641                          NODATA_value,
    1642                          (220+221+222+238+239+240+256+257+258) / 9.0,
    1643                          (224+225+226+242+243+244+260+261+262) / 9.0,
    1644                          (228+229+230+246+247+248+264+265+266) / 9.0]
    1645 
    1646         # Generate a stencil for computing the decimated values
    1647         stencil = num.ones((3,3), num.float) / 9.0
    1648 
    1649         decimate_dem(root, stencil=stencil, cellsize_new=100)
    1650 
    1651         # Open decimated NetCDF file
    1652         fid = NetCDFFile(root + '_100.dem', netcdf_mode_r)
    1653 
    1654         # Get decimated elevation
    1655         elevation = fid.variables['elevation']
    1656 
    1657         # Check values
    1658         assert num.allclose(elevation, ref_elevation)
    1659 
    1660         # Cleanup
    1661         fid.close()
    1662 
    1663         os.remove(root + '.dem')
    1664         os.remove(root + '_100.dem')     
    1665 
    1666        
    16671483    def test_file_boundary_stsIV_sinewave_ordering(self):
    16681484        """test_file_boundary_stsIV_sinewave_ordering(self):
     
    16711487        """
    16721488       
    1673         from anuga.shallow_water import Domain
     1489        from anuga.shallow_water.shallow_water_domain import Domain
    16741490        from anuga.shallow_water import Reflective_boundary
    16751491        from anuga.shallow_water import Dirichlet_boundary
     
    18671683        """
    18681684       
    1869         from anuga.shallow_water import Domain
     1685        from anuga.shallow_water.shallow_water_domain import Domain
    18701686        from anuga.shallow_water import Reflective_boundary
    18711687        from anuga.shallow_water import Dirichlet_boundary
  • trunk/anuga_core/source/anuga/shallow_water/test_forcing.py

    r7744 r7778  
    33
    44import unittest, os
    5 from anuga.shallow_water import Domain
     5from anuga.shallow_water.shallow_water_domain import Domain
    66from boundaries import Reflective_boundary
    77from anuga.coordinate_transforms.geo_reference import Geo_reference
  • trunk/anuga_core/source/anuga/shallow_water/test_shallow_water_domain.py

    r7771 r7778  
    55from math import pi, sqrt
    66import tempfile
    7 
    8 from anuga.shallow_water import Domain
    97
    108from Scientific.IO.NetCDF import NetCDFFile
     
    3129from anuga.shallow_water.sww_interrogate import get_flow_through_cross_section
    3230
     31from shallow_water_domain import Domain
     32
    3333# boundary functions
    34 from anuga.shallow_water.boundaries import Reflective_boundary, \
     34from boundaries import Reflective_boundary, \
    3535            Field_boundary, Transmissive_momentum_set_stage_boundary, \
    3636            Transmissive_stage_zero_momentum_boundary
     
    15481548        from anuga.abstract_2d_finite_volumes.mesh_factory \
    15491549                import rectangular_cross
    1550         from anuga.shallow_water import Domain
     1550        from anuga.shallow_water.shallow_water_domain import Domain
    15511551        from anuga.shallow_water import Reflective_boundary
    15521552        from anuga.shallow_water import Dirichlet_boundary
     
    53595359        verbose = False
    53605360
    5361         from anuga.shallow_water import Domain
     5361        from anuga.shallow_water.shallow_water_domain import Domain
    53625362        from anuga.pmesh.mesh_interface import create_mesh_from_regions
    53635363        from anuga.geospatial_data.geospatial_data import Geospatial_data
     
    54725472        verbose = False       
    54735473
    5474         from anuga.shallow_water import Domain
     5474        from anuga.shallow_water.shallow_water_domain import Domain
    54755475        from anuga.pmesh.mesh_interface import create_mesh_from_regions
    54765476        from anuga.geospatial_data import Geospatial_data
     
    55365536        from anuga.abstract_2d_finite_volumes.mesh_factory \
    55375537                import rectangular_cross
    5538         from anuga.shallow_water import Domain
     5538        from anuga.shallow_water.shallow_water_domain import Domain
    55395539
    55405540        #----------------------------------------------------------------------
     
    56015601        from anuga.abstract_2d_finite_volumes.mesh_factory \
    56025602                import rectangular_cross
    5603         from anuga.shallow_water import Domain
     5603        from anuga.shallow_water.shallow_water_domain import Domain
    56045604        from anuga.shallow_water.forcing import Inflow
    56055605
     
    56905690        #---------------------------------------------------------------------
    56915691        from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
    5692         from anuga.shallow_water import Domain
     5692        from anuga.shallow_water.shallow_water_domain import Domain
    56935693
    56945694        #----------------------------------------------------------------------
     
    58875887        #---------------------------------------------------------------------
    58885888        from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
    5889         from anuga.shallow_water import Domain
     5889        from anuga.shallow_water.shallow_water_domain import Domain
    58905890        from anuga.shallow_water.boundaries import Reflective_boundary
    58915891        from anuga.shallow_water.shallow_water_domain import Dirichlet_boundary
     
    60036003        #---------------------------------------------------------------------
    60046004        from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
    6005         from anuga.shallow_water import Domain
     6005        from anuga.shallow_water.shallow_water_domain import Domain
    60066006        from anuga.shallow_water import Reflective_boundary
    60076007        from anuga.shallow_water import Dirichlet_boundary
     
    62546254       
    62556255
    6256         #----------------------------------------------------------------------
    6257         # Import necessary modules
    6258         #----------------------------------------------------------------------
    6259         from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
    6260         from anuga.shallow_water.data_manager import get_flow_through_cross_section
    6261         from anuga.abstract_2d_finite_volumes.util import sww2csv_gauges, csv2timeseries_graphs
    6262 
    62636256
    62646257        #----------------------------------------------------------------------
     
    63916384        from anuga.abstract_2d_finite_volumes.mesh_factory \
    63926385                import rectangular_cross
    6393         from anuga.shallow_water import Domain
     6386        from anuga.shallow_water.shallow_water_domain import Domain
    63946387        from anuga.shallow_water.shallow_water_domain import Reflective_boundary
    63956388        from anuga.shallow_water.shallow_water_domain import Dirichlet_boundary
  • trunk/anuga_core/source/anuga/shallow_water/test_sww_interrogate.py

    r7771 r7778  
    694694            pass
    695695            #FIXME(Ole): Windows won't allow removal of this
    696        
     696
     697
     698    def test_get_maximum_inundation_from_sww(self):
     699        """test_get_maximum_inundation_from_sww(self)
     700
     701        Test of get_maximum_inundation_elevation()
     702        and get_maximum_inundation_location().
     703   
     704        This is based on test_get_maximum_inundation_3(self) but works with the
     705        stored results instead of with the internal data structure.
     706
     707        This test uses the underlying get_maximum_inundation_data for tests
     708        """
     709
     710        initial_runup_height = -0.4
     711        final_runup_height = -0.3
     712
     713        #--------------------------------------------------------------
     714        # Setup computational domain
     715        #--------------------------------------------------------------
     716        N = 10
     717        points, vertices, boundary = rectangular_cross(N, N)
     718        domain = Domain(points, vertices, boundary)
     719        domain.set_name('runup_test')
     720        domain.set_maximum_allowed_speed(1.0)
     721
     722        # FIXME: This works better with old limiters so far
     723        domain.tight_slope_limiters = 0
     724
     725        #--------------------------------------------------------------
     726        # Setup initial conditions
     727        #--------------------------------------------------------------
     728        def topography(x, y):
     729            return -x/2                             # linear bed slope
     730
     731        # Use function for elevation
     732        domain.set_quantity('elevation', topography)
     733        domain.set_quantity('friction', 0.)                # Zero friction
     734        # Constant negative initial stage
     735        domain.set_quantity('stage', initial_runup_height)
     736
     737        #--------------------------------------------------------------
     738        # Setup boundary conditions
     739        #--------------------------------------------------------------
     740        Br = Reflective_boundary(domain)                       # Reflective wall
     741        Bd = Dirichlet_boundary([final_runup_height, 0, 0])    # Constant inflow
     742
     743        # All reflective to begin with (still water)
     744        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
     745
     746        #--------------------------------------------------------------
     747        # Test initial inundation height
     748        #--------------------------------------------------------------
     749        indices = domain.get_wet_elements()
     750        z = domain.get_quantity('elevation').\
     751                get_values(location='centroids', indices=indices)
     752        assert num.alltrue(z < initial_runup_height)
     753
     754        q_ref = domain.get_maximum_inundation_elevation()
     755        # First order accuracy
     756        assert num.allclose(q_ref, initial_runup_height, rtol=1.0/N)
     757
     758        #--------------------------------------------------------------
     759        # Let triangles adjust
     760        #--------------------------------------------------------------
     761        for t in domain.evolve(yieldstep = 0.1, finaltime = 1.0):
     762            pass
     763
     764        #--------------------------------------------------------------
     765        # Test inundation height again
     766        #--------------------------------------------------------------
     767        q_ref = domain.get_maximum_inundation_elevation()
     768        q = get_maximum_inundation_elevation('runup_test.sww')
     769        msg = 'We got %f, should have been %f' % (q, q_ref)
     770        assert num.allclose(q, q_ref, rtol=1.0/N), msg
     771
     772        q = get_maximum_inundation_elevation('runup_test.sww')
     773        msg = 'We got %f, should have been %f' % (q, initial_runup_height)
     774        assert num.allclose(q, initial_runup_height, rtol = 1.0/N), msg
     775
     776        # Test error condition if time interval is out
     777        try:
     778            q = get_maximum_inundation_elevation('runup_test.sww',
     779                                                 time_interval=[2.0, 3.0])
     780        except ValueError:
     781            pass
     782        else:
     783            msg = 'should have caught wrong time interval'
     784            raise Exception, msg
     785
     786        # Check correct time interval
     787        q, loc = get_maximum_inundation_data('runup_test.sww',
     788                                             time_interval=[0.0, 3.0])
     789        msg = 'We got %f, should have been %f' % (q, initial_runup_height)
     790        assert num.allclose(q, initial_runup_height, rtol = 1.0/N), msg
     791        assert num.allclose(-loc[0]/2, q)    # From topography formula
     792
     793        #--------------------------------------------------------------
     794        # Update boundary to allow inflow
     795        #--------------------------------------------------------------
     796        domain.set_boundary({'right': Bd})
     797
     798        #--------------------------------------------------------------
     799        # Evolve system through time
     800        #--------------------------------------------------------------
     801        q_max = None
     802        for t in domain.evolve(yieldstep = 0.1, finaltime = 3.0,
     803                               skip_initial_step = True):
     804            q = domain.get_maximum_inundation_elevation()
     805            if q > q_max:
     806                q_max = q
     807
     808        #--------------------------------------------------------------
     809        # Test inundation height again
     810        #--------------------------------------------------------------
     811        indices = domain.get_wet_elements()
     812        z = domain.get_quantity('elevation').\
     813                get_values(location='centroids', indices=indices)
     814
     815        assert num.alltrue(z < final_runup_height)
     816
     817        q = domain.get_maximum_inundation_elevation()
     818        # First order accuracy
     819        assert num.allclose(q, final_runup_height, rtol=1.0/N)
     820
     821        q, loc = get_maximum_inundation_data('runup_test.sww',
     822                                             time_interval=[3.0, 3.0])
     823        msg = 'We got %f, should have been %f' % (q, final_runup_height)
     824        assert num.allclose(q, final_runup_height, rtol=1.0/N), msg
     825        assert num.allclose(-loc[0]/2, q)    # From topography formula
     826
     827        q = get_maximum_inundation_elevation('runup_test.sww')
     828        loc = get_maximum_inundation_location('runup_test.sww')
     829        msg = 'We got %f, should have been %f' % (q, q_max)
     830        assert num.allclose(q, q_max, rtol=1.0/N), msg
     831        assert num.allclose(-loc[0]/2, q)    # From topography formula
     832
     833        q = get_maximum_inundation_elevation('runup_test.sww',
     834                                             time_interval=[0, 3])
     835        msg = 'We got %f, should have been %f' % (q, q_max)
     836        assert num.allclose(q, q_max, rtol=1.0/N), msg
     837
     838        # Check polygon mode
     839        # Runup region
     840        polygon = [[0.3, 0.0], [0.9, 0.0], [0.9, 1.0], [0.3, 1.0]]
     841        q = get_maximum_inundation_elevation('runup_test.sww',
     842                                             polygon = polygon,
     843                                             time_interval=[0, 3])
     844        msg = 'We got %f, should have been %f' % (q, q_max)
     845        assert num.allclose(q, q_max, rtol=1.0/N), msg
     846
     847        # Offshore region
     848        polygon = [[0.9, 0.0], [1.0, 0.0], [1.0, 1.0], [0.9, 1.0]]
     849        q, loc = get_maximum_inundation_data('runup_test.sww',
     850                                             polygon = polygon,
     851                                             time_interval=[0, 3])
     852        msg = 'We got %f, should have been %f' % (q, -0.475)
     853        assert num.allclose(q, -0.475, rtol=1.0/N), msg
     854        assert is_inside_polygon(loc, polygon)
     855        assert num.allclose(-loc[0]/2, q)    # From topography formula
     856
     857        # Dry region
     858        polygon = [[0.0, 0.0], [0.4, 0.0], [0.4, 1.0], [0.0, 1.0]]
     859        q, loc = get_maximum_inundation_data('runup_test.sww',
     860                                             polygon = polygon,
     861                                             time_interval=[0, 3])
     862        msg = 'We got %s, should have been None' % (q)
     863        assert q is None, msg
     864        msg = 'We got %s, should have been None' % (loc)
     865        assert loc is None, msg
     866
     867        # Check what happens if no time point is within interval
     868        try:
     869            q = get_maximum_inundation_elevation('runup_test.sww',
     870                                                 time_interval=[2.75, 2.75])
     871        except AssertionError:
     872            pass
     873        else:
     874            msg = 'Time interval should have raised an exception'
     875            raise Exception, msg
     876
     877        # Cleanup
     878        try:
     879            os.remove(domain.get_name() + '.sww')
     880        except:
     881            pass
     882            #FIXME(Ole): Windows won't allow removal of this
     883
     884 
     885 
    697886 
    698887if __name__ == "__main__":
  • trunk/anuga_core/source/anuga/shallow_water/test_system.py

    r7736 r7778  
    99import numpy as num
    1010
    11 from anuga.shallow_water import Domain
    12 from anuga.shallow_water import Dirichlet_boundary, Time_boundary
    13 from anuga.shallow_water import File_boundary
     11import anuga
    1412from anuga.pmesh.mesh import Mesh
    1513from anuga.abstract_2d_finite_volumes.pmesh2domain import \
     
    4745        mesh.generate_mesh(verbose=False)
    4846       
    49         domain = pmesh_to_domain_instance(mesh, Domain)
     47        domain = pmesh_to_domain_instance(mesh, anuga.Domain)
    5048        domain.set_name(boundary_name)                 
    5149        domain.set_datadir(dir)         
     
    5755
    5856        # Setup boundary conditions
    59         Bd = Dirichlet_boundary([tide,0.,0.]) # Constant boundary values
    60         Bd = Time_boundary(domain=domain,     # Time dependent boundary 
     57        Bd = anuga.Dirichlet_boundary([tide,0.,0.]) # Constant boundary values
     58        Bd = anuga.Time_boundary(domain=domain,     # Time dependent boundary 
    6159                   f=lambda t: [t, 0.0, 0.0])
    6260        domain.set_boundary({'exterior': Bd})
     
    9694        mesh.generate_mesh(verbose=False)
    9795       
    98         domain = pmesh_to_domain_instance(mesh, Domain)
     96        domain = pmesh_to_domain_instance(mesh, anuga.Domain)
    9997        domain.set_name(senario_name)                 
    10098        domain.set_datadir(dir)
     
    103101        domain.set_quantity('elevation', 0.0)
    104102        domain.set_quantity('stage', 0.0)         
    105         Bf = File_boundary(boundary_filename,
     103        Bf = anuga.File_boundary(boundary_filename,
    106104                           domain,  use_cache=False, verbose=False)
    107105
     
    146144        mesh.generate_mesh(verbose=False)
    147145       
    148         domain = pmesh_to_domain_instance(mesh, Domain)
     146        domain = pmesh_to_domain_instance(mesh, anuga.Domain)
    149147        domain.set_name(senario_name)                 
    150148        domain.set_datadir(dir)
     
    155153        domain.set_quantity('elevation', 0.0)
    156154        domain.set_quantity('stage', 0.0)         
    157         Bf = File_boundary(boundary_filename,
     155        Bf = anuga.File_boundary(boundary_filename,
    158156                           domain,  use_cache=False, verbose=False)
    159157
  • trunk/anuga_core/source/anuga/utilities/file_utils.py

    r7765 r7778  
    99import shutil
    1010
    11 ##
    12 # @brief Convert a possible filename into a standard form.
    13 # @param s Filename to process.
    14 # @return The new filename string.
    1511def make_filename(s):
    16     """Transform argument string into a Sexsuitable filename
     12    """Transform argument string into a standard filename
     13   
     14        Convert a possible filename into a standard form.
     15        s Filename to process.
     16        The new filename string.
    1717    """
    1818
     
    2626
    2727
    28 ##
    29 # @brief Check that a specified filesystem directory path exists.
    30 # @param path The dirstory path to check.
    31 # @param verbose True if this function is to be verbose.
    32 # @note If directory path doesn't exist, it will be created.
    3328def check_dir(path, verbose=None):
    3429    """Check that specified path exists.
Note: See TracChangeset for help on using the changeset viewer.