source: inundation/ga/storm_surge/analytical solutions/Sydney.py @ 647

"""Example of shallow water wave equation analytical solution
consists of a parabolic profile in a parabolic basin. Analytical
solutiuon to this problem was derived by Carrier and Greenspan
and used by Yoon and Chou.

   Copyright 2004
   Christopher Zoppou, Stephen Roberts, Ole Nielsen, Duncan Gray
   Geoscience Australia
   
"""

######################
# Module imports 
import sys
from os import sep
sys.path.append('..'+sep+'pyvolution')

from shallow_water import Transmissive_boundary, Reflective_boundary, \
     Dirichlet_boundary, Domain
from math import sqrt, cos, sin, pi
from mesh_factory import strang_mesh
from pmesh2domain import pmesh_to_domain_instance


######################
# Domain
filename = 'Sydney_UBD.tsh'
print 'Creating domain from', filename
domain = pmesh_to_domain_instance(filename, Domain)
print 'Number of triangles = ', len(domain)


domain.default_order = 2
domain.smooth = True


# Provide file name for storing output
domain.store = True
domain.format = 'sww'
domain.filename = 'sydney'

print 'Number of triangles = ', len(domain)

#Reduction operation for get_vertex_values
from util import mean
domain.reduction = mean 
#domain.reduction = min  #Looks better near steep slopes


######################
#Initial condition
#
print 'Initial condition'

#Set bed-elevation and friction(None)
def x_slope(x,y):
    n = x.shape[0]
    z = 0*x
    return z

domain.set_quantity('elevation', x_slope)

#Set the initial water level
def level(x,y):
    z = x_slope(x,y)
    n = x.shape[0]
    h = 0*x
    return h   

domain.set_quantity('level', level)


############
#Boundary
tags = {}
#tags['input'] = Dirichlet_boundary([10.0, 0.0, 0.0])
tags['reflective'] = Reflective_boundary(domain) 
#tags['far'] = Transmissive_boundary(domain)
tags['far'] = Dirichlet_boundary([10.0, 0.0, 0.0])
domain.set_boundary(tags)
    

######################
#Evolution
import time
t0 = time.time()
for t in domain.evolve(yieldstep = 10.0, finaltime = 5):
    domain.write_time()

print 'That took %.2f seconds' %(time.time()-t0)