1 |
ahallam |
2469 |
|
2 |
|
|
######################################################## |
3 |
|
|
# |
4 |
jfenwick |
2548 |
# Copyright (c) 2003-2009 by University of Queensland |
5 |
ahallam |
2469 |
# Earth Systems Science Computational Center (ESSCC) |
6 |
|
|
# http://www.uq.edu.au/esscc |
7 |
|
|
# |
8 |
|
|
# Primary Business: Queensland, Australia |
9 |
|
|
# Licensed under the Open Software License version 3.0 |
10 |
|
|
# http://www.opensource.org/licenses/osl-3.0.php |
11 |
|
|
# |
12 |
|
|
######################################################## |
13 |
|
|
|
14 |
jfenwick |
2549 |
__copyright__="""Copyright (c) 2003-2009 by University of Queensland |
15 |
ahallam |
2469 |
Earth Systems Science Computational Center (ESSCC) |
16 |
|
|
http://www.uq.edu.au/esscc |
17 |
|
|
Primary Business: Queensland, Australia""" |
18 |
|
|
__license__="""Licensed under the Open Software License version 3.0 |
19 |
|
|
http://www.opensource.org/licenses/osl-3.0.php""" |
20 |
|
|
__url__="https://launchpad.net/escript-finley" |
21 |
|
|
|
22 |
|
|
# You can shorten the execution time by reducing variable tend from 60 to 0.5 |
23 |
|
|
|
24 |
|
|
# Importing all the necessary modules required. |
25 |
|
|
from esys.escript import * |
26 |
|
|
from esys.finley import Rectangle |
27 |
|
|
import sys |
28 |
|
|
import os |
29 |
|
|
from cblib import * |
30 |
|
|
|
31 |
|
|
# Establish a save path. |
32 |
|
|
savepath = "data/wavesolver2d002nwtest" |
33 |
|
|
# Creating a directory automatically to store the output data. |
34 |
|
|
if not os.path.isdir("data"): |
35 |
|
|
os.mkdir("data") |
36 |
|
|
if not os.path.isdir(savepath): |
37 |
|
|
os.mkdir(savepath) |
38 |
|
|
|
39 |
|
|
|
40 |
|
|
#Geometric and material property related variables. |
41 |
|
|
mx = 1000 # model lenght |
42 |
|
|
my = 200 # model width |
43 |
|
|
ndx = 200 # steps in x direction |
44 |
|
|
ndy = 40 # steps in y direction |
45 |
|
|
lam=3.462e9 #lames constant |
46 |
|
|
mu=3.462e9 #bulk modulus |
47 |
|
|
rho=1154. #density |
48 |
|
|
# Time related variables. |
49 |
|
|
tend=0.1 #end time |
50 |
|
|
#calculating )the timestep |
51 |
|
|
h=(1./5.)*sqrt(rho/(lam+2*mu))*(mx/ndx) |
52 |
|
|
#Check to make sure number of time steps is not too large. |
53 |
|
|
print "Time step size= ",h, "Expected number of outputs= ",tend/h |
54 |
|
|
proceeder = raw_input("Is this ok?(y/n)") |
55 |
|
|
#Exit if user thinks too many outputs. |
56 |
|
|
if proceeder == "n": |
57 |
|
|
sys.exit() |
58 |
|
|
|
59 |
|
|
U0=0.01 # amplitude of point source |
60 |
|
|
# spherical source at middle of bottom face |
61 |
|
|
|
62 |
|
|
xc=[300,200] |
63 |
|
|
|
64 |
|
|
mydomain=Rectangle(l0=mx,l1=my,n0=ndx, n1=ndy) |
65 |
|
|
wavesolver2df(mydomain,h,tend,lam,mu,rho,U0,xc,savepath) |
66 |
|
|
|