Brendan Smithyman | January 2015
In [1]:
import numpy as np
import networkx
from zephyr.Problem import SeisFDFDProblem
import time
In [2]:
import matplotlib.pyplot as plt
import matplotlib.cm as cm
import matplotlib
%matplotlib inline
In [3]:
cellSize = 1 # m
freqs = [2e2] # Hz
density = 2700 # units of density
Q = np.inf # can be inf
nx = 164 # count
nz = 264 # count
freeSurf = [False, False, False, False] # t r b l
dims = (nx,nz) # tuple
nPML = 32
rho = np.fliplr(np.ones(dims) * density)
nfreq = len(freqs) # number of frequencies
nky = 48 # number of y-directional plane-wave components
nsp = nfreq * nky # total number of 2D subproblems
velocity = 2500 # m/s
vanom = 500 # m/s
cPert = np.zeros(dims)
cPert[(nx/2)-20:(nx/2)+20,(nz/2)-20:(nz/2)+20] = vanom
c = np.fliplr(np.ones(dims) * velocity)
cFlat = c
c += np.fliplr(cPert)
cTrue = c
srcs = np.array([np.ones(101)*32, np.zeros(101), np.linspace(32, 232, 101)]).T
recs = np.array([np.ones(101)*132, np.zeros(101), np.linspace(32, 232, 101)]).T
nsrc = len(srcs)
nrec = len(recs)
recmode = 'fixed'
geom = {
'src': srcs,
'rec': recs,
'mode': 'fixed',
}
cache = False
cacheDir = '.'
parFac = 1
profile = 'phobos'
# Base configuration for all subproblems
systemConfig = {
'dx': cellSize, # m
'dz': cellSize, # m
'c': c.T, # m/s
'rho': rho.T, # density
'Q': Q, # can be inf
'nx': nx, # count
'nz': nz, # count
'freeSurf': freeSurf, # t r b l
'nPML': nPML,
'geom': geom,
'cache': cache,
'cacheDir': cacheDir,
'freqs': freqs,
'nky': nky,
'parFac': parFac,
'profile': profile,
}
In [4]:
sp = SeisFDFDProblem(systemConfig)
In [5]:
tick = time.time()
jobs, G = sp.forwardAccumulate()
sp.par['lview'].wait([j for i in jobs['endJobs'].values() for j in i])
tock = time.time()
print('Time elapsed: %r'%(tock - tick))
In [6]:
parFac = 2
systemConfig.update({'parFac': parFac})
In [7]:
tick = time.time()
jobs, G = sp.forwardAccumulate()
sp.par['lview'].wait([j for i in jobs['endJobs'].values() for j in i])
tock = time.time()
print('Time elapsed: %r'%(tock - tick))
In [8]:
parFac = 3
systemConfig.update({'parFac': parFac})
In [9]:
tick = time.time()
jobs, G = sp.forwardAccumulate()
sp.par['lview'].wait([j for i in jobs['endJobs'].values() for j in i])
tock = time.time()
print('Time elapsed: %r'%(tock - tick))
In [ ]: