source: production/merimbula_2005/run_cairns.py @ 3289

"""
Main meribula script using new interface
"""

#-------------------------------
# Module imports
#-------------------------------
import sys, os
from pyvolution.shallow_water import Domain, Reflective_boundary,\
     File_boundary, Transmissive_Momentum_Set_Stage_boundary
from pyvolution.mesh_factory import rectangular_cross
from pyvolution.pmesh2domain import pmesh_to_domain_instance
from Numeric import array, zeros, Float, allclose
import cairns_project
from caching import cache




#-------------------------------
# Domain
#-------------------------------
#print 'Creating domain from', cairns_project.mesh_filename

#domain = cache(pmesh_to_domain_instance,
#               (cairns_project.mesh_filename, Domain),
#               dependencies = [cairns_project.mesh_filename])

print 'Creating domain'
N = 50
M = 60
#r_earth = 5e6
#from math import sin, pi
#scale_y = r_earth*sin(74/180*pi)
#print scale_y
points, elements, boundary = rectangular_cross(N, M, len1=30.0, len2=16.0,
                                               origin = (145.0, -24.0))
domain = Domain(points, elements, boundary)

domain.check_integrity()
print 'Number of triangles = ', len(domain)
print 'The extent is ', domain.get_extent()



#-------------------------------
# Initial Conditions
#-------------------------------
print 'Initial values'

domain.set_quantity('elevation',
                    filename = cairns_project.bathymetry_filename[:-4] + '.xya',
                    alpha = 10.0,
                    verbose = True,
                    use_cache = True)

domain.set_quantity('friction', 0.0)
domain.set_quantity('stage', 0.0)

#-------------------------------
# Boundary conditions
#-------------------------------
print 'Boundaries'

#   Tidal cycle recorded at Eden as open
#print 'Open sea boundary condition from ',cairns_project.boundary_filename
#from pyvolution.util import file_function
#tide_function = file_function(cairns_project.boundary_filename[:-4] + '.tms', domain,
#                         verbose = True)
#Bts = Transmissive_Momentum_Set_Stage_boundary(domain, tide_function)

#   All other boundaries are reflective
Br = Reflective_boundary(domain)

#domain.set_boundary({'exterior': Br, 'open': Br})
domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})

#-------------------------------
# Setup domain runtime parameters
#-------------------------------
domain.visualise = True
domain.visualise_color_stage = True

base = os.path.basename(sys.argv[0])
domain.filename, _ = os.path.splitext(base)
domain.default_order = 2
domain.store = True    #Store for visualisation purposes
domain.smooth = False



#-------------------------------
# Evolve
#-------------------------------
import time
t0 = time.time()
yieldstep = 0.1
finaltime = 1


for t in domain.evolve(yieldstep = yieldstep, finaltime = finaltime):
    domain.write_time()

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