本文整理匯總了Python中anuga.shallow_water.shallow_water_domain.Domain.get_quantity方法的典型用法代碼示例。如果您正苦於以下問題：Python Domain.get_quantity方法的具體用法？Python Domain.get_quantity怎麽用？Python Domain.get_quantity使用的例子？那麽, 這裏精選的方法代碼示例或許可以為您提供幫助。您也可以進一步了解該方法所在類anuga.shallow_water.shallow_water_domain.Domain的用法示例。
在下文中一共展示了Domain.get_quantity方法的5個代碼示例，這些例子默認根據受歡迎程度排序。您可以為喜歡或者感覺有用的代碼點讚，您的評價將有助於係統推薦出更棒的Python代碼示例。
示例1: test_read_sww
# 需要導入模塊: from anuga.shallow_water.shallow_water_domain import Domain [as 別名]
# 或者: from anuga.shallow_water.shallow_water_domain.Domain import get_quantity [as 別名]
def test_read_sww(self):
"""
Save to an sww file and then read back the info.
Here we store the info "uniquely"
"""
# ---------------------------------------------------------------------
# Import necessary modules
# ---------------------------------------------------------------------
from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
from anuga.shallow_water.shallow_water_domain import Domain
from anuga import Reflective_boundary
from anuga.abstract_2d_finite_volumes.generic_boundary_conditions import Dirichlet_boundary, Time_boundary
# ---------------------------------------------------------------------
# Setup computational domain
# ---------------------------------------------------------------------
length = 8.0
width = 4.0
dx = dy = 2 # Resolution: Length of subdivisions on both axes
inc = 0.05 # Elevation increment
points, vertices, boundary = rectangular_cross(int(length / dx), int(width / dy), len1=length, len2=width)
domain = Domain(points, vertices, boundary)
domain.set_name("read_sww_test" + str(domain.processor)) # Output name
domain.set_quantities_to_be_stored({"elevation": 2, "stage": 2, "xmomentum": 2, "ymomentum": 2, "friction": 1})
domain.set_store_vertices_uniquely(True)
# ---------------------------------------------------------------------
# Setup initial conditions
# ---------------------------------------------------------------------
domain.set_quantity("elevation", 0.0) # Flat bed initially
domain.set_quantity("friction", 0.01) # Constant friction
domain.set_quantity("stage", 0.0) # Dry initial condition
# ------------------------------------------------------------------
# Setup boundary conditions
# ------------------------------------------------------------------
Bi = Dirichlet_boundary([0.4, 0, 0]) # Inflow
Br = Reflective_boundary(domain) # Solid reflective wall
Bo = Dirichlet_boundary([-5, 0, 0]) # Outflow
domain.set_boundary({"left": Bi, "right": Bo, "top": Br, "bottom": Br})
# -------------------------------------------------------------------
# Evolve system through time
# -------------------------------------------------------------------
for t in domain.evolve(yieldstep=1, finaltime=4.0):
pass
# Check that quantities have been stored correctly
source = domain.get_name() + ".sww"
# x = fid.variables['x'][:]
# y = fid.variables['y'][:]
# stage = fid.variables['stage'][:]
# elevation = fid.variables['elevation'][:]
# fid.close()
# assert len(stage.shape) == 2
# assert len(elevation.shape) == 2
# M, N = stage.shape
sww_file = sww.Read_sww(source)
# print 'last frame number',sww_file.get_last_frame_number()
assert num.allclose(sww_file.x, domain.get_vertex_coordinates()[:, 0])
assert num.allclose(sww_file.y, domain.get_vertex_coordinates()[:, 1])
assert num.allclose(sww_file.time, [0.0, 1.0, 2.0, 3.0, 4.0])
M = domain.get_number_of_triangles()
assert num.allclose(num.reshape(num.arange(3 * M), (M, 3)), sww_file.vertices)
last_frame_number = sww_file.get_last_frame_number()
assert last_frame_number == 4
assert num.allclose(sww_file.get_bounds(), [0.0, length, 0.0, width])
assert "stage" in sww_file.quantities.keys()
assert "friction" in sww_file.quantities.keys()
assert "elevation" in sww_file.quantities.keys()
assert "xmomentum" in sww_file.quantities.keys()
assert "ymomentum" in sww_file.quantities.keys()
for qname, q in sww_file.read_quantities(last_frame_number).items():
# print qname
# print num.linalg.norm(num.abs((domain.get_quantity(qname).get_values()-q).flatten()), ord=1)
assert num.allclose(domain.get_quantity(qname).get_values(), q)
# -----------------------------------------
# Start the evolution off again at frame 3
#.........這裏部分代碼省略.........
示例2: test_get_mesh_and_quantities_from_unique_vertices_DE0_sww_file
# 需要導入模塊: from anuga.shallow_water.shallow_water_domain import Domain [as 別名]
# 或者: from anuga.shallow_water.shallow_water_domain.Domain import get_quantity [as 別名]
def test_get_mesh_and_quantities_from_unique_vertices_DE0_sww_file(self):
"""test_get_mesh_and_quantities_from_unique_vertices_sww_file(self):
"""
# Generate a test sww file with non trivial georeference
import time, os
# Setup
#from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular
# Create basic mesh (100m x 5m)
width = 5
length = 50
t_end = 10
points, vertices, boundary = rectangular(10, 1, length, width)
# Create shallow water domain
domain = Domain(points, vertices, boundary,
geo_reference = Geo_reference(56,308500,6189000))
domain.set_name('test_get_mesh_and_quantities_from_unique_vertices_sww_file')
swwfile = domain.get_name() + '.sww'
domain.set_datadir('.')
domain.set_flow_algorithm('DE0')
domain.set_store_vertices_uniquely()
Br = Reflective_boundary(domain) # Side walls
Bd = Dirichlet_boundary([1, 0, 0]) # inflow
domain.set_boundary( {'left': Bd, 'right': Bd, 'top': Br, 'bottom': Br})
for t in domain.evolve(yieldstep=1, finaltime = t_end):
pass
# Read it
# Get mesh and quantities from sww file
X = get_mesh_and_quantities_from_file(swwfile,
quantities=['elevation',
'stage',
'xmomentum',
'ymomentum'],
verbose=False)
mesh, quantities, time = X
#print quantities
#print time
dhash = domain.get_nodes()[:,0]*10+domain.get_nodes()[:,1]
mhash = mesh.nodes[:,0]*10+mesh.nodes[:,1]
#print 'd_nodes',len(dhash)
#print 'm_nodes',len(mhash)
di = num.argsort(dhash)
mi = num.argsort(mhash)
minv = num.argsort(mi)
dinv = num.argsort(di)
#print 'd_tri',len(domain.get_triangles())
#print 'm_tri',len(mesh.triangles)
# Check that mesh has been recovered
# triangle order should be ok
assert num.allclose(mesh.nodes[mi,:],domain.get_nodes()[di,:])
assert num.alltrue(minv[mesh.triangles] == dinv[domain.get_triangles()])
# Check that time has been recovered
assert num.allclose(time, range(t_end+1))
z=domain.get_quantity('elevation').get_values(location='vertices').flatten()
assert num.allclose(quantities['elevation'], z)
for q in ['stage', 'xmomentum', 'ymomentum']:
# Get quantity at last timestep
q_ref=domain.get_quantity(q).get_values(location='vertices').flatten()
#print q,quantities[q]
q_sww=quantities[q][-1,:]
msg = 'Quantity %s failed to be recovered' %q
assert num.allclose(q_ref, q_sww, atol=1.0e-6), msg示例3: test_get_maximum_inundation_from_sww
# 需要導入模塊: from anuga.shallow_water.shallow_water_domain import Domain [as 別名]
# 或者: from anuga.shallow_water.shallow_water_domain.Domain import get_quantity [as 別名]
def test_get_maximum_inundation_from_sww(self):
"""test_get_maximum_inundation_from_sww(self)
Test of get_maximum_inundation_elevation()
and get_maximum_inundation_location().
This is based on test_get_maximum_inundation_3(self) but works with the
stored results instead of with the internal data structure.
This test uses the underlying get_maximum_inundation_data for tests
"""
verbose = False
from anuga.config import minimum_storable_height
initial_runup_height = -0.4
final_runup_height = -0.3
filename = 'runup_test_2'
#--------------------------------------------------------------
# Setup computational domain
#--------------------------------------------------------------
N = 10
points, vertices, boundary = rectangular_cross(N, N)
domain = Domain(points, vertices, boundary)
domain.set_name(filename)
domain.set_maximum_allowed_speed(1.0)
#domain.set_minimum_storable_height(1.0e-5)
domain.set_store_vertices_uniquely()
# FIXME: This works better with old limiters so far
domain.tight_slope_limiters = 0
#--------------------------------------------------------------
# Setup initial conditions
#--------------------------------------------------------------
def topography(x, y):
return -x/2 # linear bed slope
# Use function for elevation
domain.set_quantity('elevation', topography)
domain.set_quantity('friction', 0.) # Zero friction
# Constant negative initial stage
domain.set_quantity('stage', initial_runup_height)
#--------------------------------------------------------------
# Setup boundary conditions
#--------------------------------------------------------------
Br = Reflective_boundary(domain) # Reflective wall
Bd = Dirichlet_boundary([final_runup_height, 0, 0]) # Constant inflow
# All reflective to begin with (still water)
domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
#--------------------------------------------------------------
# Test initial inundation height
#--------------------------------------------------------------
indices = domain.get_wet_elements()
z = domain.get_quantity('elevation').\
get_values(location='centroids', indices=indices)
assert num.alltrue(z < initial_runup_height)
q_ref = domain.get_maximum_inundation_elevation(minimum_height=minimum_storable_height)
# First order accuracy
assert num.allclose(q_ref, initial_runup_height, rtol=1.0/N)
#--------------------------------------------------------------
# Let triangles adjust
#--------------------------------------------------------------
q_max = None
for t in domain.evolve(yieldstep = 0.1, finaltime = 1.0):
q = domain.get_maximum_inundation_elevation(minimum_height=minimum_storable_height)
if verbose:
domain.write_time()
print q
if q > q_max:
q_max = q
#--------------------------------------------------------------
# Test inundation height again
#--------------------------------------------------------------
#q_ref = domain.get_maximum_inundation_elevation()
q = get_maximum_inundation_elevation(filename+'.sww')
msg = 'We got %f, should have been %f' % (q, q_max)
assert num.allclose(q, q_max, rtol=2.0/N), msg
msg = 'We got %f, should have been %f' % (q, initial_runup_height)
assert num.allclose(q, initial_runup_height, rtol = 1.0/N), msg
# Test error condition if time interval is out
try:
q = get_maximum_inundation_elevation(filename+'.sww',
time_interval=[2.0, 3.0])
except ValueError:
pass
else:
msg = 'should have caught wrong time interval'
raise Exception, msg
#.........這裏部分代碼省略.........
示例4: test_get_mesh_and_quantities_from_de0_sww_file
# 需要導入模塊: from anuga.shallow_water.shallow_water_domain import Domain [as 別名]
# 或者: from anuga.shallow_water.shallow_water_domain.Domain import get_quantity [as 別名]
def test_get_mesh_and_quantities_from_de0_sww_file(self):
"""test_get_mesh_and_quantities_from_sww_file(self):
"""
# Generate a test sww file with non trivial georeference
import time, os
# Setup
#from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular
# Create basic mesh (100m x 5m)
width = 5
length = 50
t_end = 10
points, vertices, boundary = rectangular(length, width, 50, 5)
# Create shallow water domain
domain = Domain(points, vertices, boundary,
geo_reference = Geo_reference(56,308500,6189000))
domain.set_name('test_get_mesh_and_quantities_from_sww_file')
swwfile = domain.get_name() + '.sww'
domain.set_datadir('.')
domain.set_flow_algorithm('DE0')
Br = Reflective_boundary(domain) # Side walls
Bd = Dirichlet_boundary([1, 0, 0]) # inflow
domain.set_boundary( {'left': Bd, 'right': Bd, 'top': Br, 'bottom': Br})
for t in domain.evolve(yieldstep=1, finaltime = t_end):
pass
# Read it
# Get mesh and quantities from sww file
X = get_mesh_and_quantities_from_file(swwfile,
quantities=['elevation',
'stage',
'xmomentum',
'ymomentum'],
verbose=False)
mesh, quantities, time = X
# Check that mesh has been recovered
assert num.alltrue(mesh.triangles == domain.get_triangles())
assert num.allclose(mesh.nodes, domain.get_nodes())
# Check that time has been recovered
assert num.allclose(time, range(t_end+1))
# Check that quantities have been recovered
# (sww files use single precision)
z=domain.get_quantity('elevation').get_values(location='unique vertices')
assert num.allclose(quantities['elevation'], z)
for q in ['stage', 'xmomentum', 'ymomentum']:
# Get quantity at last timestep
q_ref=domain.get_quantity(q).get_values(location='unique vertices')
#print q,quantities[q]
q_sww=quantities[q][-1,:]
msg = 'Quantity %s failed to be recovered' %q
assert num.allclose(q_ref, q_sww, atol=1.0e-2), msg示例5: test_get_flow_through_cross_section_with_geo
# 需要導入模塊: from anuga.shallow_water.shallow_water_domain import Domain [as 別名]
# 或者: from anuga.shallow_water.shallow_water_domain.Domain import get_quantity [as 別名]
def test_get_flow_through_cross_section_with_geo(self):
"""test_get_flow_through_cross_section(self):
Test that the total flow through a cross section can be
correctly obtained at run-time from the ANUGA domain.
This test creates a flat bed with a known flow through it and tests
that the function correctly returns the expected flow.
The specifics are
e = -1 m
u = 2 m/s
h = 2 m
w = 3 m (width of channel)
q = u*h*w = 12 m^3/s
This run tries it with georeferencing and with elevation = -1
"""
# Create basic mesh (20m x 3m)
width = 3
length = 20
t_end = 1
points, vertices, boundary = rectangular(length, width, length, width)
# Create shallow water domain
domain = Domain(points, vertices, boundary,
geo_reference=Geo_reference(56, 308500, 6189000))
domain.default_order = 2
domain.set_quantities_to_be_stored(None)
e = -1.0
w = 1.0
h = w-e
u = 2.0
uh = u*h
Br = Reflective_boundary(domain) # Side walls
Bd = Dirichlet_boundary([w, uh, 0]) # 2 m/s across the 3 m inlet:
# Initial conditions
domain.set_quantity('elevation', e)
domain.set_quantity('stage', w)
domain.set_quantity('xmomentum', uh)
domain.set_boundary({'left': Bd, 'right': Bd, 'top': Br, 'bottom': Br})
# Interpolation points down the middle
I = [[0, width/2.],
[length/2., width/2.],
[length, width/2.]]
interpolation_points = domain.geo_reference.get_absolute(I)
for t in domain.evolve(yieldstep=0.1, finaltime=0.5):
# Shortcuts to quantites
stage = domain.get_quantity('stage')
xmomentum = domain.get_quantity('xmomentum')
ymomentum = domain.get_quantity('ymomentum')
# Check that quantities are they should be in the interior
w_t = stage.get_values(interpolation_points)
uh_t = xmomentum.get_values(interpolation_points)
vh_t = ymomentum.get_values(interpolation_points)
assert num.allclose(w_t, w)
assert num.allclose(uh_t, uh)
assert num.allclose(vh_t, 0.0, atol=1.0e-6)
# Check flows through the middle
for i in range(5):
x = length/2. + i*0.23674563 # Arbitrary
cross_section = [[x, 0], [x, width]]
cross_section = domain.geo_reference.get_absolute(cross_section)
Q = domain.get_flow_through_cross_section(cross_section,
verbose=False)
assert num.allclose(Q, uh*width)
import cPickle
cPickle.dump(domain, open('domain_pickle.pickle', 'w'))
domain_restored = cPickle.load(open('domain_pickle.pickle'))
for t in domain_restored.evolve(yieldstep=0.1, finaltime=1.0):
# Shortcuts to quantites
stage = domain_restored.get_quantity('stage')
xmomentum = domain_restored.get_quantity('xmomentum')
ymomentum = domain_restored.get_quantity('ymomentum')
# Check that quantities are they should be in the interior
w_t = stage.get_values(interpolation_points)
uh_t = xmomentum.get_values(interpolation_points)
vh_t = ymomentum.get_values(interpolation_points)
assert num.allclose(w_t, w)
assert num.allclose(uh_t, uh)
assert num.allclose(vh_t, 0.0, atol=1.0e-6)
#.........這裏部分代碼省略.........