diff --git a/gprMax/exceptions.py b/gprMax/exceptions.py
index e6165032..e4d84d58 100644
--- a/gprMax/exceptions.py
+++ b/gprMax/exceptions.py
@@ -16,10 +16,12 @@
# You should have received a copy of the GNU General Public License
# along with gprMax. If not, see .
+
class GeneralError(ValueError):
"""Handles general errors. Subclasses the ValueError class."""
pass
+
class CmdInputError(ValueError):
"""Handles errors in user specified commands. Subclasses the ValueError class."""
pass
diff --git a/gprMax/fractals.py b/gprMax/fractals.py
index 9f36ce5e..bf9a5bcb 100644
--- a/gprMax/fractals.py
+++ b/gprMax/fractals.py
@@ -76,7 +76,7 @@ class FractalSurface(object):
self.fractalsurface = np.zeros(surfacedims, dtype=complextype)
# Positional vector at centre of array, scaled by weighting
- v1 = np.array([self.weighting[0]*(surfacedims[0])/2, self.weighting[1]*(surfacedims[1])/2])
+ v1 = np.array([self.weighting[0] * (surfacedims[0]) / 2, self.weighting[1] * (surfacedims[1]) / 2])
# 2D array of random numbers to be convolved with the fractal function
R = np.random.RandomState(self.seed)
@@ -88,13 +88,13 @@ class FractalSurface(object):
for i in range(surfacedims[0]):
for j in range(surfacedims[1]):
# Positional vector for current position
- v2 = np.array([self.weighting[0]*i, self.weighting[1]*j])
+ v2 = np.array([self.weighting[0] * i, self.weighting[1] * j])
rr = np.linalg.norm(v2 - v1)
try:
- self.fractalsurface[i, j] = A[i, j] * 1/(rr**self.b)
+ self.fractalsurface[i, j] = A[i, j] * 1 / (rr**self.b)
except FloatingPointError:
rr = 0.9
- self.fractalsurface[i, j] = A[i, j] * 1/(rr**self.b)
+ self.fractalsurface[i, j] = A[i, j] * 1 / (rr**self.b)
# Shift the zero frequency component to the centre of the spectrum
self.fractalsurface = np.fft.ifftshift(self.fractalsurface)
@@ -104,7 +104,7 @@ class FractalSurface(object):
fractalmin = np.amin(self.fractalsurface)
fractalmax = np.amax(self.fractalsurface)
fractalrange = fractalmax - fractalmin
- self.fractalsurface = self.fractalsurface * ((self.fractalrange[1] - self.fractalrange[0])/fractalrange) + self.fractalrange[0] - ((self.fractalrange[1] - self.fractalrange[0])/fractalrange) * fractalmin
+ self.fractalsurface = self.fractalsurface * ((self.fractalrange[1] - self.fractalrange[0]) / fractalrange) + self.fractalrange[0] - ((self.fractalrange[1] - self.fractalrange[0]) / fractalrange) * fractalmin
class FractalVolume(object):
@@ -118,7 +118,7 @@ class FractalVolume(object):
"""
self.ID = None
- self.operatingonID= None
+ self.operatingonID = None
self.xs = xs
self.xf = xf
self.ys = ys
@@ -146,7 +146,7 @@ class FractalVolume(object):
self.fractalvolume = np.zeros((self.nx + 1, self.ny + 1, self.nz + 1), dtype=complextype)
# Positional vector at centre of array, scaled by weighting
- v1 = np.array([self.weighting[0]*(self.nx + 1)/2, self.weighting[1]*(self.ny + 1)/2, self.weighting[2]*(self.nz + 1)/2])
+ v1 = np.array([self.weighting[0] * (self.nx + 1) / 2, self.weighting[1] * (self.ny + 1) / 2, self.weighting[2] * (self.nz + 1) / 2])
# 3D array of random numbers to be convolved with the fractal function
R = np.random.RandomState(self.seed)
@@ -159,13 +159,13 @@ class FractalVolume(object):
for j in range(self.ny + 1):
for k in range(self.nz + 1):
# Positional vector for current position
- v2 = np.array([self.weighting[0]*i, self.weighting[1]*j, self.weighting[2]*k])
+ v2 = np.array([self.weighting[0] * i, self.weighting[1] * j, self.weighting[2] * k])
rr = np.linalg.norm(v2 - v1)
try:
- self.fractalvolume[i, j, k] = A[i, j, k] * 1/(rr**self.b)
+ self.fractalvolume[i, j, k] = A[i, j, k] * 1 / (rr**self.b)
except FloatingPointError:
rr = 0.9
- self.fractalvolume[i, j, k] = A[i, j, k] * 1/(rr**self.b)
+ self.fractalvolume[i, j, k] = A[i, j, k] * 1 / (rr**self.b)
# Shift the zero frequency component to the centre of the spectrum
self.fractalvolume = np.fft.ifftshift(self.fractalvolume)
diff --git a/gprMax/geometry_views.py b/gprMax/geometry_views.py
index 3b137ff2..4c06e4d3 100644
--- a/gprMax/geometry_views.py
+++ b/gprMax/geometry_views.py
@@ -153,7 +153,7 @@ class GeometryView(object):
vtk_numpoints = (self.nx + 1) * (self.ny + 1) * (self.nz + 1)
vtk_numpoint_components = 3
vtk_numlines = 2 * self.nx * self.ny + 2 * self.ny * self.nz + 2 * self.nx * self.nz + 3 * self.nx * self.ny * self.nz + self.nx + self.ny + self.nz
- vtk_numline_components = 2
+ vtk_numline_components = 2
vtk_connectivity_offset = round_value((vtk_numpoints * vtk_numpoint_components * np.dtype(np.float32).itemsize) + np.dtype(np.uint32).itemsize)
vtk_offsets_offset = round_value(vtk_connectivity_offset + (vtk_numlines * vtk_numline_components * np.dtype(np.uint32).itemsize) + np.dtype(np.uint32).itemsize)
@@ -244,7 +244,6 @@ class GeometryView(object):
self.write_gprmax_info(f, G, materialsonly=True)
-
def write_gprmax_info(self, f, G, materialsonly=False):
"""Writes gprMax specific information relating material, source, and receiver names to numeric identifiers.
diff --git a/gprMax/grid.py b/gprMax/grid.py
index 670ae354..a9172c38 100644
--- a/gprMax/grid.py
+++ b/gprMax/grid.py
@@ -63,6 +63,7 @@ class Grid(object):
co = round_value(float(val) / getattr(self, 'd' + coord))
return co
+
class FDTDGrid(Grid):
"""Holds attributes associated with the entire grid. A convenient way for accessing regularly used parameters."""
@@ -155,7 +156,7 @@ def dispersion_check(G):
if waveform.type == 'sine' or waveform.type == 'contsine':
maxfreqs.append(4 * waveform.freq)
- elif waveform.type =='impulse':
+ elif waveform.type == 'impulse':
pass
else:
@@ -239,6 +240,7 @@ def Ix(x, y, z, Hy, Hz, G):
Ix = G.dy * (Hy[x, y, z - 1] - Hy[x, y, z]) + G.dz * (Hz[x, y, z] - Hz[x, y - 1, z])
return Ix
+
def Iy(x, y, z, Hx, Hz, G):
"""Calculates the y-component of current at a grid position.
@@ -256,6 +258,7 @@ def Iy(x, y, z, Hx, Hz, G):
Iy = G.dx * (Hx[x, y, z] - Hx[x, y, z - 1]) + G.dz * (Hz[x - 1, y, z] - Hz[x, y, z])
return Iy
+
def Iz(x, y, z, Hx, Hy, G):
"""Calculates the z-component of current at a grid position.
diff --git a/gprMax/input_cmd_funcs.py b/gprMax/input_cmd_funcs.py
index 019e67b7..c4f51f8b 100644
--- a/gprMax/input_cmd_funcs.py
+++ b/gprMax/input_cmd_funcs.py
@@ -35,6 +35,7 @@ Coordinate(x=0.1, y=0.2, z=0.3)
from collections import namedtuple
Coordinate_tuple = namedtuple('Coordinate', ['x', 'y', 'z'])
+
class Coordinate(Coordinate_tuple):
"""Subclass of a namedtuple where __str__ outputs 'x y z'"""
def __str__(self):
@@ -62,7 +63,7 @@ def command(cmd, *parameters):
except TypeError as e:
# append info about cmd and parameters to the exception:
if not e.args:
- e.args=('', )
+ e.args = ('', )
additional_info = "Creating cmd = #{} with parameters {} -> {} failed".format(cmd, parameters, filtered_list)
e.args = e.args + (additional_info,)
raise e
diff --git a/gprMax/input_cmds_geometry.py b/gprMax/input_cmds_geometry.py
index 2cbe5247..be7ebec3 100644
--- a/gprMax/input_cmds_geometry.py
+++ b/gprMax/input_cmds_geometry.py
@@ -46,9 +46,9 @@ def process_geometrycmds(geometry, G):
if len(tmp) != 6:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires exactly five parameters')
- xs = round_value(float(tmp[1])/G.dx)
- ys = round_value(float(tmp[2])/G.dy)
- zs = round_value(float(tmp[3])/G.dz)
+ xs = round_value(float(tmp[1]) / G.dx)
+ ys = round_value(float(tmp[2]) / G.dy)
+ zs = round_value(float(tmp[3]) / G.dz)
geofile = tmp[4]
matfile = tmp[5]
@@ -90,17 +90,16 @@ def process_geometrycmds(geometry, G):
if G.messages:
print('Geometry objects from file {} inserted at {:g}m, {:g}m, {:g}m, with corresponding materials file {}.'.format(geofile, xs * G.dx, ys * G.dy, zs * G.dz, matfile))
-
elif tmp[0] == '#edge:':
if len(tmp) != 8:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires exactly seven parameters')
- xs = round_value(float(tmp[1])/G.dx)
- xf = round_value(float(tmp[4])/G.dx)
- ys = round_value(float(tmp[2])/G.dy)
- yf = round_value(float(tmp[5])/G.dy)
- zs = round_value(float(tmp[3])/G.dz)
- zf = round_value(float(tmp[6])/G.dz)
+ xs = round_value(float(tmp[1]) / G.dx)
+ xf = round_value(float(tmp[4]) / G.dx)
+ ys = round_value(float(tmp[2]) / G.dy)
+ yf = round_value(float(tmp[5]) / G.dy)
+ zs = round_value(float(tmp[3]) / G.dz)
+ zf = round_value(float(tmp[6]) / G.dz)
if xs < 0 or xs > G.nx:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' the lower x-coordinate {:g}m is not within the model domain'.format(xs * G.dx))
@@ -150,7 +149,6 @@ def process_geometrycmds(geometry, G):
if G.messages:
print('Edge from {:g}m, {:g}m, {:g}m, to {:g}m, {:g}m, {:g}m of material {} created.'.format(xs * G.dx, ys * G.dy, zs * G.dz, xf * G.dx, yf * G.dy, zf * G.dz, tmp[7]))
-
elif tmp[0] == '#plate:':
if len(tmp) < 8:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires at least seven parameters')
@@ -166,12 +164,12 @@ def process_geometrycmds(geometry, G):
else:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' too many parameters have been given')
- xs = round_value(float(tmp[1])/G.dx)
- xf = round_value(float(tmp[4])/G.dx)
- ys = round_value(float(tmp[2])/G.dy)
- yf = round_value(float(tmp[5])/G.dy)
- zs = round_value(float(tmp[3])/G.dz)
- zf = round_value(float(tmp[6])/G.dz)
+ xs = round_value(float(tmp[1]) / G.dx)
+ xf = round_value(float(tmp[4]) / G.dx)
+ ys = round_value(float(tmp[2]) / G.dy)
+ yf = round_value(float(tmp[5]) / G.dy)
+ zs = round_value(float(tmp[3]) / G.dz)
+ zf = round_value(float(tmp[6]) / G.dz)
if xs < 0 or xs > G.nx:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' the lower x-coordinate {:g}m is not within the model domain'.format(xs * G.dx))
@@ -259,7 +257,6 @@ def process_geometrycmds(geometry, G):
if G.messages:
print('Plate from {:g}m, {:g}m, {:g}m, to {:g}m, {:g}m, {:g}m of material(s) {} created.'.format(xs * G.dx, ys * G.dy, zs * G.dz, xf * G.dx, yf * G.dy, zf * G.dz, ', '.join(materialsrequested)))
-
elif tmp[0] == '#triangle:':
if len(tmp) < 12:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires at least eleven parameters')
@@ -286,15 +283,15 @@ def process_geometrycmds(geometry, G):
else:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' too many parameters have been given')
- x1 = round_value(float(tmp[1])/G.dx) * G.dx
- y1 = round_value(float(tmp[2])/G.dy) * G.dy
- z1 = round_value(float(tmp[3])/G.dz) * G.dz
- x2 = round_value(float(tmp[4])/G.dx) * G.dx
- y2 = round_value(float(tmp[5])/G.dy) * G.dy
- z2 = round_value(float(tmp[6])/G.dz) * G.dz
- x3 = round_value(float(tmp[7])/G.dx) * G.dx
- y3 = round_value(float(tmp[8])/G.dy) * G.dy
- z3 = round_value(float(tmp[9])/G.dz) * G.dz
+ x1 = round_value(float(tmp[1]) / G.dx) * G.dx
+ y1 = round_value(float(tmp[2]) / G.dy) * G.dy
+ z1 = round_value(float(tmp[3]) / G.dz) * G.dz
+ x2 = round_value(float(tmp[4]) / G.dx) * G.dx
+ y2 = round_value(float(tmp[5]) / G.dy) * G.dy
+ z2 = round_value(float(tmp[6]) / G.dz) * G.dz
+ x3 = round_value(float(tmp[7]) / G.dx) * G.dx
+ y3 = round_value(float(tmp[8]) / G.dy) * G.dy
+ z3 = round_value(float(tmp[9]) / G.dz) * G.dz
thickness = float(tmp[10])
if x1 < 0 or x2 < 0 or x3 < 0 or x1 > G.nx or x2 > G.nx or x3 > G.nx:
@@ -379,7 +376,6 @@ def process_geometrycmds(geometry, G):
else:
print('Triangle with coordinates {:g}m {:g}m {:g}m, {:g}m {:g}m {:g}m, {:g}m {:g}m {:g}m of material(s) {} created.'.format(x1, y1, z1, x2, y2, z2, x3, y3, z3, ', '.join(materialsrequested)))
-
elif tmp[0] == '#box:':
if len(tmp) < 8:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires at least seven parameters')
@@ -406,12 +402,12 @@ def process_geometrycmds(geometry, G):
else:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' too many parameters have been given')
- xs = round_value(float(tmp[1])/G.dx)
- xf = round_value(float(tmp[4])/G.dx)
- ys = round_value(float(tmp[2])/G.dy)
- yf = round_value(float(tmp[5])/G.dy)
- zs = round_value(float(tmp[3])/G.dz)
- zf = round_value(float(tmp[6])/G.dz)
+ xs = round_value(float(tmp[1]) / G.dx)
+ xf = round_value(float(tmp[4]) / G.dx)
+ ys = round_value(float(tmp[2]) / G.dy)
+ yf = round_value(float(tmp[5]) / G.dy)
+ zs = round_value(float(tmp[3]) / G.dz)
+ zf = round_value(float(tmp[6]) / G.dz)
if xs < 0 or xs > G.nx:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' the lower x-coordinate {:g}m is not within the model domain'.format(xs * G.dx))
@@ -471,7 +467,6 @@ def process_geometrycmds(geometry, G):
dielectricsmoothing = 'off'
print('Box from {:g}m, {:g}m, {:g}m, to {:g}m, {:g}m, {:g}m of material(s) {} created, dielectric smoothing is {}.'.format(xs * G.dx, ys * G.dy, zs * G.dz, xf * G.dx, yf * G.dy, zf * G.dz, ', '.join(materialsrequested), dielectricsmoothing))
-
elif tmp[0] == '#cylinder:':
if len(tmp) < 9:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires at least eight parameters')
@@ -498,12 +493,12 @@ def process_geometrycmds(geometry, G):
else:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' too many parameters have been given')
- x1 = round_value(float(tmp[1])/G.dx) * G.dx
- y1 = round_value(float(tmp[2])/G.dy) * G.dy
- z1 = round_value(float(tmp[3])/G.dz) * G.dz
- x2 = round_value(float(tmp[4])/G.dx) * G.dx
- y2 = round_value(float(tmp[5])/G.dy) * G.dy
- z2 = round_value(float(tmp[6])/G.dz) * G.dz
+ x1 = round_value(float(tmp[1]) / G.dx) * G.dx
+ y1 = round_value(float(tmp[2]) / G.dy) * G.dy
+ z1 = round_value(float(tmp[3]) / G.dz) * G.dz
+ x2 = round_value(float(tmp[4]) / G.dx) * G.dx
+ y2 = round_value(float(tmp[5]) / G.dy) * G.dy
+ z2 = round_value(float(tmp[6]) / G.dz) * G.dz
r = float(tmp[7])
if r <= 0:
@@ -552,7 +547,6 @@ def process_geometrycmds(geometry, G):
dielectricsmoothing = 'off'
print('Cylinder with face centres {:g}m, {:g}m, {:g}m and {:g}m, {:g}m, {:g}m, with radius {:g}m, of material(s) {} created, dielectric smoothing is {}.'.format(x1, y1, z1, x2, y2, z2, r, ', '.join(materialsrequested), dielectricsmoothing))
-
elif tmp[0] == '#cylindrical_sector:':
if len(tmp) < 10:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires at least nine parameters')
@@ -586,8 +580,8 @@ def process_geometrycmds(geometry, G):
extent2 = float(tmp[5])
thickness = extent2 - extent1
r = float(tmp[6])
- sectorstartangle = 2*np.pi*(float(tmp[7])/360)
- sectorangle = 2*np.pi*(float(tmp[8])/360)
+ sectorstartangle = 2 * np.pi * (float(tmp[7]) / 360)
+ sectorangle = 2 * np.pi * (float(tmp[8]) / 360)
if normal != 'x' and normal != 'y' and normal != 'z':
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' the normal direction must be either x, y or z.')
@@ -648,21 +642,21 @@ def process_geometrycmds(geometry, G):
# yz-plane cylindrical sector
if normal == 'x':
- ctr1 = round_value(ctr1/G.dy) * G.dy
- ctr2 = round_value(ctr2/G.dz) * G.dz
- level = round_value(extent1/G.dx)
+ ctr1 = round_value(ctr1 / G.dy) * G.dy
+ ctr2 = round_value(ctr2 / G.dz) * G.dz
+ level = round_value(extent1 / G.dx)
# xz-plane cylindrical sector
elif normal == 'y':
- ctr1 = round_value(ctr1/G.dx) * G.dx
- ctr2 = round_value(ctr2/G.dz) * G.dz
- level = round_value(extent1/G.dy)
+ ctr1 = round_value(ctr1 / G.dx) * G.dx
+ ctr2 = round_value(ctr2 / G.dz) * G.dz
+ level = round_value(extent1 / G.dy)
# xy-plane cylindrical sector
elif normal == 'z':
- ctr1 = round_value(ctr1/G.dx) * G.dx
- ctr2 = round_value(ctr2/G.dy) * G.dy
- level = round_value(extent1/G.dz)
+ ctr1 = round_value(ctr1 / G.dx) * G.dx
+ ctr2 = round_value(ctr2 / G.dy) * G.dy
+ level = round_value(extent1 / G.dz)
build_cylindrical_sector(ctr1, ctr2, level, sectorstartangle, sectorangle, r, normal, thickness, G.dx, G.dy, G.dz, numID, numIDx, numIDy, numIDz, averaging, G.solid, G.rigidE, G.rigidH, G.ID)
@@ -672,10 +666,9 @@ def process_geometrycmds(geometry, G):
dielectricsmoothing = 'on'
else:
dielectricsmoothing = 'off'
- print('Cylindrical sector with centre {:g}m, {:g}m, radius {:g}m, starting angle {:.1f} degrees, sector angle {:.1f} degrees, thickness {:g}m, of material(s) {} created, dielectric smoothing is {}.'.format(ctr1, ctr2, r, (sectorstartangle/(2*np.pi))*360, (sectorangle/(2*np.pi))*360, thickness, ', '.join(materialsrequested), dielectricsmoothing))
+ print('Cylindrical sector with centre {:g}m, {:g}m, radius {:g}m, starting angle {:.1f} degrees, sector angle {:.1f} degrees, thickness {:g}m, of material(s) {} created, dielectric smoothing is {}.'.format(ctr1, ctr2, r, (sectorstartangle / (2 * np.pi)) * 360, (sectorangle / (2 * np.pi)) * 360, thickness, ', '.join(materialsrequested), dielectricsmoothing))
else:
- print('Cylindrical sector with centre {:g}m, {:g}m, radius {:g}m, starting angle {:.1f} degrees, sector angle {:.1f} degrees, of material(s) {} created.'.format(ctr1, ctr2, r, (sectorstartangle/(2*np.pi))*360, (sectorangle/(2*np.pi))*360, ', '.join(materialsrequested)))
-
+ print('Cylindrical sector with centre {:g}m, {:g}m, radius {:g}m, starting angle {:.1f} degrees, sector angle {:.1f} degrees, of material(s) {} created.'.format(ctr1, ctr2, r, (sectorstartangle / (2 * np.pi)) * 360, (sectorangle / (2 * np.pi)) * 360, ', '.join(materialsrequested)))
elif tmp[0] == '#sphere:':
if len(tmp) < 6:
@@ -704,9 +697,9 @@ def process_geometrycmds(geometry, G):
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' too many parameters have been given')
# Centre of sphere
- xc = round_value(float(tmp[1])/G.dx)
- yc = round_value(float(tmp[2])/G.dy)
- zc = round_value(float(tmp[3])/G.dz)
+ xc = round_value(float(tmp[1]) / G.dx)
+ yc = round_value(float(tmp[2]) / G.dy)
+ zc = round_value(float(tmp[3]) / G.dz)
r = float(tmp[4])
# Look up requested materials in existing list of material instances
@@ -752,7 +745,6 @@ def process_geometrycmds(geometry, G):
dielectricsmoothing = 'off'
print('Sphere with centre {:g}m, {:g}m, {:g}m, radius {:g}m, of material(s) {} created, dielectric smoothing is {}.'.format(xc * G.dx, yc * G.dy, zc * G.dz, r, ', '.join(materialsrequested), dielectricsmoothing))
-
elif tmp[0] == '#fractal_box:':
if len(tmp) < 14:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires at least thirteen parameters')
@@ -763,12 +755,12 @@ def process_geometrycmds(geometry, G):
else:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' too many parameters have been given')
- xs = round_value(float(tmp[1])/G.dx)
- xf = round_value(float(tmp[4])/G.dx)
- ys = round_value(float(tmp[2])/G.dy)
- yf = round_value(float(tmp[5])/G.dy)
- zs = round_value(float(tmp[3])/G.dz)
- zf = round_value(float(tmp[6])/G.dz)
+ xs = round_value(float(tmp[1]) / G.dx)
+ xf = round_value(float(tmp[4]) / G.dx)
+ ys = round_value(float(tmp[2]) / G.dy)
+ yf = round_value(float(tmp[5]) / G.dy)
+ zs = round_value(float(tmp[3]) / G.dz)
+ zf = round_value(float(tmp[6]) / G.dz)
if xs < 0 or xs > G.nx:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' the lower x-coordinate {:g}m is not within the model domain'.format(xs * G.dx))
@@ -836,12 +828,12 @@ def process_geometrycmds(geometry, G):
# Only process rough surfaces for this fractal volume
if tmp[12] == volume.ID:
- xs = round_value(float(tmp[1])/G.dx)
- xf = round_value(float(tmp[4])/G.dx)
- ys = round_value(float(tmp[2])/G.dy)
- yf = round_value(float(tmp[5])/G.dy)
- zs = round_value(float(tmp[3])/G.dz)
- zf = round_value(float(tmp[6])/G.dz)
+ xs = round_value(float(tmp[1]) / G.dx)
+ xf = round_value(float(tmp[4]) / G.dx)
+ ys = round_value(float(tmp[2]) / G.dy)
+ yf = round_value(float(tmp[5]) / G.dy)
+ zs = round_value(float(tmp[3]) / G.dz)
+ zf = round_value(float(tmp[6]) / G.dz)
if xs < 0 or xs > G.nx:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' the lower x-coordinate {:g}m is not within the model domain'.format(xs * G.dx))
@@ -943,12 +935,12 @@ def process_geometrycmds(geometry, G):
# Only process surfaces for this fractal volume
if tmp[8] == volume.ID:
- xs = round_value(float(tmp[1])/G.dx)
- xf = round_value(float(tmp[4])/G.dx)
- ys = round_value(float(tmp[2])/G.dy)
- yf = round_value(float(tmp[5])/G.dy)
- zs = round_value(float(tmp[3])/G.dz)
- zf = round_value(float(tmp[6])/G.dz)
+ xs = round_value(float(tmp[1]) / G.dx)
+ xf = round_value(float(tmp[4]) / G.dx)
+ ys = round_value(float(tmp[2]) / G.dy)
+ yf = round_value(float(tmp[5]) / G.dy)
+ zs = round_value(float(tmp[3]) / G.dz)
+ zf = round_value(float(tmp[6]) / G.dz)
depth = float(tmp[7])
if xs < 0 or xs > G.nx:
@@ -1044,7 +1036,6 @@ def process_geometrycmds(geometry, G):
if G.messages:
print('Water on surface from {:g}m, {:g}m, {:g}m, to {:g}m, {:g}m, {:g}m with depth {:g}m, added to {}.'.format(xs * G.dx, ys * G.dy, zs * G.dz, xf * G.dx, yf * G.dy, zf * G.dz, filldepth, surface.operatingonID))
-
if tmp[0] == '#add_grass:':
if len(tmp) < 12:
raise CmdInputError("'" + ' '.join(tmp) + "'" + ' requires at least eleven parameters')
@@ -1057,12 +1048,12 @@ def process_geometrycmds(geometry, G):
# Only process grass for this fractal volume
if tmp[11] == volume.ID:
- xs = round_value(float(tmp[1])/G.dx)
- xf = round_value(float(tmp[4])/G.dx)
- ys = round_value(float(tmp[2])/G.dy)
- yf = round_value(float(tmp[5])/G.dy)
- zs = round_value(float(tmp[3])/G.dz)
- zf = round_value(float(tmp[6])/G.dz)
+ xs = round_value(float(tmp[1]) / G.dx)
+ xf = round_value(float(tmp[4]) / G.dx)
+ ys = round_value(float(tmp[2]) / G.dy)
+ yf = round_value(float(tmp[5]) / G.dy)
+ zs = round_value(float(tmp[3]) / G.dz)
+ zf = round_value(float(tmp[6]) / G.dz)
numblades = int(tmp[10])
if xs < 0 or xs > G.nx:
@@ -1194,7 +1185,6 @@ def process_geometrycmds(geometry, G):
if G.messages:
print('{} blades of grass on surface from {:g}m, {:g}m, {:g}m, to {:g}m, {:g}m, {:g}m with fractal dimension {:g}, fractal seeding {}, and range {:g}m to {:g}m, added to {}.'.format(numblades, xs * G.dx, ys * G.dy, zs * G.dz, xf * G.dx, yf * G.dy, zf * G.dz, surface.dimension, surface.seed, float(tmp[8]), float(tmp[9]), surface.operatingonID))
-
# Process any modifications to the original fractal box then generate it
if volume.fractalsurfaces:
volume.originalxs = volume.xs
@@ -1439,7 +1429,6 @@ def process_geometrycmds(geometry, G):
depth += 1
k -= 1
-
# Build voxels from any true values of the 3D mask array
waternumID = next((x.numID for x in G.materials if x.ID == 'water'), 0)
grassnumID = next((x.numID for x in G.materials if x.ID == 'grass'), 0)
diff --git a/gprMax/input_cmds_multiuse.py b/gprMax/input_cmds_multiuse.py
index 20495e33..7feb5288 100644
--- a/gprMax/input_cmds_multiuse.py
+++ b/gprMax/input_cmds_multiuse.py
@@ -71,7 +71,6 @@ def process_multicmds(multicmds, G):
G.waveforms.append(w)
-
# Voltage source
cmdname = '#voltage_source'
if multicmds[cmdname] != 'None':
@@ -98,7 +97,7 @@ def process_multicmds(multicmds, G):
raise CmdInputError("'" + cmdname + ': ' + ' '.join(tmp) + "'" + ' there is no waveform with the identifier {}'.format(tmp[5]))
v = VoltageSource()
- v.polarisation= tmp[0]
+ v.polarisation = tmp[0]
v.xcoord = xcoord
v.ycoord = ycoord
v.zcoord = zcoord
@@ -132,7 +131,6 @@ def process_multicmds(multicmds, G):
G.voltagesources.append(v)
-
# Hertzian dipole
cmdname = '#hertzian_dipole'
if multicmds[cmdname] != 'None':
@@ -192,7 +190,6 @@ def process_multicmds(multicmds, G):
G.hertziandipoles.append(h)
-
# Magnetic dipole
cmdname = '#magnetic_dipole'
if multicmds[cmdname] != 'None':
@@ -325,9 +322,9 @@ def process_multicmds(multicmds, G):
raise CmdInputError("'" + cmdname + ': ' + ' '.join(tmp) + "'" + ' has an incorrect number of parameters')
# Check position parameters
- xcoord = round_value(float(tmp[0])/G.dx)
- ycoord = round_value(float(tmp[1])/G.dy)
- zcoord = round_value(float(tmp[2])/G.dz)
+ xcoord = round_value(float(tmp[0]) / G.dx)
+ ycoord = round_value(float(tmp[1]) / G.dy)
+ zcoord = round_value(float(tmp[2]) / G.dz)
check_coordinates(xcoord, ycoord, zcoord)
if xcoord < G.pmlthickness[0] or xcoord > G.nx - G.pmlthickness[3] or ycoord < G.pmlthickness[1] or ycoord > G.ny - G.pmlthickness[4] or zcoord < G.pmlthickness[2] or zcoord > G.nz - G.pmlthickness[5]:
print("WARNING: '" + cmdname + ': ' + ' '.join(tmp) + "'" + ' sources and receivers should not normally be positioned within the PML.\n')
@@ -359,7 +356,6 @@ def process_multicmds(multicmds, G):
G.rxs.append(r)
-
# Receiver array
cmdname = '#rx_array'
if multicmds[cmdname] != 'None':
@@ -425,7 +421,6 @@ def process_multicmds(multicmds, G):
print(' Receiver at {:g}m, {:g}m, {:g}m with output(s) {} created.'.format(r.xcoord * G.dx, r.ycoord * G.dy, r.zcoord * G.dz, ', '.join(r.outputs)))
G.rxs.append(r)
-
# Snapshot
cmdname = '#snapshot'
if multicmds[cmdname] != 'None':
@@ -511,7 +506,6 @@ def process_multicmds(multicmds, G):
# Append the new material object to the materials list
G.materials.append(m)
-
cmdname = '#add_dispersion_debye'
if multicmds[cmdname] != 'None':
for cmdinstance in multicmds[cmdname]:
@@ -583,7 +577,6 @@ def process_multicmds(multicmds, G):
if G.messages:
print('Lorentz disperion added to {} with delta_epsr={}, omega={} secs, and gamma={} created.'.format(material.ID, ', '.join('%4.2f' % deltaer for deltaer in material.deltaer), ', '.join('%4.3e' % tau for tau in material.tau), ', '.join('%4.3e' % alpha for alpha in material.alpha)))
-
cmdname = '#add_dispersion_drude'
if multicmds[cmdname] != 'None':
for cmdinstance in multicmds[cmdname]:
@@ -619,7 +612,6 @@ def process_multicmds(multicmds, G):
if G.messages:
print('Drude disperion added to {} with omega={} secs, and gamma={} secs created.'.format(material.ID, ', '.join('%4.3e' % tau for tau in material.tau), ', '.join('%4.3e' % alpha for alpha in material.alpha)))
-
cmdname = '#soil_peplinski'
if multicmds[cmdname] != 'None':
for cmdinstance in multicmds[cmdname]:
@@ -650,7 +642,6 @@ def process_multicmds(multicmds, G):
# Append the new material object to the materials list
G.mixingmodels.append(s)
-
# Geometry views (creates VTK-based geometry files)
cmdname = '#geometry_view'
if multicmds[cmdname] != 'None':
@@ -695,7 +686,6 @@ def process_multicmds(multicmds, G):
# Append the new GeometryView object to the geometry views list
G.geometryviews.append(g)
-
# Complex frequency shifted (CFS) PML parameter
cmdname = '#pml_cfs'
if multicmds[cmdname] != 'None':
diff --git a/gprMax/input_cmds_singleuse.py b/gprMax/input_cmds_singleuse.py
index 57a5f058..9466a8e1 100644
--- a/gprMax/input_cmds_singleuse.py
+++ b/gprMax/input_cmds_singleuse.py
@@ -50,7 +50,6 @@ def process_singlecmds(singlecmds, G):
else:
raise CmdInputError(cmd + ' requires input values of either y or n')
-
# Title
cmd = '#title'
if singlecmds[cmd] != 'None':
@@ -58,7 +57,6 @@ def process_singlecmds(singlecmds, G):
if G.messages:
print('Model title: {}'.format(G.title))
-
# Number of processors to run on (OpenMP)
cmd = '#num_threads'
os.environ['OMP_WAIT_POLICY'] = 'active'
@@ -85,7 +83,6 @@ def process_singlecmds(singlecmds, G):
if G.nthreads > psutil.cpu_count(logical=False):
print('\nWARNING: You have specified more threads ({}) than available physical CPU cores ({}). This may lead to degraded performance.'.format(G.nthreads, psutil.cpu_count(logical=False)))
-
# Spatial discretisation
cmd = '#dx_dy_dz'
tmp = [float(x) for x in singlecmds[cmd].split()]
@@ -103,15 +100,14 @@ def process_singlecmds(singlecmds, G):
if G.messages:
print('Spatial discretisation: {:g} x {:g} x {:g}m'.format(G.dx, G.dy, G.dz))
-
# Domain
cmd = '#domain'
tmp = [float(x) for x in singlecmds[cmd].split()]
if len(tmp) != 3:
raise CmdInputError(cmd + ' requires exactly three parameters')
- G.nx = round_value(tmp[0]/G.dx)
- G.ny = round_value(tmp[1]/G.dy)
- G.nz = round_value(tmp[2]/G.dz)
+ G.nx = round_value(tmp[0] / G.dx)
+ G.ny = round_value(tmp[1] / G.dy)
+ G.nz = round_value(tmp[2] / G.dz)
if G.nx == 0 or G.ny == 0 or G.nz == 0:
raise CmdInputError(cmd + ' requires at least one cell in every dimension')
if G.messages:
@@ -120,7 +116,6 @@ def process_singlecmds(singlecmds, G):
mem = (((G.nx + 1) * (G.ny + 1) * (G.nz + 1) * 13 * np.dtype(floattype).itemsize + (G.nx + 1) * (G.ny + 1) * (G.nz + 1) * 18) * 1.1) + 30e6
print('Memory (RAM) usage: ~{} required, {} available'.format(human_size(mem), human_size(psutil.virtual_memory().total)))
-
# Time step CFL limit (use either 2D or 3D) and default PML thickness
if G.nx == 1:
G.dt = 1 / (c * np.sqrt((1 / G.dy) * (1 / G.dy) + (1 / G.dz) * (1 / G.dz)))
@@ -145,7 +140,6 @@ def process_singlecmds(singlecmds, G):
if G.messages:
print('Time step (at {} CFL limit): {:g} secs'.format(G.dtlimit, G.dt))
-
# Time step stability factor
cmd = '#time_step_stability_factor'
if singlecmds[cmd] != 'None':
@@ -158,7 +152,6 @@ def process_singlecmds(singlecmds, G):
if G.messages:
print('Time step (modified): {:g} secs'.format(G.dt))
-
# Time window
cmd = '#time_window'
tmp = singlecmds[cmd].split()
@@ -182,7 +175,6 @@ def process_singlecmds(singlecmds, G):
if G.messages:
print('Time window: {:g} secs ({} iterations)'.format(G.timewindow, G.iterations))
-
# PML
cmd = '#pml_cells'
if singlecmds[cmd] != 'None':
@@ -193,36 +185,33 @@ def process_singlecmds(singlecmds, G):
G.pmlthickness = (int(tmp[0]), int(tmp[0]), int(tmp[0]), int(tmp[0]), int(tmp[0]), int(tmp[0]))
else:
G.pmlthickness = (int(tmp[0]), int(tmp[1]), int(tmp[2]), int(tmp[3]), int(tmp[4]), int(tmp[5]))
- if 2*G.pmlthickness[0] >= G.nx or 2*G.pmlthickness[1] >= G.ny or 2*G.pmlthickness[2] >= G.nz or 2*G.pmlthickness[3] >= G.nx or 2*G.pmlthickness[4] >= G.ny or 2*G.pmlthickness[5] >= G.nz:
+ if 2 * G.pmlthickness[0] >= G.nx or 2 * G.pmlthickness[1] >= G.ny or 2 * G.pmlthickness[2] >= G.nz or 2 * G.pmlthickness[3] >= G.nx or 2 * G.pmlthickness[4] >= G.ny or 2 * G.pmlthickness[5] >= G.nz:
raise CmdInputError(cmd + ' has too many cells for the domain size')
-
# src_steps
cmd = '#src_steps'
if singlecmds[cmd] != 'None':
tmp = singlecmds[cmd].split()
if len(tmp) != 3:
raise CmdInputError(cmd + ' requires exactly three parameters')
- G.srcstepx = round_value(float(tmp[0])/G.dx)
- G.srcstepy = round_value(float(tmp[1])/G.dy)
- G.srcstepz = round_value(float(tmp[2])/G.dz)
+ G.srcstepx = round_value(float(tmp[0]) / G.dx)
+ G.srcstepy = round_value(float(tmp[1]) / G.dy)
+ G.srcstepz = round_value(float(tmp[2]) / G.dz)
if G.messages:
print('Simple sources will step {:g}m, {:g}m, {:g}m for each model run.'.format(G.srcstepx * G.dx, G.srcstepy * G.dy, G.srcstepz * G.dz))
-
# rx_steps
cmd = '#rx_steps'
if singlecmds[cmd] != 'None':
tmp = singlecmds[cmd].split()
if len(tmp) != 3:
raise CmdInputError(cmd + ' requires exactly three parameters')
- G.rxstepx = round_value(float(tmp[0])/G.dx)
- G.rxstepy = round_value(float(tmp[1])/G.dy)
- G.rxstepz = round_value(float(tmp[2])/G.dz)
+ G.rxstepx = round_value(float(tmp[0]) / G.dx)
+ G.rxstepy = round_value(float(tmp[1]) / G.dy)
+ G.rxstepz = round_value(float(tmp[2]) / G.dz)
if G.messages:
print('All receivers will step {:g}m, {:g}m, {:g}m for each model run.'.format(G.rxstepx * G.dx, G.rxstepy * G.dy, G.rxstepz * G.dz))
-
# Excitation file for user-defined source waveforms
cmd = '#excitation_file'
if singlecmds[cmd] != 'None':
@@ -251,7 +240,7 @@ def process_singlecmds(singlecmds, G):
if len(waveformvalues.shape) == 1:
w.uservalues = waveformvalues[:]
else:
- w.uservalues = waveformvalues[:,waveform]
+ w.uservalues = waveformvalues[:, waveform]
if G.messages:
print('User waveform {} created.'.format(w.ID))
diff --git a/gprMax/materials.py b/gprMax/materials.py
index fa6efc68..45cf6bfa 100644
--- a/gprMax/materials.py
+++ b/gprMax/materials.py
@@ -74,8 +74,8 @@ class Material(object):
G (class): Grid class instance - holds essential parameters describing the model.
"""
- HA = (m0*self.mr / G.dt) + 0.5*self.sm
- HB = (m0*self.mr / G.dt) - 0.5*self.sm
+ HA = (m0 * self.mr / G.dt) + 0.5 * self.sm
+ HB = (m0 * self.mr / G.dt) - 0.5 * self.sm
self.DA = HB / HA
self.DBx = (1 / G.dx) * 1 / HA
self.DBy = (1 / G.dy) * 1 / HA
@@ -122,12 +122,12 @@ class Material(object):
self.zt[x] = (self.w[x] / self.q[x]) * (1 - self.eqt[x]) / G.dt
self.zt2[x] = (self.w[x] / self.q[x]) * (1 - self.eqt2[x])
- EA = (e0*self.er / G.dt) + 0.5*self.se - (e0 / G.dt) * np.sum(self.zt2.real)
- EB = (e0*self.er / G.dt) - 0.5*self.se - (e0 / G.dt) * np.sum(self.zt2.real)
+ EA = (e0 * self.er / G.dt) + 0.5 * self.se - (e0 / G.dt) * np.sum(self.zt2.real)
+ EB = (e0 * self.er / G.dt) - 0.5 * self.se - (e0 / G.dt) * np.sum(self.zt2.real)
else:
- EA = (e0*self.er / G.dt) + 0.5*self.se
- EB = (e0*self.er / G.dt) - 0.5*self.se
+ EA = (e0 * self.er / G.dt) + 0.5 * self.se
+ EB = (e0 * self.er / G.dt) - 0.5 * self.se
if self.ID == 'pec' or self.se == float('inf'):
self.CA = 0
@@ -198,15 +198,15 @@ class PeplinskiSoil(object):
muiter = np.nditer(mumaterials, flags=['c_index'])
while not muiter.finished:
# Real part for frequencies in the range 1.4GHz to 18GHz
- er1 = (1 + (self.rb/self.rs) * ((es**a) - 1) + (muiter[0]**b1 * erealw**a) - muiter[0]) ** (1/a)
+ er1 = (1 + (self.rb / self.rs) * ((es**a) - 1) + (muiter[0]**b1 * erealw**a) - muiter[0]) ** (1 / a)
# Real part for frequencies in the range 0.3GHz to 1.3GHz
er2 = 1.15 * er1 - 0.68
# Imaginary part for frequencies in the range 0.3GHz to 1.3GHz
- eri = er2 - (muiter[0]**(b2/a) * Material.waterdeltaer)
+ eri = er2 - (muiter[0]**(b2 / a) * Material.waterdeltaer)
# Effective conductivity
- sig = muiter[0]**(b2/a) * ((sigf1 * (self.rs - self.rb)) / (self.rs * muiter[0]))
+ sig = muiter[0]**(b2 / a) * ((sigf1 * (self.rs - self.rb)) / (self.rs * muiter[0]))
# Check to see if the material already exists before creating a new one
requiredID = '|{:.4f}|'.format(float(muiter[0]))
diff --git a/gprMax/optimisation_taguchi.py b/gprMax/optimisation_taguchi.py
index 783f2e1e..102eb7fc 100644
--- a/gprMax/optimisation_taguchi.py
+++ b/gprMax/optimisation_taguchi.py
@@ -68,7 +68,7 @@ def run_opt_sim(args, numbermodelruns, inputfile, usernamespace):
# Select OA
OA, N, cols, k, s, t = construct_OA(optparams)
- print('\n{}\nTaguchi optimisation...\n'.format(68*'*'))
+ print('\n{}\nTaguchi optimisation...\n'.format(68 * '*'))
print('\tOrthogonal array: {:g} experiments per iteration, {:g} parameters ({:g} will be used), {:g} levels, and strength {:g}'.format(N, cols, k, s, t))
tmp = [(k, v) for k, v in optparams.items()]
print('\tParameters to optimise with ranges: {}'.format(str(tmp).strip('[]')))
@@ -131,7 +131,7 @@ def run_opt_sim(args, numbermodelruns, inputfile, usernamespace):
# Rename confirmation experiment output file so that it is retained for each iteraction
os.rename(outputfile, os.path.splitext(outputfile)[0] + '_final' + str(iteration + 1) + '.out')
- print('\nTaguchi optimisation, iteration {} completed. History of optimal parameter values {} and of fitness values {}'.format(iteration + 1, dict(optparamshist), fitnessvalueshist, 68*'*'))
+ print('\nTaguchi optimisation, iteration {} completed. History of optimal parameter values {} and of fitness values {}'.format(iteration + 1, dict(optparamshist), fitnessvalueshist, 68 * '*'))
iteration += 1
# Stop optimisation if stopping criterion has been reached
@@ -154,7 +154,7 @@ def run_opt_sim(args, numbermodelruns, inputfile, usernamespace):
pickle.dump(fitnessvalueshist, f)
pickle.dump(optparamsinit, f)
- print('\n{}\nTaguchi optimisation completed after {} iteration(s).\nHistory of optimal parameter values {} and of fitness values {}\n{}\n'.format(68*'*', iteration, dict(optparamshist), fitnessvalueshist, 68*'*'))
+ print('\n{}\nTaguchi optimisation completed after {} iteration(s).\nHistory of optimal parameter values {} and of fitness values {}\n{}\n'.format(68 * '*', iteration, dict(optparamshist), fitnessvalueshist, 68 * '*'))
def taguchi_code_blocks(inputfile, taguchinamespace):
@@ -313,7 +313,7 @@ def calculate_ranges_experiments(optparams, optparamsinit, levels, levelsopt, le
# Gaussian reduction function used for calculating levels
T = 18 # Usually values between 15 - 20
- RR = np.exp(-(i/T)**2)
+ RR = np.exp(-(i / T)**2)
# Calculate levels for each parameter
for p in range(k):
@@ -341,7 +341,6 @@ def calculate_ranges_experiments(optparams, optparamsinit, levels, levelsopt, le
levels[0, p] = levels[1, p] - levelsdiff[p]
levels[2, p] = levels[1, p] + levelsdiff[p]
-
# Update dictionary of parameters to optimise with lists of new values; clear dictionary first
optparams = OrderedDict((key, list()) for key in optparams)
p = 0
diff --git a/gprMax/pml.py b/gprMax/pml.py
index 52ea5b36..4c94e476 100644
--- a/gprMax/pml.py
+++ b/gprMax/pml.py
@@ -98,7 +98,7 @@ class CFS(object):
Hvalues (float): numpy array holding scaling profile values for magnetic PML update.
"""
- tmp = (np.linspace(0, (len(Evalues) - 1) + 0.5, num=2*len(Evalues)) / (len(Evalues) - 1)) ** order
+ tmp = (np.linspace(0, (len(Evalues) - 1) + 0.5, num=2 * len(Evalues)) / (len(Evalues) - 1)) ** order
Evalues = tmp[0:-1:2]
Hvalues = tmp[1::2]
return Evalues, Hvalues
@@ -115,7 +115,7 @@ class CFS(object):
Hvalues (float): numpy array holding profile value for magnetic PML update.
"""
- Evalues= np.zeros(thickness + 1, dtype=floattype)
+ Evalues = np.zeros(thickness + 1, dtype=floattype)
Hvalues = np.zeros(thickness + 1, dtype=floattype)
if parameter.scalingprofile == 'constant':
@@ -218,18 +218,18 @@ class PML(object):
Esigma, Hsigma = cfs.calculate_values(self.thickness, cfs.sigma)
# Electric PML update coefficients
- tmp = (2*e0*Ekappa) + G.dt * (Ealpha * Ekappa + Esigma)
- self.ERA[x, :] = (2*e0 + G.dt*Ealpha) / tmp
- self.ERB[x, :] = (2*e0*Ekappa) / tmp
- self.ERE[x, :] = ((2*e0*Ekappa) - G.dt * (Ealpha * Ekappa + Esigma)) / tmp
- self.ERF[x, :] = (2*Esigma*G.dt) / (Ekappa * tmp)
+ tmp = (2 * e0 * Ekappa) + G.dt * (Ealpha * Ekappa + Esigma)
+ self.ERA[x, :] = (2 * e0 + G.dt * Ealpha) / tmp
+ self.ERB[x, :] = (2 * e0 * Ekappa) / tmp
+ self.ERE[x, :] = ((2 * e0 * Ekappa) - G.dt * (Ealpha * Ekappa + Esigma)) / tmp
+ self.ERF[x, :] = (2 * Esigma * G.dt) / (Ekappa * tmp)
# Magnetic PML update coefficients
- tmp = (2*e0*Hkappa) + G.dt * (Halpha * Hkappa + Hsigma)
- self.HRA[x, :] = (2*e0 + G.dt*Halpha) / tmp
- self.HRB[x, :] = (2*e0*Hkappa) / tmp
- self.HRE[x, :] = ((2*e0*Hkappa) - G.dt * (Halpha * Hkappa + Hsigma)) / tmp
- self.HRF[x, :] = (2*Hsigma*G.dt) / (Hkappa * tmp)
+ tmp = (2 * e0 * Hkappa) + G.dt * (Halpha * Hkappa + Hsigma)
+ self.HRA[x, :] = (2 * e0 + G.dt * Halpha) / tmp
+ self.HRB[x, :] = (2 * e0 * Hkappa) / tmp
+ self.HRE[x, :] = ((2 * e0 * Hkappa) - G.dt * (Halpha * Hkappa + Hsigma)) / tmp
+ self.HRF[x, :] = (2 * Hsigma * G.dt) / (Hkappa * tmp)
def build_pmls(G):
@@ -244,7 +244,7 @@ def build_pmls(G):
if pmlthickness > 0:
sumer = 0 # Sum of relative permittivities in PML slab
summr = 0 # Sum of relative permeabilities in PML slab
- pmldirection=PML.directions[index]
+ pmldirection = PML.directions[index]
if pmldirection[0] == 'x':
if pmldirection == 'xminus':
diff --git a/gprMax/sources.py b/gprMax/sources.py
index 66b3e63d..7e3fabbe 100644
--- a/gprMax/sources.py
+++ b/gprMax/sources.py
@@ -185,13 +185,13 @@ class MagneticDipole(Source):
waveform = next(x for x in G.waveforms if x.ID == self.waveformID)
if self.polarisation is 'x':
- Hx[i, j, k] -= waveform.amp * waveform.calculate_value(time, G.dt) * (G.dt / (G.dx * G.dy * G.dz))
+ Hx[i, j, k] -= waveform.amp * waveform.calculate_value(time, G.dt) * (G.dt / (G.dx * G.dy * G.dz))
elif self.polarisation is 'y':
- Hy[i, j, k] -= waveform.amp * waveform.calculate_value(time, G.dt) * (G.dt / (G.dx * G.dy * G.dz))
+ Hy[i, j, k] -= waveform.amp * waveform.calculate_value(time, G.dt) * (G.dt / (G.dx * G.dy * G.dz))
elif self.polarisation is 'z':
- Hz[i, j, k] -= waveform.amp * waveform.calculate_value(time, G.dt) * (G.dt / (G.dx * G.dy * G.dz))
+ Hz[i, j, k] -= waveform.amp * waveform.calculate_value(time, G.dt) * (G.dt / (G.dx * G.dy * G.dz))
class TransmissionLine(Source):
diff --git a/gprMax/utilities.py b/gprMax/utilities.py
index 5ba3f364..d47e431d 100644
--- a/gprMax/utilities.py
+++ b/gprMax/utilities.py
@@ -57,9 +57,9 @@ along with gprMax. If not, see ."""
url = 'www.gprmax.com'
print('\nElectromagnetic modelling software based on the Finite-Difference \nTime-Domain (FDTD) method')
- print('\n{} {} {}'.format('*'*round((width - len(url))/2), url, '*'*round((width - len(url))/2)))
+ print('\n{} {} {}'.format('*' * round((width - len(url)) / 2), url, '*' * round((width - len(url)) / 2)))
print('{}'.format(gprMaxlogo))
- print('{} v{} {}'.format('*'*round((width - len(version))/2), (version), '*'*round((width - len(version))/2)))
+ print('{} v{} {}'.format('*' * round((width - len(version)) / 2), (version), '*' * round((width - len(version)) / 2)))
print(licenseinfo)
@@ -91,12 +91,12 @@ def round_value(value, decimalplaces=0):
# Rounds to nearest integer (half values are rounded downwards)
if decimalplaces == 0:
- rounded = int(d.Decimal(value).quantize(d.Decimal('1'),rounding=d.ROUND_HALF_DOWN))
+ rounded = int(d.Decimal(value).quantize(d.Decimal('1'), rounding=d.ROUND_HALF_DOWN))
# Rounds down to nearest float represented by number of decimal places
else:
precision = '1.{places}'.format(places='0' * decimalplaces)
- rounded = float(d.Decimal(value).quantize(d.Decimal(precision),rounding=d.ROUND_FLOOR))
+ rounded = float(d.Decimal(value).quantize(d.Decimal(precision), rounding=d.ROUND_FLOOR))
return rounded
diff --git a/gprMax/writer_hdf5.py b/gprMax/writer_hdf5.py
index 2ce2ffa0..26206c22 100644
--- a/gprMax/writer_hdf5.py
+++ b/gprMax/writer_hdf5.py
@@ -43,7 +43,6 @@ def write_hdf5(outputfile, Ex, Ey, Ez, Hx, Hy, Hz, G):
f.attrs['srcsteps'] = (G.srcstepx, G.srcstepy, G.srcstepz)
f.attrs['rxsteps'] = (G.rxstepx, G.rxstepy, G.rxstepz)
-
# Create group for sources (except transmission lines); add type and positional data attributes
srclist = G.voltagesources + G.hertziandipoles + G.magneticdipoles
for srcindex, src in enumerate(srclist):
@@ -65,7 +64,6 @@ def write_hdf5(outputfile, Ex, Ey, Ez, Hx, Hy, Hz, G):
f['/tls/tl' + str(tlindex + 1) + '/Vtotal'] = tl.Vtotal
f['/tls/tl' + str(tlindex + 1) + '/Itotal'] = tl.Itotal
-
# Create group, add positional data and write field component arrays for receivers
for rxindex, rx in enumerate(G.rxs):
grp = f.create_group('/rxs/rx' + str(rxindex + 1))
diff --git a/gprMax/xdmf.py b/gprMax/xdmf.py
index 92841401..07ef71b2 100644
--- a/gprMax/xdmf.py
+++ b/gprMax/xdmf.py
@@ -6,6 +6,7 @@ import numpy as np
from gprMax.grid import Grid
+
class ListCounter():
def __init__(self, npArray):
@@ -16,10 +17,10 @@ class ListCounter():
self.array[self.count] = item
self.count += 1
+
class EdgeLabels:
def __init__(self, grid):
-
"""
Class to define some connectivity of for an n x l x m
grid
@@ -37,6 +38,7 @@ class EdgeLabels:
edge = np.array([in_label, out_label])
self.edge_counter.add(edge)
+
class EdgeMaterials:
def __init__(self, fdtd_grid):
@@ -85,7 +87,6 @@ class SolidLabels():
self.label_counter = ListCounter(self.solid_labels)
def hexCellPicker(self, grid, i, j, k):
-
"""
This is the ordering of nodes in the hexahedron cell.
@@ -124,6 +125,7 @@ class SolidLabels():
solid_labels = self.hexCellPicker(self.label_grid.grid, i, j, k)
self.label_counter.add(solid_labels)
+
class SolidManager(Grid):
def __init__(self, label_grid, fdtd_grid):