publicImage/gprMax
1
0
派生 0
镜像自地址 https://gitee.com/sunhf/gprMax.git 已同步 2025-08-07 23:14:03 +08:00
代码 工单 软件包 项目 版本发布 百科 动态

Merge branch 'devel' of https://github.com/gprMax/gprMax into devel

浏览代码
这个提交包含在：
Antonis Giannopoulos
2023-04-21 17:04:04 +01:00
父节点 577aa6445a 19bb9ee474
当前提交 5df218601b
共有 56 个文件被更改，包括 773 次插入 和 823 次删除
下载 Patch 文件 下载 Diff 文件 展开所有文件 折叠所有文件

13
.github/FUNDING.yml vendored 普通文件
查看文件

@@ -0,0 +1,13 @@
# These are supported funding model platforms
github: # Replace with up to 4 GitHub Sponsors-enabled usernames e.g., [user1, user2]
patreon: # Replace with a single Patreon username
open_collective: gprMax
ko_fi: # Replace with a single Ko-fi username
tidelift: # Replace with a single Tidelift platform-name/package-name e.g., npm/babel
community_bridge: # Replace with a single Community Bridge project-name e.g., cloud-foundry
liberapay: # Replace with a single Liberapay username
issuehunt: # Replace with a single IssueHunt username
otechie: # Replace with a single Otechie username
lfx_crowdfunding: # Replace with a single LFX Crowdfunding project-name e.g., cloud-foundry
custom: # Replace with up to 4 custom sponsorship URLs e.g., ['link1', 'link2']

二进制
examples/antenna_like_GSSI_1500.vti
查看文件

二进制文件未显示。

二进制
examples/antenna_like_GSSI_1500_fs.h5
查看文件

二进制文件未显示。

2
examples/antenna_like_GSSI_1500_fs.py
查看文件

@@ -41,7 +41,7 @@ gv2 = gprMax.GeometryView(p1=(ant_pos[0] - 0.170/2, ant_pos[1] - 0.108/2, ant_po
p2=(ant_pos[0] + 0.170/2, ant_pos[1] + 0.108/2, ant_pos[2] + 0.010), p2=(ant_pos[0] + 0.170/2, ant_pos[1] + 0.108/2, ant_pos[2] + 0.010),
dl=(dl, dl, dl), filename='antenna_like_GSSI_1500_pcb', dl=(dl, dl, dl), filename='antenna_like_GSSI_1500_pcb',
output_type='f') output_type='f')
scene.add(gv1) #scene.add(gv1)
#scene.add(gv2) #scene.add(gv2)
# Run model # Run model

二进制
examples/heterogeneous_soil_v2.h5
查看文件

二进制文件未显示。

13
examples/subgrids/cylinder_fs.py
查看文件

@@ -89,17 +89,10 @@ gv1 = gprMax.GeometryView(p1=(0, 0, 0), p2=(x, y, z), dl=(dl, dl, dl),
scene.add(gv1) scene.add(gv1)
# Create some snapshots of entire domain # Create some snapshots of entire domain
# for i in range(5):
# s = gprMax.Snapshot(p1=(0, 0, 0), p2=(x, y, z), dl=(dl, dl, dl),
# time=(i + 0.5) * 1e-9,
# filename=fn.with_suffix('').parts[-1] + '_' + str(i + 1))
# scene.add(s)
for i in range(5): for i in range(5):
s = gprMax.Snapshot(p1=sg1, p2=sg2, dl=(dl_sg, dl_sg, dl_sg), s = gprMax.Snapshot(p1=(0, 0, 0), p2=(x, y, z), dl=(dl, dl, dl),
time=(i + 0.5) * 1e-9, time=(i + 0.5) * 1e-9,
filename=fn.with_suffix('').parts[-1] + '_' + str(i + 1), filename=fn.with_suffix('').parts[-1] + '_' + str(i + 1))
outputs=['Ez']) scene.add(s)
subgrid.add(s)
gprMax.run(scenes=[scene], n=1, geometry_only=False, outputfile=fn, subgrid=True, autotranslate=True) gprMax.run(scenes=[scene], n=1, geometry_only=False, outputfile=fn, subgrid=True, autotranslate=True)

8
examples/subgrids/gssi_400_over_fractal_subsurface.py
查看文件

@@ -102,9 +102,11 @@ halfspace = gprMax.Box(p1=(0, 0, 0), p2=(x, y, antenna_p[2]), material_id='soil'
scene.add(halfspace) scene.add(halfspace)
for i in range(1, 51): for i in range(1, 51):
snap = gprMax.Snapshot(p1=(0, y / 2, 0), p2=(x, y / 2 + dl, z), dl=(dl, dl, dl), snap = gprMax.Snapshot(p1=(0, y / 2, 0),
filename=Path(*parts[:-1], parts[-1] + '_' + str(i)).name, p2=(x, y / 2 + dl, z),
time=i * tw / 50) dl=(dl, dl, dl),
filename=Path(*parts[:-1], f'{parts[-1]}_{str(i)}').name,
time=i * tw / 50,)
scene.add(snap) scene.add(snap)
# create a geometry view of the main grid and the sub grid stitched together # create a geometry view of the main grid and the sub grid stitched together

5
gprMax/__init__.py
查看文件

@@ -25,8 +25,9 @@ from .cmds_geometry.triangle import Triangle
from .cmds_multiuse import (PMLCFS, AddDebyeDispersion, AddDrudeDispersion, from .cmds_multiuse import (PMLCFS, AddDebyeDispersion, AddDrudeDispersion,
AddLorentzDispersion, GeometryObjectsWrite, AddLorentzDispersion, GeometryObjectsWrite,
GeometryView, HertzianDipole, MagneticDipole, GeometryView, HertzianDipole, MagneticDipole,
Material, Rx, RxArray, Snapshot, SoilPeplinski, Material, MaterialList, MaterialRange, Rx, RxArray,
TransmissionLine, VoltageSource, Waveform, MaterialRange, MaterialList) Snapshot, SoilPeplinski, TransmissionLine,
VoltageSource, Waveform)
from .cmds_singleuse import (Discretisation, Domain, ExcitationFile, from .cmds_singleuse import (Discretisation, Domain, ExcitationFile,
OMPThreads, PMLProps, RxSteps, SrcSteps, OMPThreads, PMLProps, RxSteps, SrcSteps,
TimeStepStabilityFactor, TimeWindow, Title) TimeStepStabilityFactor, TimeWindow, Title)

40
gprMax/cmds_geometry/add_grass.py
查看文件

@@ -74,7 +74,7 @@ class AddGrass(UserObjectGeometry):
limits = self.kwargs['limits'] limits = self.kwargs['limits']
n_blades = self.kwargs['n_blades'] n_blades = self.kwargs['n_blades']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires at least eleven parameters') logger.exception(f'{self.__str__()} requires at least eleven parameters')
raise raise
try: try:
@@ -90,7 +90,7 @@ class AddGrass(UserObjectGeometry):
try: try:
volume = volumes[0] volume = volumes[0]
except NameError: except NameError:
logger.exception(self.__str__() + f' cannot find FractalBox {fractal_box_id}') logger.exception(f'{self.__str__()} cannot find FractalBox {fractal_box_id}')
raise raise
p1, p2 = uip.check_box_points(p1, p2, self.__str__()) p1, p2 = uip.check_box_points(p1, p2, self.__str__())
@@ -98,32 +98,32 @@ class AddGrass(UserObjectGeometry):
xf, yf, zf = p2 xf, yf, zf = p2
if frac_dim < 0: if frac_dim < 0:
logger.exception(self.__str__() + ' requires a positive value for ' + logger.exception(f'{self.__str__()} requires a positive value for ' +
'the fractal dimension') 'the fractal dimension')
raise ValueError raise ValueError
if limits[0] < 0 or limits[1] < 0: if limits[0] < 0 or limits[1] < 0:
logger.exception(self.__str__() + ' requires a positive value for ' + logger.exception(f'{self.__str__()} requires a positive value for ' +
'the minimum and maximum heights for grass blades') 'the minimum and maximum heights for grass blades')
raise ValueError raise ValueError
# Check for valid orientations # Check for valid orientations
if xs == xf: if xs == xf:
if ys == yf or zs == zf: if ys == yf or zs == zf:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
if xs != volume.xs and xs != volume.xf: if xs != volume.xs and xs != volume.xf:
logger.exception(self.__str__() + ' must specify external surfaces on a fractal box') logger.exception(f'{self.__str__()} must specify external surfaces on a fractal box')
raise ValueError raise ValueError
fractalrange = (round_value(limits[0] / grid.dx), round_value(limits[1] / grid.dx)) fractalrange = (round_value(limits[0] / grid.dx), round_value(limits[1] / grid.dx))
# xminus surface # xminus surface
if xs == volume.xs: if xs == volume.xs:
logger.exception(self.__str__() + ' grass can only be specified ' + logger.exception(f'{self.__str__()} grass can only be specified ' +
'on surfaces in the positive axis direction') 'on surfaces in the positive axis direction')
raise ValueError raise ValueError
# xplus surface # xplus surface
elif xf == volume.xf: elif xf == volume.xf:
if fractalrange[1] > grid.nx: if fractalrange[1] > grid.nx:
logger.exception(self.__str__() + ' cannot apply grass to ' + logger.exception(f'{self.__str__()} cannot apply grass to ' +
'fractal box as it would exceed the domain ' + 'fractal box as it would exceed the domain ' +
'size in the x direction') 'size in the x direction')
raise ValueError raise ValueError
@@ -131,21 +131,21 @@ class AddGrass(UserObjectGeometry):
elif ys == yf: elif ys == yf:
if xs == xf or zs == zf: if xs == xf or zs == zf:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
if ys != volume.ys and ys != volume.yf: if ys != volume.ys and ys != volume.yf:
logger.exception(self.__str__() + ' must specify external surfaces on a fractal box') logger.exception(f'{self.__str__()} must specify external surfaces on a fractal box')
raise ValueError raise ValueError
fractalrange = (round_value(limits[0] / grid.dy), round_value(limits[1] / grid.dy)) fractalrange = (round_value(limits[0] / grid.dy), round_value(limits[1] / grid.dy))
# yminus surface # yminus surface
if ys == volume.ys: if ys == volume.ys:
logger.exception(self.__str__() + ' grass can only be specified ' + logger.exception(f'{self.__str__()} grass can only be specified ' +
'on surfaces in the positive axis direction') 'on surfaces in the positive axis direction')
raise ValueError raise ValueError
# yplus surface # yplus surface
elif yf == volume.yf: elif yf == volume.yf:
if fractalrange[1] > grid.ny: if fractalrange[1] > grid.ny:
logger.exception(self.__str__() + ' cannot apply grass to ' + logger.exception(f'{self.__str__()} cannot apply grass to ' +
'fractal box as it would exceed the domain ' + 'fractal box as it would exceed the domain ' +
'size in the y direction') 'size in the y direction')
raise ValueError raise ValueError
@@ -153,28 +153,28 @@ class AddGrass(UserObjectGeometry):
elif zs == zf: elif zs == zf:
if xs == xf or ys == yf: if xs == xf or ys == yf:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
if zs != volume.zs and zs != volume.zf: if zs != volume.zs and zs != volume.zf:
logger.exception(self.__str__() + ' must specify external surfaces on a fractal box') logger.exception(f'{self.__str__()} must specify external surfaces on a fractal box')
raise ValueError raise ValueError
fractalrange = (round_value(limits[0] / grid.dz), round_value(limits[1] / grid.dz)) fractalrange = (round_value(limits[0] / grid.dz), round_value(limits[1] / grid.dz))
# zminus surface # zminus surface
if zs == volume.zs: if zs == volume.zs:
logger.exception(self.__str__() + ' grass can only be specified ' + logger.exception(f'{self.__str__()} grass can only be specified ' +
'on surfaces in the positive axis direction') 'on surfaces in the positive axis direction')
raise ValueError raise ValueError
# zplus surface # zplus surface
elif zf == volume.zf: elif zf == volume.zf:
if fractalrange[1] > grid.nz: if fractalrange[1] > grid.nz:
logger.exception(self.__str__() + ' cannot apply grass to ' + logger.exception(f'{self.__str__()} cannot apply grass to ' +
'fractal box as it would exceed the domain ' + 'fractal box as it would exceed the domain ' +
'size in the z direction') 'size in the z direction')
raise ValueError raise ValueError
requestedsurface = 'zplus' requestedsurface = 'zplus'
else: else:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
surface = FractalSurface(xs, xf, ys, yf, zs, zf, frac_dim) surface = FractalSurface(xs, xf, ys, yf, zs, zf, frac_dim)
@@ -187,7 +187,7 @@ class AddGrass(UserObjectGeometry):
surface.operatingonID = volume.ID surface.operatingonID = volume.ID
surface.generate_fractal_surface() surface.generate_fractal_surface()
if n_blades > surface.fractalsurface.shape[0] * surface.fractalsurface.shape[1]: if n_blades > surface.fractalsurface.shape[0] * surface.fractalsurface.shape[1]:
logger.exception(self.__str__() + ' the specified surface is not large ' + logger.exception(f'{self.__str__()} the specified surface is not large ' +
'enough for the number of grass blades/roots specified') 'enough for the number of grass blades/roots specified')
raise ValueError raise ValueError
@@ -234,13 +234,13 @@ class AddGrass(UserObjectGeometry):
grass = next((x for x in grid.materials if x.ID == 'grass')) grass = next((x for x in grid.materials if x.ID == 'grass'))
testgrass = next((x for x in grass.tau if x < grid.dt), None) testgrass = next((x for x in grass.tau if x < grid.dt), None)
if testgrass: if testgrass:
logger.exception(self.__str__() + ' requires the time step for the ' + logger.exception(f'{self.__str__()} requires the time step for the ' +
'model to be less than the relaxation time required to model grass.') 'model to be less than the relaxation time required to model grass.')
raise ValueError raise ValueError
volume.fractalsurfaces.append(surface) volume.fractalsurfaces.append(surface)
logger.info(self.grid_name(grid) + f'{n_blades} blades of grass on surface from ' + logger.info(f'{self.grid_name(grid)}{n_blades} blades of grass on surface from ' +
f'{xs * grid.dx:g}m, {ys * grid.dy:g}m, {zs * grid.dz:g}m, ' + f'{xs * grid.dx:g}m, {ys * grid.dy:g}m, {zs * grid.dz:g}m, ' +
f'to {xf * grid.dx:g}m, {yf * grid.dy:g}m, {zf * grid.dz:g}m ' + f'to {xf * grid.dx:g}m, {yf * grid.dy:g}m, {zf * grid.dz:g}m ' +
f'with fractal dimension {surface.dimension:g}, fractal seeding ' + f'with fractal dimension {surface.dimension:g}, fractal seeding ' +

42
gprMax/cmds_geometry/add_surface_roughness.py
查看文件

@@ -74,13 +74,13 @@ class AddSurfaceRoughness(UserObjectGeometry):
limits = np.array(self.kwargs['limits']) limits = np.array(self.kwargs['limits'])
fractal_box_id = self.kwargs['fractal_box_id'] fractal_box_id = self.kwargs['fractal_box_id']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' incorrect parameters') logger.exception(f'{self.__str__()} incorrect parameters')
raise raise
try: try:
seed = self.kwargs['seed'] seed = self.kwargs['seed']
except KeyError: except KeyError:
logger.warning(self.__str__() + ' no value for seed detected. This ' + logger.warning(f'{self.__str__()} no value for seed detected. This ' +
'means you will get a different fractal distribution ' + 'means you will get a different fractal distribution ' +
'every time the model runs.') 'every time the model runs.')
seed = None seed = None
@@ -93,7 +93,7 @@ class AddSurfaceRoughness(UserObjectGeometry):
if volumes: if volumes:
volume = volumes[0] volume = volumes[0]
else: else:
logger.exception(self.__str__() + f' cannot find FractalBox {fractal_box_id}') logger.exception(f'{self.__str__()} cannot find FractalBox {fractal_box_id}')
raise ValueError raise ValueError
p1, p2 = uip.check_box_points(p1, p2, self.__str__()) p1, p2 = uip.check_box_points(p1, p2, self.__str__())
@@ -101,25 +101,25 @@ class AddSurfaceRoughness(UserObjectGeometry):
xf, yf, zf = p2 xf, yf, zf = p2
if frac_dim < 0: if frac_dim < 0:
logger.exception(self.__str__() + ' requires a positive value for the ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'fractal dimension') 'fractal dimension')
raise ValueError raise ValueError
if weighting[0] < 0: if weighting[0] < 0:
logger.exception(self.__str__() + ' requires a positive value for the ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'fractal weighting in the first direction of the surface') 'fractal weighting in the first direction of the surface')
raise ValueError raise ValueError
if weighting[1] < 0: if weighting[1] < 0:
logger.exception(self.__str__() + ' requires a positive value for the ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'fractal weighting in the second direction of the surface') 'fractal weighting in the second direction of the surface')
raise ValueError raise ValueError
# Check for valid orientations # Check for valid orientations
if xs == xf: if xs == xf:
if ys == yf or zs == zf: if ys == yf or zs == zf:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
if xs != volume.xs and xs != volume.xf: if xs != volume.xs and xs != volume.xf:
logger.exception(self.__str__() + ' can only be used on the external ' + logger.exception(f'{self.__str__()} can only be used on the external ' +
'surfaces of a fractal box') 'surfaces of a fractal box')
raise ValueError raise ValueError
fractalrange = (round_value(limits[0] / grid.dx), fractalrange = (round_value(limits[0] / grid.dx),
@@ -127,7 +127,7 @@ class AddSurfaceRoughness(UserObjectGeometry):
# xminus surface # xminus surface
if xs == volume.xs: if xs == volume.xs:
if fractalrange[0] < 0 or fractalrange[1] > volume.xf: if fractalrange[0] < 0 or fractalrange[1] > volume.xf:
logger.exception(self.__str__() + ' cannot apply fractal surface ' + logger.exception(f'{self.__str__()} cannot apply fractal surface ' +
'to fractal box as it would exceed either the ' + 'to fractal box as it would exceed either the ' +
'upper coordinates of the fractal box or the ' + 'upper coordinates of the fractal box or the ' +
'domain in the x direction') 'domain in the x direction')
@@ -136,7 +136,7 @@ class AddSurfaceRoughness(UserObjectGeometry):
# xplus surface # xplus surface
elif xf == volume.xf: elif xf == volume.xf:
if fractalrange[0] < volume.xs or fractalrange[1] > grid.nx: if fractalrange[0] < volume.xs or fractalrange[1] > grid.nx:
logger.exception(self.__str__() + ' cannot apply fractal surface ' + logger.exception(f'{self.__str__()} cannot apply fractal surface ' +
'to fractal box as it would exceed either the ' + 'to fractal box as it would exceed either the ' +
'lower coordinates of the fractal box or the ' + 'lower coordinates of the fractal box or the ' +
'domain in the x direction') 'domain in the x direction')
@@ -145,10 +145,10 @@ class AddSurfaceRoughness(UserObjectGeometry):
elif ys == yf: elif ys == yf:
if xs == xf or zs == zf: if xs == xf or zs == zf:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
if ys != volume.ys and ys != volume.yf: if ys != volume.ys and ys != volume.yf:
logger.exception(self.__str__() + ' can only be used on the external ' + logger.exception(f'{self.__str__()} can only be used on the external ' +
'surfaces of a fractal box') 'surfaces of a fractal box')
raise ValueError raise ValueError
fractalrange = (round_value(limits[0] / grid.dy), fractalrange = (round_value(limits[0] / grid.dy),
@@ -156,7 +156,7 @@ class AddSurfaceRoughness(UserObjectGeometry):
# yminus surface # yminus surface
if ys == volume.ys: if ys == volume.ys:
if fractalrange[0] < 0 or fractalrange[1] > volume.yf: if fractalrange[0] < 0 or fractalrange[1] > volume.yf:
logger.exception(self.__str__() + ' cannot apply fractal surface ' + logger.exception(f'{self.__str__()} cannot apply fractal surface ' +
'to fractal box as it would exceed either the ' + 'to fractal box as it would exceed either the ' +
'upper coordinates of the fractal box or the ' + 'upper coordinates of the fractal box or the ' +
'domain in the y direction') 'domain in the y direction')
@@ -165,7 +165,7 @@ class AddSurfaceRoughness(UserObjectGeometry):
# yplus surface # yplus surface
elif yf == volume.yf: elif yf == volume.yf:
if fractalrange[0] < volume.ys or fractalrange[1] > grid.ny: if fractalrange[0] < volume.ys or fractalrange[1] > grid.ny:
logger.exception(self.__str__() + ' cannot apply fractal surface ' + logger.exception(f'{self.__str__()} cannot apply fractal surface ' +
'to fractal box as it would exceed either the ' + 'to fractal box as it would exceed either the ' +
'lower coordinates of the fractal box or the ' + 'lower coordinates of the fractal box or the ' +
'domain in the y direction') 'domain in the y direction')
@@ -174,10 +174,10 @@ class AddSurfaceRoughness(UserObjectGeometry):
elif zs == zf: elif zs == zf:
if xs == xf or ys == yf: if xs == xf or ys == yf:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
if zs != volume.zs and zs != volume.zf: if zs != volume.zs and zs != volume.zf:
logger.exception(self.__str__() + ' can only be used on the external ' + logger.exception(f'{self.__str__()} can only be used on the external ' +
'surfaces of a fractal box') 'surfaces of a fractal box')
raise ValueError raise ValueError
fractalrange = (round_value(limits[0] / grid.dz), fractalrange = (round_value(limits[0] / grid.dz),
@@ -185,7 +185,7 @@ class AddSurfaceRoughness(UserObjectGeometry):
# zminus surface # zminus surface
if zs == volume.zs: if zs == volume.zs:
if fractalrange[0] < 0 or fractalrange[1] > volume.zf: if fractalrange[0] < 0 or fractalrange[1] > volume.zf:
logger.exception(self.__str__() + ' cannot apply fractal surface ' + logger.exception(f'{self.__str__()} cannot apply fractal surface ' +
'to fractal box as it would exceed either the ' + 'to fractal box as it would exceed either the ' +
'upper coordinates of the fractal box or the ' + 'upper coordinates of the fractal box or the ' +
'domain in the x direction') 'domain in the x direction')
@@ -194,7 +194,7 @@ class AddSurfaceRoughness(UserObjectGeometry):
# zplus surface # zplus surface
elif zf == volume.zf: elif zf == volume.zf:
if fractalrange[0] < volume.zs or fractalrange[1] > grid.nz: if fractalrange[0] < volume.zs or fractalrange[1] > grid.nz:
logger.exception(self.__str__() + ' cannot apply fractal surface ' + logger.exception(f'{self.__str__()} cannot apply fractal surface ' +
'to fractal box as it would exceed either the ' + 'to fractal box as it would exceed either the ' +
'lower coordinates of the fractal box or the ' + 'lower coordinates of the fractal box or the ' +
'domain in the z direction') 'domain in the z direction')
@@ -202,7 +202,7 @@ class AddSurfaceRoughness(UserObjectGeometry):
requestedsurface = 'zplus' requestedsurface = 'zplus'
else: else:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
surface = FractalSurface(xs, xf, ys, yf, zs, zf, frac_dim) surface = FractalSurface(xs, xf, ys, yf, zs, zf, frac_dim)
@@ -218,14 +218,14 @@ class AddSurfaceRoughness(UserObjectGeometry):
# List of existing surfaces IDs # List of existing surfaces IDs
existingsurfaceIDs = [x.surfaceID for x in volume.fractalsurfaces] existingsurfaceIDs = [x.surfaceID for x in volume.fractalsurfaces]
if surface.surfaceID in existingsurfaceIDs: if surface.surfaceID in existingsurfaceIDs:
logger.exception(self.__str__() + f' has already been used on the ' + logger.exception(f'{self.__str__()} has already been used on the ' +
f'{surface.surfaceID} surface') f'{surface.surfaceID} surface')
raise ValueError raise ValueError
surface.generate_fractal_surface() surface.generate_fractal_surface()
volume.fractalsurfaces.append(surface) volume.fractalsurfaces.append(surface)
logger.info(self.grid_name(grid) + f'Fractal surface from {xs * grid.dx:g}m, ' + logger.info(f'{self.grid_name(grid)}Fractal surface from {xs * grid.dx:g}m, ' +
f'{ys * grid.dy:g}m, {zs * grid.dz:g}m, to {xf * grid.dx:g}m, ' + f'{ys * grid.dy:g}m, {zs * grid.dz:g}m, to {xf * grid.dx:g}m, ' +
f'{yf * grid.dy:g}m, {zf * grid.dz:g}m with fractal dimension ' + f'{yf * grid.dy:g}m, {zf * grid.dz:g}m with fractal dimension ' +
f'{surface.dimension:g}, fractal weightings {surface.weighting[0]:g}, ' + f'{surface.dimension:g}, fractal weightings {surface.weighting[0]:g}, ' +

62
gprMax/cmds_geometry/add_surface_water.py
查看文件

@@ -68,18 +68,16 @@ class AddSurfaceWater(UserObjectGeometry):
fractal_box_id = self.kwargs['fractal_box_id'] fractal_box_id = self.kwargs['fractal_box_id']
depth = self.kwargs['depth'] depth = self.kwargs['depth']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires exactly eight parameters') logger.exception(f'{self.__str__()} requires exactly eight parameters')
raise raise
if self.do_rotate: if self.do_rotate:
self._do_rotate() self._do_rotate()
# Get the correct fractal volume if volumes := [volume for volume in grid.fractalvolumes if volume.ID == fractal_box_id]:
volumes = [volume for volume in grid.fractalvolumes if volume.ID == fractal_box_id]
if volumes:
volume = volumes[0] volume = volumes[0]
else: else:
logger.exception(self.__str__() + f' cannot find FractalBox {fractal_box_id}') logger.exception(f'{self.__str__()} cannot find FractalBox {fractal_box_id}')
raise ValueError raise ValueError
p1, p2 = uip.check_box_points(p1, p2, self.__str__()) p1, p2 = uip.check_box_points(p1, p2, self.__str__())
@@ -87,18 +85,18 @@ class AddSurfaceWater(UserObjectGeometry):
xf, yf, zf = p2 xf, yf, zf = p2
if depth <= 0: if depth <= 0:
logger.exception(self.__str__() + ' requires a positive value for ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'the depth of water') f'depth of water')
raise ValueError raise ValueError
# Check for valid orientations # Check for valid orientations
if xs == xf: if xs == xf:
if ys == yf or zs == zf: if ys == yf or zs == zf:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
if xs != volume.xs and xs != volume.xf: if xs not in [volume.xs, volume.xf]:
logger.exception(self.__str__() + ' can only be used on the ' + logger.exception(f'{self.__str__()} can only be used on the external surfaces '
'external surfaces of a fractal box') f'of a fractal box')
raise ValueError raise ValueError
# xminus surface # xminus surface
if xs == volume.xs: if xs == volume.xs:
@@ -110,12 +108,12 @@ class AddSurfaceWater(UserObjectGeometry):
filldepth = filldepthcells * grid.dx filldepth = filldepthcells * grid.dx
elif ys == yf: elif ys == yf:
if xs == xf or zs == zf: if zs == zf:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
if ys != volume.ys and ys != volume.yf: if ys not in [volume.ys, volume.yf]:
logger.exception(self.__str__() + ' can only be used on the ' + logger.exception(f'{self.__str__()} can only be used on the external surfaces ' +
'external surfaces of a fractal box') f'of a fractal box')
raise ValueError raise ValueError
# yminus surface # yminus surface
if ys == volume.ys: if ys == volume.ys:
@@ -127,12 +125,9 @@ class AddSurfaceWater(UserObjectGeometry):
filldepth = filldepthcells * grid.dy filldepth = filldepthcells * grid.dy
elif zs == zf: elif zs == zf:
if xs == xf or ys == yf: if zs not in [volume.zs, volume.zf]:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} can only be used on the external surfaces '
raise ValueError f'of a fractal box')
if zs != volume.zs and zs != volume.zf:
logger.exception(self.__str__() + ' can only be used on the ' +
'external surfaces of a fractal box')
raise ValueError raise ValueError
# zminus surface # zminus surface
if zs == volume.zs: if zs == volume.zs:
@@ -144,38 +139,35 @@ class AddSurfaceWater(UserObjectGeometry):
filldepth = filldepthcells * grid.dz filldepth = filldepthcells * grid.dz
else: else:
logger.exception(self.__str__() + ' dimensions are not specified correctly') logger.exception(f'{self.__str__()} dimensions are not specified correctly')
raise ValueError raise ValueError
surface = next((x for x in volume.fractalsurfaces if x.surfaceID == requestedsurface), None) surface = next((x for x in volume.fractalsurfaces if x.surfaceID == requestedsurface), None)
if not surface: if not surface:
logger.exception(self.__str__() + f' specified surface {requestedsurface} ' + logger.exception(f'{self.__str__()} specified surface {requestedsurface} ' +
'does not have a rough surface applied') f'does not have a rough surface applied')
raise ValueError raise ValueError
surface.filldepth = filldepthcells surface.filldepth = filldepthcells
# Check that requested fill depth falls within range of surface roughness # Check that requested fill depth falls within range of surface roughness
if surface.filldepth < surface.fractalrange[0] or surface.filldepth > surface.fractalrange[1]: if surface.filldepth < surface.fractalrange[0] or surface.filldepth > surface.fractalrange[1]:
logger.exception(self.__str__() + ' requires a value for the depth ' + logger.exception(f'{self.__str__()} requires a value for the depth of water that lies with the ' +
'of water that lies with the range of the requested ' + f'range of the requested surface roughness')
'surface roughness')
raise ValueError raise ValueError
# Check to see if water has been already defined as a material # Check to see if water has been already defined as a material
if not any(x.ID == 'water' for x in grid.materials): if all(x.ID != 'water' for x in grid.materials):
create_water(grid) create_water(grid)
# Check if time step for model is suitable for using water # Check if time step for model is suitable for using water
water = next((x for x in grid.materials if x.ID == 'water')) water = next((x for x in grid.materials if x.ID == 'water'))
testwater = next((x for x in water.tau if x < grid.dt), None) if testwater := next((x for x in water.tau if x < grid.dt), None):
if testwater: logger.exception(f'{self.__str__()} requires the time step for the model '
logger.exception(self.__str__() + ' requires the time step for the ' + f'to be less than the relaxation time required to model water.')
'model to be less than the relaxation time required ' +
'to model water.')
raise ValueError raise ValueError
logger.info(self.grid_name(grid) + f'Water on surface from {xs * grid.dx:g}m, ' + logger.info(f'{self.grid_name(grid)}Water on surface from {xs * grid.dx:g}m, ' +
f'{ys * grid.dy:g}m, {zs * grid.dz:g}m, to {xf * grid.dx:g}m, ' + f'{ys * grid.dy:g}m, {zs * grid.dz:g}m, to {xf * grid.dx:g}m, ' +
f'{yf * grid.dy:g}m, {zf * grid.dz:g}m with depth {filldepth:g}m, ' + f'{yf * grid.dy:g}m, {zf * grid.dz:g}m with depth {filldepth:g}m, ' +
f'added to {surface.operatingonID}.') f'added to {surface.operatingonID}.')

8
gprMax/cmds_geometry/box.py
查看文件

@@ -65,7 +65,7 @@ class Box(UserObjectGeometry):
p1 = self.kwargs['p1'] p1 = self.kwargs['p1']
p2 = self.kwargs['p2'] p2 = self.kwargs['p2']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' Please specify two points.') logger.exception(f'{self.__str__()} Please specify two points.')
raise raise
if self.do_rotate: if self.do_rotate:
@@ -80,7 +80,7 @@ class Box(UserObjectGeometry):
try: try:
materialsrequested = self.kwargs['material_ids'] materialsrequested = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' No materials have been specified') logger.exception(f'{self.__str__()} No materials have been specified')
raise raise
# Check averaging # Check averaging
@@ -103,7 +103,7 @@ class Box(UserObjectGeometry):
if len(materials) != len(materialsrequested): if len(materials) != len(materialsrequested):
notfound = [x for x in materialsrequested if x not in materials] notfound = [x for x in materialsrequested if x not in materials]
logger.exception(self.__str__() + f' material(s) {notfound} do not exist') logger.exception(f'{self.__str__()} material(s) {notfound} do not exist')
raise ValueError raise ValueError
# Isotropic case # Isotropic case
@@ -144,7 +144,7 @@ class Box(UserObjectGeometry):
dielectricsmoothing = 'on' if averaging else 'off' dielectricsmoothing = 'on' if averaging else 'off'
logger.info(self.grid_name(grid) + f"Box from {p5[0]:g}m, {p5[1]:g}m, " + logger.info(f"{self.grid_name(grid)}Box from {p5[0]:g}m, {p5[1]:g}m, " +
f"{p5[2]:g}m, to {p6[0]:g}m, {p6[1]:g}m, {p6[2]:g}m of " + f"{p5[2]:g}m, to {p6[0]:g}m, {p6[1]:g}m, {p6[2]:g}m of " +
f"material(s) {', '.join(materialsrequested)} created, " + f"material(s) {', '.join(materialsrequested)} created, " +
f"dielectric smoothing is {dielectricsmoothing}.") f"dielectric smoothing is {dielectricsmoothing}.")

3
gprMax/cmds_geometry/build_templates.py
查看文件

@@ -66,6 +66,5 @@ r = template.render(
}] }]
) )
f = open('cython/dispersive_updates_test.pyx', 'w') with open('cython/dispersive_updates_test.pyx', 'w') as f:
f.write(r) f.write(r)
f.close()

6
gprMax/cmds_geometry/cmds_geometry.py
查看文件

@@ -38,9 +38,9 @@ class UserObjectGeometry:
"""Readable string of parameters given to object.""" """Readable string of parameters given to object."""
s = '' s = ''
for _, v in self.kwargs.items(): for _, v in self.kwargs.items():
if isinstance(v, tuple) or isinstance(v, list): if isinstance(v, (tuple, list)):
v = ' '.join([str(el) for el in v]) v = ' '.join([str(el) for el in v])
s += str(v) + ' ' s += f'{str(v)} '
return f'{self.hash}: {s[:-1]}' return f'{self.hash}: {s[:-1]}'
@@ -120,7 +120,7 @@ def rotate_2point_object(pts, axis, angle, origin=None):
raise ValueError raise ValueError
# Check axis is valid # Check axis is valid
if axis != 'x' and axis != 'y' and axis != 'z': if axis not in ['x', 'y', 'z']:
logger.exception('Axis of rotation must be x, y, or z') logger.exception('Axis of rotation must be x, y, or z')
raise ValueError raise ValueError

14
gprMax/cmds_geometry/cone.py
查看文件

@@ -55,7 +55,7 @@ class Cone(UserObjectGeometry):
r1 = self.kwargs['r1'] r1 = self.kwargs['r1']
r2 = self.kwargs['r2'] r2 = self.kwargs['r2']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' please specify 2 points and two radii') logger.exception(f'{self.__str__()} please specify 2 points and two radii')
raise raise
# Check averaging # Check averaging
@@ -75,7 +75,7 @@ class Cone(UserObjectGeometry):
try: try:
materialsrequested = self.kwargs['material_ids'] materialsrequested = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' no materials have been specified') logger.exception(f'{self.__str__()} no materials have been specified')
raise raise
p3 = uip.round_to_grid_static_point(p1) p3 = uip.round_to_grid_static_point(p1)
@@ -85,15 +85,15 @@ class Cone(UserObjectGeometry):
x2, y2, z2 = uip.round_to_grid(p2) x2, y2, z2 = uip.round_to_grid(p2)
if r1 < 0: if r1 < 0:
logger.exception(self.__str__() + f' the radius of the first face {r1:g} should be a positive value.') logger.exception(f'{self.__str__()} the radius of the first face {r1:g} should be a positive value.')
raise ValueError raise ValueError
if r2 < 0: if r2 < 0:
logger.exception(self.__str__() + f' the radius of the second face {r2:g} should be a positive value.') logger.exception(f'{self.__str__()} the radius of the second face {r2:g} should be a positive value.')
raise ValueError raise ValueError
if r1 == 0 and r2 == 0: if r1 == 0 and r2 == 0:
logger.exception(self.__str__() + f' not both radii can be zero.') logger.exception(f'{self.__str__()} not both radii can be zero.')
raise ValueError raise ValueError
# Look up requested materials in existing list of material instances # Look up requested materials in existing list of material instances
@@ -101,7 +101,7 @@ class Cone(UserObjectGeometry):
if len(materials) != len(materialsrequested): if len(materials) != len(materialsrequested):
notfound = [x for x in materialsrequested if x not in materials] notfound = [x for x in materialsrequested if x not in materials]
logger.exception(self.__str__() + f' material(s) {notfound} do not exist') logger.exception(f'{self.__str__()} material(s) {notfound} do not exist')
raise ValueError raise ValueError
# Isotropic case # Isotropic case
@@ -141,7 +141,7 @@ class Cone(UserObjectGeometry):
grid.rigidE, grid.rigidH, grid.ID) grid.rigidE, grid.rigidH, grid.ID)
dielectricsmoothing = 'on' if averaging else 'off' dielectricsmoothing = 'on' if averaging else 'off'
logger.info(self.grid_name(grid) + f"Cone with face centres {p3[0]:g}m, " + logger.info(f"{self.grid_name(grid)}Cone with face centres {p3[0]:g}m, " +
f"{p3[1]:g}m, {p3[2]:g}m and {p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m, " + f"{p3[1]:g}m, {p3[2]:g}m and {p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m, " +
f"with radii {r1:g}m and {r2:g}, of material(s) {', '.join(materialsrequested)} " + f"with radii {r1:g}m and {r2:g}, of material(s) {', '.join(materialsrequested)} " +
f"created, dielectric smoothing is {dielectricsmoothing}.") f"created, dielectric smoothing is {dielectricsmoothing}.")

10
gprMax/cmds_geometry/cylinder.py
查看文件

@@ -52,7 +52,7 @@ class Cylinder(UserObjectGeometry):
p2 = self.kwargs['p2'] p2 = self.kwargs['p2']
r = self.kwargs['r'] r = self.kwargs['r']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' please specify 2 points and a radius') logger.exception(f'{self.__str__()} please specify 2 points and a radius')
raise raise
# Check averaging # Check averaging
@@ -72,7 +72,7 @@ class Cylinder(UserObjectGeometry):
try: try:
materialsrequested = self.kwargs['material_ids'] materialsrequested = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' no materials have been specified') logger.exception(f'{self.__str__()} no materials have been specified')
raise raise
p3 = uip.round_to_grid_static_point(p1) p3 = uip.round_to_grid_static_point(p1)
@@ -82,7 +82,7 @@ class Cylinder(UserObjectGeometry):
x2, y2, z2 = uip.round_to_grid(p2) x2, y2, z2 = uip.round_to_grid(p2)
if r <= 0: if r <= 0:
logger.exception(self.__str__() + f' the radius {r:g} should be a positive value.') logger.exception(f'{self.__str__()} the radius {r:g} should be a positive value.')
raise ValueError raise ValueError
# Look up requested materials in existing list of material instances # Look up requested materials in existing list of material instances
@@ -90,7 +90,7 @@ class Cylinder(UserObjectGeometry):
if len(materials) != len(materialsrequested): if len(materials) != len(materialsrequested):
notfound = [x for x in materialsrequested if x not in materials] notfound = [x for x in materialsrequested if x not in materials]
logger.exception(self.__str__() + f' material(s) {notfound} do not exist') logger.exception(f'{self.__str__()} material(s) {notfound} do not exist')
raise ValueError raise ValueError
# Isotropic case # Isotropic case
@@ -130,7 +130,7 @@ class Cylinder(UserObjectGeometry):
grid.rigidE, grid.rigidH, grid.ID) grid.rigidE, grid.rigidH, grid.ID)
dielectricsmoothing = 'on' if averaging else 'off' dielectricsmoothing = 'on' if averaging else 'off'
logger.info(self.grid_name(grid) + f"Cylinder with face centres {p3[0]:g}m, " + logger.info(f"{self.grid_name(grid)}Cylinder with face centres {p3[0]:g}m, " +
f"{p3[1]:g}m, {p3[2]:g}m and {p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m, " + f"{p3[1]:g}m, {p3[2]:g}m and {p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m, " +
f"with radius {r:g}m, of material(s) {', '.join(materialsrequested)} " + f"with radius {r:g}m, of material(s) {', '.join(materialsrequested)} " +
f"created, dielectric smoothing is {dielectricsmoothing}.") f"created, dielectric smoothing is {dielectricsmoothing}.")

26
gprMax/cmds_geometry/cylindrical_sector.py
查看文件

@@ -88,24 +88,25 @@ class CylindricalSector(UserObjectGeometry):
try: try:
materialsrequested = self.kwargs['material_ids'] materialsrequested = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' No materials have been specified') logger.exception(f'{self.__str__()} No materials have been specified')
raise raise
sectorstartangle = 2 * np.pi * (start / 360) sectorstartangle = 2 * np.pi * (start / 360)
sectorangle = 2 * np.pi * (end / 360) sectorangle = 2 * np.pi * (end / 360)
if normal != 'x' and normal != 'y' and normal != 'z': if normal not in ['x', 'y', 'z']:
logger.exception(self.__str__() + ' the normal direction must be either x, y or z.') logger.exception(f'{self.__str__()} the normal direction must be either ' +
f'x, y or z.')
raise ValueError raise ValueError
if r <= 0: if r <= 0:
logger.exception(self.__str__() + f' the radius {r:g} should be a positive value.') logger.exception(f'{self.__str__()} the radius {r:g} should be a positive value.')
if sectorstartangle < 0 or sectorangle <= 0: if sectorstartangle < 0 or sectorangle <= 0:
logger.exception(self.__str__() + ' the starting angle and sector ' + logger.exception(f'{self.__str__()} the starting angle and sector angle should be ' +
'angle should be a positive values.') f'a positive values.')
raise ValueError raise ValueError
if sectorstartangle >= 2 * np.pi or sectorangle >= 2 * np.pi: if sectorstartangle >= 2 * np.pi or sectorangle >= 2 * np.pi:
logger.exception(self.__str__() + ' the starting angle and sector ' + logger.exception(f'{self.__str__()} the starting angle and sector angle must be ' +
'angle must be less than 360 degrees.') f'less than 360 degrees.')
raise ValueError raise ValueError
# Look up requested materials in existing list of material instances # Look up requested materials in existing list of material instances
@@ -113,7 +114,7 @@ class CylindricalSector(UserObjectGeometry):
if len(materials) != len(materialsrequested): if len(materials) != len(materialsrequested):
notfound = [x for x in materialsrequested if x not in materials] notfound = [x for x in materialsrequested if x not in materials]
logger.exception(self.__str__() + f' material(s) {notfound} do not exist') logger.exception(f'{self.__str__()} material(s) {notfound} do not exist')
raise ValueError raise ValueError
if thickness > 0: if thickness > 0:
@@ -122,13 +123,12 @@ class CylindricalSector(UserObjectGeometry):
averaging = materials[0].averagable and averagecylindricalsector averaging = materials[0].averagable and averagecylindricalsector
numID = numIDx = numIDy = numIDz = materials[0].numID numID = numIDx = numIDy = numIDz = materials[0].numID
# Uniaxial anisotropic case
elif len(materials) == 3: elif len(materials) == 3:
averaging = False averaging = False
numIDx = materials[0].numID numIDx = materials[0].numID
numIDy = materials[1].numID numIDy = materials[1].numID
numIDz = materials[2].numID numIDz = materials[2].numID
requiredID = materials[0].ID + '+' + materials[1].ID + '+' + materials[2].ID requiredID = f'{materials[0].ID}+{materials[1].ID}+{materials[2].ID}'
averagedmaterial = [x for x in grid.materials if x.ID == requiredID] averagedmaterial = [x for x in grid.materials if x.ID == requiredID]
if averagedmaterial: if averagedmaterial:
numID = averagedmaterial.numID numID = averagedmaterial.numID
@@ -181,14 +181,14 @@ class CylindricalSector(UserObjectGeometry):
if thickness > 0: if thickness > 0:
dielectricsmoothing = 'on' if averaging else 'off' dielectricsmoothing = 'on' if averaging else 'off'
logger.info(self.grid_name(grid) + f"Cylindrical sector with centre " + logger.info(f"{self.grid_name(grid)}Cylindrical sector with centre " +
f"{ctr1:g}m, {ctr2:g}m, radius {r:g}m, starting angle " + f"{ctr1:g}m, {ctr2:g}m, radius {r:g}m, starting angle " +
f"{(sectorstartangle / (2 * np.pi)) * 360:.1f} degrees, " + f"{(sectorstartangle / (2 * np.pi)) * 360:.1f} degrees, " +
f"sector angle {(sectorangle / (2 * np.pi)) * 360:.1f} degrees, " + f"sector angle {(sectorangle / (2 * np.pi)) * 360:.1f} degrees, " +
f"thickness {thickness:g}m, of material(s) {', '.join(materialsrequested)} " + f"thickness {thickness:g}m, of material(s) {', '.join(materialsrequested)} " +
f"created, dielectric smoothing is {dielectricsmoothing}.") f"created, dielectric smoothing is {dielectricsmoothing}.")
else: else:
logger.info(self.grid_name(grid) + f"Cylindrical sector with centre " + logger.info(f"{self.grid_name(grid)}Cylindrical sector with centre " +
f"{ctr1:g}m, {ctr2:g}m, radius {r:g}m, starting angle " + f"{ctr1:g}m, {ctr2:g}m, radius {r:g}m, starting angle " +
f"{(sectorstartangle / (2 * np.pi)) * 360:.1f} degrees, " + f"{(sectorstartangle / (2 * np.pi)) * 360:.1f} degrees, " +
f"sector angle {(sectorangle / (2 * np.pi)) * 360:.1f} " + f"sector angle {(sectorangle / (2 * np.pi)) * 360:.1f} " +

23
gprMax/cmds_geometry/edge.py
查看文件

@@ -62,7 +62,7 @@ class Edge(UserObjectGeometry):
p2 = self.kwargs['p2'] p2 = self.kwargs['p2']
material_id = self.kwargs['material_id'] material_id = self.kwargs['material_id']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires exactly 3 parameters') logger.exception(f'{self.__str__()} requires exactly 3 parameters')
raise raise
if self.do_rotate: if self.do_rotate:
@@ -83,32 +83,23 @@ class Edge(UserObjectGeometry):
# Check for valid orientations # Check for valid orientations
# x-orientated edge # x-orientated edge
if xs != xf: if ((xs != xf and (ys != yf or zs != zf))
if ys != yf or zs != zf: or (ys != yf and (xs != xf or zs != zf))
logger.exception(self.__str__() + ' the edge is not specified correctly') or (zs != zf and (xs != xf or ys != yf))):
logger.exception(f'{self.__str__()} the edge is not specified correctly')
raise ValueError raise ValueError
else: elif xs != xf:
for i in range(xs, xf): for i in range(xs, xf):
build_edge_x(i, ys, zs, material.numID, grid.rigidE, grid.rigidH, grid.ID) build_edge_x(i, ys, zs, material.numID, grid.rigidE, grid.rigidH, grid.ID)
# y-orientated edge
elif ys != yf: elif ys != yf:
if xs != xf or zs != zf:
logger.exception(self.__str__() + ' the edge is not specified correctly')
raise ValueError
else:
for j in range(ys, yf): for j in range(ys, yf):
build_edge_y(xs, j, zs, material.numID, grid.rigidE, grid.rigidH, grid.ID) build_edge_y(xs, j, zs, material.numID, grid.rigidE, grid.rigidH, grid.ID)
# z-orientated edge
elif zs != zf: elif zs != zf:
if xs != xf or ys != yf:
logger.exception(self.__str__() + ' the edge is not specified correctly')
raise ValueError
else:
for k in range(zs, zf): for k in range(zs, zf):
build_edge_z(xs, ys, k, material.numID, grid.rigidE, grid.rigidH, grid.ID) build_edge_z(xs, ys, k, material.numID, grid.rigidE, grid.rigidH, grid.ID)
logger.info(self.grid_name(grid) + f'Edge from {p3[0]:g}m, {p3[1]:g}m, ' + logger.info(f'{self.grid_name(grid)}Edge from {p3[0]:g}m, {p3[1]:g}m, ' +
f'{p3[2]:g}m, to {p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m of ' + f'{p3[2]:g}m, to {p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m of ' +
f'material {material_id} created.') f'material {material_id} created.')

8
gprMax/cmds_geometry/ellipsoid.py
查看文件

@@ -53,7 +53,7 @@ class Ellipsoid(UserObjectGeometry):
zr = self.kwargs['zr'] zr = self.kwargs['zr']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' please specify a point and the three semiaxes.') logger.exception(f'{self.__str__()} please specify a point and the three semiaxes.')
raise raise
# Check averaging # Check averaging
@@ -73,7 +73,7 @@ class Ellipsoid(UserObjectGeometry):
try: try:
materialsrequested = self.kwargs['material_ids'] materialsrequested = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' no materials have been specified') logger.exception(f'{self.__str__()} no materials have been specified')
raise raise
# Centre of sphere # Centre of sphere
@@ -86,7 +86,7 @@ class Ellipsoid(UserObjectGeometry):
if len(materials) != len(materialsrequested): if len(materials) != len(materialsrequested):
notfound = [x for x in materialsrequested if x not in materials] notfound = [x for x in materialsrequested if x not in materials]
logger.exception(self.__str__() + f' material(s) {notfound} do not exist') logger.exception(f'{self.__str__()} material(s) {notfound} do not exist')
raise ValueError raise ValueError
# Isotropic case # Isotropic case
@@ -126,7 +126,7 @@ class Ellipsoid(UserObjectGeometry):
grid.rigidE, grid.rigidH, grid.ID) grid.rigidE, grid.rigidH, grid.ID)
dielectricsmoothing = 'on' if averaging else 'off' dielectricsmoothing = 'on' if averaging else 'off'
logger.info(self.grid_name(grid) + f"Ellipsoid with centre {p2[0]:g}m, " + logger.info(f"{self.grid_name(grid)}Ellipsoid with centre {p2[0]:g}m, " +
f"{p2[1]:g}m, {p2[2]:g}m, x-semiaxis {xr:g}m, y-semiaxis {yr:g}m and z-semiaxis {zr:g}m of material(s) " + f"{p2[1]:g}m, {p2[2]:g}m, x-semiaxis {xr:g}m, y-semiaxis {yr:g}m and z-semiaxis {zr:g}m of material(s) " +
f"{', '.join(materialsrequested)} created, dielectric " + f"{', '.join(materialsrequested)} created, dielectric " +
f"smoothing is {dielectricsmoothing}.") f"smoothing is {dielectricsmoothing}.")

20
gprMax/cmds_geometry/fractal_box.py
查看文件

@@ -78,13 +78,13 @@ class FractalBox(UserObjectGeometry):
mixing_model_id = self.kwargs['mixing_model_id'] mixing_model_id = self.kwargs['mixing_model_id']
ID = self.kwargs['id'] ID = self.kwargs['id']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' Incorrect parameters') logger.exception(f'{self.__str__()} Incorrect parameters')
raise raise
try: try:
seed = self.kwargs['seed'] seed = self.kwargs['seed']
except KeyError: except KeyError:
logger.warning(self.__str__() + ' no value for seed detected. This ' + logger.warning(f'{self.__str__()} no value for seed detected. This ' +
'means you will get a different fractal distribution ' + 'means you will get a different fractal distribution ' +
'every time the model runs.') 'every time the model runs.')
seed = None seed = None
@@ -109,22 +109,22 @@ class FractalBox(UserObjectGeometry):
xf, yf, zf = p2 xf, yf, zf = p2
if frac_dim < 0: if frac_dim < 0:
logger.exception(self.__str__() + ' requires a positive value for the ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'fractal dimension') 'fractal dimension')
raise ValueError raise ValueError
if weighting[0] < 0: if weighting[0] < 0:
logger.exception(self.__str__() + ' requires a positive value for the ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'fractal weighting in the x direction') 'fractal weighting in the x direction')
raise ValueError raise ValueError
if weighting[1] < 0: if weighting[1] < 0:
logger.exception(self.__str__() + ' requires a positive value for the ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'fractal weighting in the y direction') 'fractal weighting in the y direction')
raise ValueError raise ValueError
if weighting[2] < 0: if weighting[2] < 0:
logger.exception(self.__str__() + ' requires a positive value for the ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'fractal weighting in the z direction') 'fractal weighting in the z direction')
if n_materials < 0: if n_materials < 0:
logger.exception(self.__str__() + ' requires a positive value for the ' + logger.exception(f'{self.__str__()} requires a positive value for the ' +
'number of bins') 'number of bins')
raise ValueError raise ValueError
@@ -136,14 +136,14 @@ class FractalBox(UserObjectGeometry):
if mixingmodel: if mixingmodel:
if nbins == 1: if nbins == 1:
logger.exception(self.__str__() + ' must be used with more than ' + logger.exception(f'{self.__str__()} must be used with more than ' +
'one material from the mixing model.') 'one material from the mixing model.')
raise ValueError raise ValueError
# Create materials from mixing model as number of bins now known # Create materials from mixing model as number of bins now known
# from fractal_box command. # from fractal_box command.
mixingmodel.calculate_properties(nbins, grid) mixingmodel.calculate_properties(nbins, grid)
elif not material: elif not material:
logger.exception(self.__str__() + f' mixing model or material with ' + logger.exception(f'{self.__str__()} mixing model or material with ' +
'ID {mixing_model_id} does not exist') 'ID {mixing_model_id} does not exist')
raise ValueError raise ValueError
@@ -160,7 +160,7 @@ class FractalBox(UserObjectGeometry):
volume.mixingmodel = mixingmodel volume.mixingmodel = mixingmodel
dielectricsmoothing = 'on' if volume.averaging else 'off' dielectricsmoothing = 'on' if volume.averaging else 'off'
logger.info(self.grid_name(grid) + f'Fractal box {volume.ID} from ' + logger.info(f'{self.grid_name(grid)}Fractal box {volume.ID} from ' +
f'{p3[0]:g}m, {p3[1]:g}m, {p3[2]:g}m, to {p4[0]:g}m, ' + f'{p3[0]:g}m, {p3[1]:g}m, {p3[2]:g}m, to {p4[0]:g}m, ' +
f'{p4[1]:g}m, {p4[2]:g}m with {volume.operatingonID}, ' + f'{p4[1]:g}m, {p4[2]:g}m with {volume.operatingonID}, ' +
f'fractal dimension {volume.dimension:g}, fractal weightings ' + f'fractal dimension {volume.dimension:g}, fractal weightings ' +

2
gprMax/cmds_geometry/fractal_box_builder.py
查看文件

@@ -319,7 +319,7 @@ class FractalBoxBuilder(UserObjectGeometry):
else: else:
if volume.nbins == 1: if volume.nbins == 1:
logger.exception(self.__str__() + ' is being used with a ' + logger.exception(f'{self.__str__()} is being used with a ' +
'single material and no modifications, ' + 'single material and no modifications, ' +
'therefore please use a #box command instead.') 'therefore please use a #box command instead.')
raise ValueError raise ValueError

8
gprMax/cmds_geometry/geometry_objects_read.py
查看文件

@@ -46,7 +46,7 @@ class GeometryObjectsRead(UserObjectGeometry):
geofile = self.kwargs['geofile'] geofile = self.kwargs['geofile']
matfile = self.kwargs['matfile'] matfile = self.kwargs['matfile']
except KeyError: except KeyError:
logger.exception(self.__str__() + 'requires exactly five parameters') logger.exception(f'{self.__str__()} requires exactly five parameters')
raise raise
# Discretise the point using uip object. This has different behaviour # Discretise the point using uip object. This has different behaviour
@@ -100,7 +100,7 @@ class GeometryObjectsRead(UserObjectGeometry):
if round_value((dx_dy_dz[0] / grid.dx) != 1 or if round_value((dx_dy_dz[0] / grid.dx) != 1 or
round_value(dx_dy_dz[1] / grid.dy) != 1 or round_value(dx_dy_dz[1] / grid.dy) != 1 or
round_value(dx_dy_dz[2] / grid.dz) != 1): round_value(dx_dy_dz[2] / grid.dz) != 1):
logger.exception(self.__str__() + ' requires the spatial resolution ' + logger.exception(f'{self.__str__()} requires the spatial resolution ' +
'of the geometry objects file to match the spatial ' + 'of the geometry objects file to match the spatial ' +
'resolution of the model') 'resolution of the model')
raise ValueError raise ValueError
@@ -122,7 +122,7 @@ class GeometryObjectsRead(UserObjectGeometry):
grid.rigidE[:, xs:xs + rigidE.shape[1], ys:ys + rigidE.shape[2], zs:zs + rigidE.shape[3]] = rigidE grid.rigidE[:, xs:xs + rigidE.shape[1], ys:ys + rigidE.shape[2], zs:zs + rigidE.shape[3]] = rigidE
grid.rigidH[:, xs:xs + rigidH.shape[1], ys:ys + rigidH.shape[2], zs:zs + rigidH.shape[3]] = rigidH grid.rigidH[:, xs:xs + rigidH.shape[1], ys:ys + rigidH.shape[2], zs:zs + rigidH.shape[3]] = rigidH
grid.ID[:, xs:xs + ID.shape[1], ys:ys + ID.shape[2], zs:zs + ID.shape[3]] = ID + numexistmaterials grid.ID[:, xs:xs + ID.shape[1], ys:ys + ID.shape[2], zs:zs + ID.shape[3]] = ID + numexistmaterials
logger.info(self.grid_name(grid) + f'Geometry objects from file {geofile} ' + logger.info(f'{self.grid_name(grid)}Geometry objects from file {geofile} ' +
f'inserted at {xs * grid.dx:g}m, {ys * grid.dy:g}m, ' + f'inserted at {xs * grid.dx:g}m, {ys * grid.dy:g}m, ' +
f'{zs * grid.dz:g}m, with corresponding materials file ' + f'{zs * grid.dz:g}m, with corresponding materials file ' +
f'{matfile}.') f'{matfile}.')
@@ -131,7 +131,7 @@ class GeometryObjectsRead(UserObjectGeometry):
build_voxels_from_array(xs, ys, zs, config.get_model_config().ompthreads, build_voxels_from_array(xs, ys, zs, config.get_model_config().ompthreads,
numexistmaterials, averaging, data, numexistmaterials, averaging, data,
grid.solid, grid.rigidE, grid.rigidH, grid.ID) grid.solid, grid.rigidE, grid.rigidH, grid.ID)
logger.info(self.grid_name(grid) + f'Geometry objects from file ' + logger.info(f'{self.grid_name(grid)}Geometry objects from file ' +
f'(voxels only){geofile} inserted at {xs * grid.dx:g}m, ' + f'(voxels only){geofile} inserted at {xs * grid.dx:g}m, ' +
f'{ys * grid.dy:g}m, {zs * grid.dz:g}m, with corresponding ' + f'{ys * grid.dy:g}m, {zs * grid.dz:g}m, with corresponding ' +
f'materials file {matfile}.') f'materials file {matfile}.')

29
gprMax/cmds_geometry/plate.py
查看文件

@@ -61,7 +61,7 @@ class Plate(UserObjectGeometry):
p1 = self.kwargs['p1'] p1 = self.kwargs['p1']
p2 = self.kwargs['p2'] p2 = self.kwargs['p2']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' 2 points must be specified') logger.exception(f'{self.__str__()} 2 points must be specified')
raise raise
# isotropic # isotropic
@@ -72,7 +72,7 @@ class Plate(UserObjectGeometry):
try: try:
materialsrequested = self.kwargs['material_ids'] materialsrequested = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' No materials have been specified') logger.exception(f'{self.__str__()} No materials have been specified')
raise raise
if self.do_rotate: if self.do_rotate:
@@ -86,23 +86,10 @@ class Plate(UserObjectGeometry):
xf, yf, zf = p2 xf, yf, zf = p2
# Check for valid orientations # Check for valid orientations
if xs == xf: if ((xs == xf and (ys == yf or zs == zf))
if ys == yf or zs == zf: or (ys == yf and (xs == xf or zs == zf))
logger.exception(self.__str__() + ' the plate is not specified correctly') or (zs == zf and (ys != yf and xs != xf))):
raise ValueError logger.exception(f'{self.__str__()} the plate is not specified correctly')
elif ys == yf:
if xs == xf or zs == zf:
logger.exception(self.__str__() + ' the plate is not specified correctly')
raise ValueError
elif zs == zf:
if xs == xf or ys == yf:
logger.exception(self.__str__() + ' the plate is not specified correctly')
raise ValueError
else:
logger.exception(self.__str__() + ' the plate is not specified correctly')
raise ValueError raise ValueError
# Look up requested materials in existing list of material instances # Look up requested materials in existing list of material instances
@@ -110,7 +97,7 @@ class Plate(UserObjectGeometry):
if len(materials) != len(materialsrequested): if len(materials) != len(materialsrequested):
notfound = [x for x in materialsrequested if x not in materials] notfound = [x for x in materialsrequested if x not in materials]
logger.exception(self.__str__() + f' material(s) {notfound} do not exist') logger.exception(f'{self.__str__()} material(s) {notfound} do not exist')
raise ValueError raise ValueError
# yz-plane plate # yz-plane plate
@@ -161,6 +148,6 @@ class Plate(UserObjectGeometry):
build_face_xy(i, j, zs, numIDx, numIDy, grid.rigidE, build_face_xy(i, j, zs, numIDx, numIDy, grid.rigidE,
grid.rigidH, grid.ID) grid.rigidH, grid.ID)
logger.info(self.grid_name(grid) + f"Plate from {p3[0]:g}m, {p3[1]:g}m, " + logger.info(f"{self.grid_name(grid)}Plate from {p3[0]:g}m, {p3[1]:g}m, " +
f"{p3[2]:g}m, to {p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m of " + f"{p3[2]:g}m, to {p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m of " +
f"material(s) {', '.join(materialsrequested)} created.") f"material(s) {', '.join(materialsrequested)} created.")

8
gprMax/cmds_geometry/sphere.py
查看文件

@@ -48,7 +48,7 @@ class Sphere(UserObjectGeometry):
p1 = self.kwargs['p1'] p1 = self.kwargs['p1']
r = self.kwargs['r'] r = self.kwargs['r']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' please specify a point and a radius.') logger.exception(f'{self.__str__()} please specify a point and a radius.')
raise raise
# Check averaging # Check averaging
@@ -68,7 +68,7 @@ class Sphere(UserObjectGeometry):
try: try:
materialsrequested = self.kwargs['material_ids'] materialsrequested = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' no materials have been specified') logger.exception(f'{self.__str__()} no materials have been specified')
raise raise
# Centre of sphere # Centre of sphere
@@ -85,7 +85,7 @@ class Sphere(UserObjectGeometry):
if len(materials) != len(materialsrequested): if len(materials) != len(materialsrequested):
notfound = [x for x in materialsrequested if x not in materials] notfound = [x for x in materialsrequested if x not in materials]
logger.exception(self.__str__() + f' material(s) {notfound} do not exist') logger.exception(f'{self.__str__()} material(s) {notfound} do not exist')
raise ValueError raise ValueError
# Isotropic case # Isotropic case
@@ -125,7 +125,7 @@ class Sphere(UserObjectGeometry):
grid.rigidE, grid.rigidH, grid.ID) grid.rigidE, grid.rigidH, grid.ID)
dielectricsmoothing = 'on' if averaging else 'off' dielectricsmoothing = 'on' if averaging else 'off'
logger.info(self.grid_name(grid) + f"Sphere with centre {p2[0]:g}m, " + logger.info(f"{self.grid_name(grid)}Sphere with centre {p2[0]:g}m, " +
f"{p2[1]:g}m, {p2[2]:g}m, radius {r:g}m, of material(s) " + f"{p2[1]:g}m, {p2[2]:g}m, radius {r:g}m, of material(s) " +
f"{', '.join(materialsrequested)} created, dielectric " + f"{', '.join(materialsrequested)} created, dielectric " +
f"smoothing is {dielectricsmoothing}.") f"smoothing is {dielectricsmoothing}.")

20
gprMax/cmds_geometry/triangle.py
查看文件

@@ -70,7 +70,7 @@ class Triangle(UserObjectGeometry):
up3 = self.kwargs['p3'] up3 = self.kwargs['p3']
thickness = self.kwargs['thickness'] thickness = self.kwargs['thickness']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' specify 3 points and a thickness') logger.exception(f'{self.__str__()} specify 3 points and a thickness')
raise raise
if self.do_rotate: if self.do_rotate:
@@ -93,7 +93,7 @@ class Triangle(UserObjectGeometry):
try: try:
materialsrequested = self.kwargs['material_ids'] materialsrequested = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' no materials have been specified') logger.exception(f'{self.__str__()} no materials have been specified')
raise raise
p4 = uip.round_to_grid_static_point(up1) p4 = uip.round_to_grid_static_point(up1)
@@ -108,21 +108,21 @@ class Triangle(UserObjectGeometry):
x3, y3, z3 = uip.round_to_grid(up3) x3, y3, z3 = uip.round_to_grid(up3)
if thickness < 0: if thickness < 0:
logger.exception(self.__str__() + ' requires a positive value for thickness') logger.exception(f'{self.__str__()} requires a positive value for thickness')
raise ValueError raise ValueError
# Check for valid orientations # Check for valid orientations
# yz-plane triangle # yz-plane triangle
if x1 == x2 and x2 == x3: if x1 == x2 == x3:
normal = 'x' normal = 'x'
# xz-plane triangle # xz-plane triangle
elif y1 == y2 and y2 == y3: elif y1 == y2 == y3:
normal = 'y' normal = 'y'
# xy-plane triangle # xy-plane triangle
elif z1 == z2 and z2 == z3: elif z1 == z2 == z3:
normal = 'z' normal = 'z'
else: else:
logger.exception(self.__str__() + ' the triangle is not specified correctly') logger.exception(f'{self.__str__()} the triangle is not specified correctly')
raise ValueError raise ValueError
# Look up requested materials in existing list of material instances # Look up requested materials in existing list of material instances
@@ -130,7 +130,7 @@ class Triangle(UserObjectGeometry):
if len(materials) != len(materialsrequested): if len(materials) != len(materialsrequested):
notfound = [x for x in materialsrequested if x not in materials] notfound = [x for x in materialsrequested if x not in materials]
logger.exception(self.__str__() + f' material(s) {notfound} do not exist') logger.exception(f'{self.__str__()} material(s) {notfound} do not exist')
raise ValueError raise ValueError
if thickness > 0: if thickness > 0:
@@ -185,13 +185,13 @@ class Triangle(UserObjectGeometry):
if thickness > 0: if thickness > 0:
dielectricsmoothing = 'on' if averaging else 'off' dielectricsmoothing = 'on' if averaging else 'off'
logger.info(self.grid_name(grid) + f"Triangle with coordinates " + logger.info(f"{self.grid_name(grid)}Triangle with coordinates " +
f"{p4[0]:g}m {p4[1]:g}m {p4[2]:g}m, {p5[0]:g}m {p5[1]:g}m " + f"{p4[0]:g}m {p4[1]:g}m {p4[2]:g}m, {p5[0]:g}m {p5[1]:g}m " +
f"{p5[2]:g}m, {p6[0]:g}m {p6[1]:g}m {p6[2]:g}m and thickness " + f"{p5[2]:g}m, {p6[0]:g}m {p6[1]:g}m {p6[2]:g}m and thickness " +
f"{thickness:g}m of material(s) {', '.join(materialsrequested)} " + f"{thickness:g}m of material(s) {', '.join(materialsrequested)} " +
f"created, dielectric smoothing is {dielectricsmoothing}.") f"created, dielectric smoothing is {dielectricsmoothing}.")
else: else:
logger.info(self.grid_name(grid) + f"Triangle with coordinates " + logger.info(f"{self.grid_name(grid)}Triangle with coordinates " +
f"{p4[0]:g}m {p4[1]:g}m {p4[2]:g}m, {p5[0]:g}m {p5[1]:g}m " + f"{p4[0]:g}m {p4[1]:g}m {p4[2]:g}m, {p5[0]:g}m {p5[1]:g}m " +
f"{p5[2]:g}m, {p6[0]:g}m {p6[1]:g}m {p6[2]:g}m of material(s) " + f"{p5[2]:g}m, {p6[0]:g}m {p6[1]:g}m {p6[2]:g}m of material(s) " +
f"{', '.join(materialsrequested)} created.") f"{', '.join(materialsrequested)} created.")

373
gprMax/cmds_multiuse.py
查看文件

@@ -61,9 +61,9 @@ class UserObjectMulti:
"""Readable user string as per hash commands.""" """Readable user string as per hash commands."""
s = '' s = ''
for _, v in self.kwargs.items(): for _, v in self.kwargs.items():
if isinstance(v, tuple) or isinstance(v, list): if isinstance(v, (tuple, list)):
v = ' '.join([str(el) for el in v]) v = ' '.join([str(el) for el in v])
s += str(v) + ' ' s += f'{str(v)} '
return f'{self.hash}: {s[:-1]}' return f'{self.hash}: {s[:-1]}'
@@ -77,7 +77,7 @@ class UserObjectMulti:
def params_str(self): def params_str(self):
"""Readable string of parameters given to object.""" """Readable string of parameters given to object."""
return self.hash + ': ' + str(self.kwargs) return f'{self.hash}: {str(self.kwargs)}'
def grid_name(self, grid): def grid_name(self, grid):
"""Returns subgrid name for use with logging info. Returns an empty """Returns subgrid name for use with logging info. Returns an empty
@@ -114,12 +114,12 @@ class Waveform(UserObjectMulti):
try: try:
wavetype = self.kwargs['wave_type'].lower() wavetype = self.kwargs['wave_type'].lower()
except KeyError: except KeyError:
logger.exception(self.params_str() + (f" must have one of the " logger.exception(f"{self.params_str()} must have one of the " +
f"following types {','.join(WaveformUser.types)}.")) f"following types {','.join(WaveformUser.types)}.")
raise raise
if wavetype not in WaveformUser.types: if wavetype not in WaveformUser.types:
logger.exception(self.params_str() + (f" must have one of the " logger.exception(f"{self.params_str()} must have one of the " +
f"following types {','.join(WaveformUser.types)}.")) f"following types {','.join(WaveformUser.types)}.")
raise ValueError raise ValueError
if wavetype != 'user': if wavetype != 'user':
@@ -243,18 +243,15 @@ class VoltageSource(UserObjectMulti):
logger.exception(self.params_str() + (' polarisation must be ' logger.exception(self.params_str() + (' polarisation must be '
'x, y, or z.')) 'x, y, or z.'))
raise ValueError raise ValueError
if ('2D TMx' in config.get_model_config().mode and if '2D TMx' in config.get_model_config().mode and polarisation in ['y','z',]:
(polarisation == 'y' or polarisation == 'z')):
logger.exception(self.params_str() + (' polarisation must be x in ' logger.exception(self.params_str() + (' polarisation must be x in '
'2D TMx mode.')) '2D TMx mode.'))
raise ValueError raise ValueError
elif ('2D TMy' in config.get_model_config().mode and elif '2D TMy' in config.get_model_config().mode and polarisation in ['x','z',]:
(polarisation == 'x' or polarisation == 'z')):
logger.exception(self.params_str() + (' polarisation must be y in ' logger.exception(self.params_str() + (' polarisation must be y in '
'2D TMy mode.')) '2D TMy mode.'))
raise ValueError raise ValueError
elif ('2D TMz' in config.get_model_config().mode and elif '2D TMz' in config.get_model_config().mode and polarisation in ['x','y',]:
(polarisation == 'x' or polarisation == 'y')):
logger.exception(self.params_str() + (' polarisation must be z in ' logger.exception(self.params_str() + (' polarisation must be z in '
'2D TMz mode.')) '2D TMz mode.'))
raise ValueError raise ValueError
@@ -305,10 +302,7 @@ class VoltageSource(UserObjectMulti):
'less.')) 'less.'))
raise ValueError raise ValueError
v.start = start v.start = start
if stop > grid.timewindow: v.stop = min(stop, grid.timewindow)
v.stop = grid.timewindow
else:
v.stop = stop
startstop = (f' start time {v.start:g} secs, finish time ' startstop = (f' start time {v.start:g} secs, finish time '
f'{v.stop:g} secs ') f'{v.stop:g} secs ')
except KeyError: except KeyError:
@@ -318,7 +312,7 @@ class VoltageSource(UserObjectMulti):
v.calculate_waveform_values(grid) v.calculate_waveform_values(grid)
logger.info(self.grid_name(grid) + f'Voltage source with polarity ' logger.info(f'{self.grid_name(grid)}Voltage source with polarity '
f'{v.polarisation} at {p2[0]:g}m, {p2[1]:g}m, {p2[2]:g}m, ' f'{v.polarisation} at {p2[0]:g}m, {p2[1]:g}m, {p2[2]:g}m, '
f'resistance {v.resistance:.1f} Ohms,' + startstop + f'resistance {v.resistance:.1f} Ohms,' + startstop +
f'using waveform {v.waveformID} created.') f'using waveform {v.waveformID} created.')
@@ -367,8 +361,7 @@ class HertzianDipole(UserObjectMulti):
p1 = self.kwargs['p1'] p1 = self.kwargs['p1']
waveform_id = self.kwargs['waveform_id'] waveform_id = self.kwargs['waveform_id']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at least 3 ' + logger.exception(f'{self.params_str()} requires at least 3 parameters.')
'parameters.')
raise raise
if self.do_rotate: if self.do_rotate:
@@ -379,18 +372,15 @@ class HertzianDipole(UserObjectMulti):
logger.exception(self.params_str() + (' polarisation must be ' logger.exception(self.params_str() + (' polarisation must be '
'x, y, or z.')) 'x, y, or z.'))
raise ValueError raise ValueError
if ('2D TMx' in config.get_model_config().mode and if '2D TMx' in config.get_model_config().mode and polarisation in ['y','z',]:
(polarisation == 'y' or polarisation == 'z')):
logger.exception(self.params_str() + (' polarisation must be x in ' logger.exception(self.params_str() + (' polarisation must be x in '
'2D TMx mode.')) '2D TMx mode.'))
raise ValueError raise ValueError
elif ('2D TMy' in config.get_model_config().mode and elif '2D TMy' in config.get_model_config().mode and polarisation in ['x','z',]:
(polarisation == 'x' or polarisation == 'z')):
logger.exception(self.params_str() + (' polarisation must be y in ' logger.exception(self.params_str() + (' polarisation must be y in '
'2D TMy mode.')) '2D TMy mode.'))
raise ValueError raise ValueError
elif ('2D TMz' in config.get_model_config().mode and elif '2D TMz' in config.get_model_config().mode and polarisation in ['x','y',]:
(polarisation == 'x' or polarisation == 'y')):
logger.exception(self.params_str() + (' polarisation must be z in ' logger.exception(self.params_str() + (' polarisation must be z in '
'2D TMz mode.')) '2D TMz mode.'))
raise ValueError raise ValueError
@@ -401,7 +391,7 @@ class HertzianDipole(UserObjectMulti):
# Check if there is a waveformID in the waveforms list # Check if there is a waveformID in the waveforms list
if not any(x.ID == waveform_id for x in grid.waveforms): if not any(x.ID == waveform_id for x in grid.waveforms):
logger.exception(self.params_str() + ' there is no waveform ' + logger.exception(f'{self.params_str()} there is no waveform ' +
f'with the identifier {waveform_id}.') f'with the identifier {waveform_id}.')
raise ValueError raise ValueError
@@ -422,8 +412,7 @@ class HertzianDipole(UserObjectMulti):
h.xcoordorigin = xcoord h.xcoordorigin = xcoord
h.ycoordorigin = ycoord h.ycoordorigin = ycoord
h.zcoordorigin = zcoord h.zcoordorigin = zcoord
h.ID = (h.__class__.__name__ + '(' + str(h.xcoord) + ',' + h.ID = f'{h.__class__.__name__}({str(h.xcoord)},{str(h.ycoord)},{str(h.zcoord)})'
str(h.ycoord) + ',' + str(h.zcoord) + ')')
h.waveformID = waveform_id h.waveformID = waveform_id
try: try:
@@ -431,25 +420,19 @@ class HertzianDipole(UserObjectMulti):
start = self.kwargs['start'] start = self.kwargs['start']
stop = self.kwargs['stop'] stop = self.kwargs['stop']
if start < 0: if start < 0:
logger.exception(self.params_str() + (' delay of the initiation ' logger.exception(f'{self.params_str()} delay of the initiation of the ' +
'of the source should not ' f'source should not be less than zero.')
'be less than zero.'))
raise ValueError raise ValueError
if stop < 0: if stop < 0:
logger.exception(self.params_str() + (' time to remove the ' logger.exception(f'{self.params_str()} time to remove the source ' +
'source should not be ' f'should not be less than zero.')
'less than zero.'))
raise ValueError raise ValueError
if stop - start <= 0: if stop - start <= 0:
logger.exception(self.params_str() + (' duration of the source ' logger.exception(
'should not be zero or ' f'{self.params_str()} duration of the source should not be zero or less.')
'less.'))
raise ValueError raise ValueError
h.start = start h.start = start
if stop > grid.timewindow: h.stop = min(stop, grid.timewindow)
h.stop = grid.timewindow
else:
h.stop = stop
startstop = (f' start time {h.start:g} secs, finish time ' startstop = (f' start time {h.start:g} secs, finish time '
f'{h.stop:g} secs ') f'{h.stop:g} secs ')
except KeyError: except KeyError:
@@ -460,12 +443,12 @@ class HertzianDipole(UserObjectMulti):
h.calculate_waveform_values(grid) h.calculate_waveform_values(grid)
if config.get_model_config().mode == '2D': if config.get_model_config().mode == '2D':
logger.info(self.grid_name(grid) + f'Hertzian dipole is a line ' + logger.info(f'{self.grid_name(grid)}Hertzian dipole is a line ' +
f'source in 2D with polarity {h.polarisation} at ' + f'source in 2D with polarity {h.polarisation} at ' +
f'{p2[0]:g}m, {p2[1]:g}m, {p2[2]:g}m,' + startstop + f'{p2[0]:g}m, {p2[1]:g}m, {p2[2]:g}m,' + startstop +
f'using waveform {h.waveformID} created.') f'using waveform {h.waveformID} created.')
else: else:
logger.info(self.grid_name(grid) + f'Hertzian dipole with ' + logger.info(f'{self.grid_name(grid)}Hertzian dipole with ' +
f'polarity {h.polarisation} at {p2[0]:g}m, ' + f'polarity {h.polarisation} at {p2[0]:g}m, ' +
f'{p2[1]:g}m, {p2[2]:g}m,' + startstop + f'{p2[1]:g}m, {p2[2]:g}m,' + startstop +
f'using waveform {h.waveformID} created.') f'using waveform {h.waveformID} created.')
@@ -514,8 +497,7 @@ class MagneticDipole(UserObjectMulti):
p1 = self.kwargs['p1'] p1 = self.kwargs['p1']
waveform_id = self.kwargs['waveform_id'] waveform_id = self.kwargs['waveform_id']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at least five ' logger.exception(f'{self.params_str()} requires at least five parameters.')
'parameters.')
raise raise
if self.do_rotate: if self.do_rotate:
@@ -526,18 +508,15 @@ class MagneticDipole(UserObjectMulti):
logger.exception(self.params_str() + (' polarisation must be ' logger.exception(self.params_str() + (' polarisation must be '
'x, y, or z.')) 'x, y, or z.'))
raise ValueError raise ValueError
if ('2D TMx' in config.get_model_config().mode and if '2D TMx' in config.get_model_config().mode and polarisation in ['y','z',]:
(polarisation == 'y' or polarisation == 'z')):
logger.exception(self.params_str() + (' polarisation must be x in ' logger.exception(self.params_str() + (' polarisation must be x in '
'2D TMx mode.')) '2D TMx mode.'))
raise ValueError raise ValueError
elif ('2D TMy' in config.get_model_config().mode and elif '2D TMy' in config.get_model_config().mode and polarisation in ['x','z',]:
(polarisation == 'x' or polarisation == 'z')):
logger.exception(self.params_str() + (' polarisation must be y in ' logger.exception(self.params_str() + (' polarisation must be y in '
'2D TMy mode.')) '2D TMy mode.'))
raise ValueError raise ValueError
elif ('2D TMz' in config.get_model_config().mode and elif '2D TMz' in config.get_model_config().mode and polarisation in ['x','y',]:
(polarisation == 'x' or polarisation == 'y')):
logger.exception(self.params_str() + (' polarisation must be z in ' logger.exception(self.params_str() + (' polarisation must be z in '
'2D TMz mode.')) '2D TMz mode.'))
raise ValueError raise ValueError
@@ -548,7 +527,7 @@ class MagneticDipole(UserObjectMulti):
# Check if there is a waveformID in the waveforms list # Check if there is a waveformID in the waveforms list
if not any(x.ID == waveform_id for x in grid.waveforms): if not any(x.ID == waveform_id for x in grid.waveforms):
logger.exception(self.params_str() + ' there is no waveform ' + logger.exception(f'{self.params_str()} there is no waveform ' +
f'with the identifier {waveform_id}.') f'with the identifier {waveform_id}.')
raise ValueError raise ValueError
@@ -584,10 +563,7 @@ class MagneticDipole(UserObjectMulti):
'less.')) 'less.'))
raise ValueError raise ValueError
m.start = start m.start = start
if stop > grid.timewindow: m.stop = min(stop, grid.timewindow)
m.stop = grid.timewindow
else:
m.stop = stop
startstop = (f' start time {m.start:g} secs, ' startstop = (f' start time {m.start:g} secs, '
f'finish time {m.stop:g} secs ') f'finish time {m.stop:g} secs ')
except KeyError: except KeyError:
@@ -597,7 +573,7 @@ class MagneticDipole(UserObjectMulti):
m.calculate_waveform_values(grid) m.calculate_waveform_values(grid)
logger.info(self.grid_name(grid) + f'Magnetic dipole with polarity ' + logger.info(f'{self.grid_name(grid)}Magnetic dipole with polarity ' +
f'{m.polarisation} at {p2[0]:g}m, {p2[1]:g}m, {p2[2]:g}m,' + f'{m.polarisation} at {p2[0]:g}m, {p2[1]:g}m, {p2[2]:g}m,' +
startstop + f'using waveform {m.waveformID} created.') startstop + f'using waveform {m.waveformID} created.')
@@ -647,17 +623,15 @@ class TransmissionLine(UserObjectMulti):
waveform_id = self.kwargs['waveform_id'] waveform_id = self.kwargs['waveform_id']
resistance = self.kwargs['resistance'] resistance = self.kwargs['resistance']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at least six ' logger.exception(f'{self.params_str()} requires at least six parameters.')
'parameters.')
raise raise
if self.do_rotate: if self.do_rotate:
self._do_rotate(grid) self._do_rotate(grid)
# Warn about using a transmission line on GPU # Warn about using a transmission line on GPU
if (config.sim_config.general['solver'] == 'cuda' or if (config.sim_config.general['solver'] in ['cuda', 'opencl']):
config.sim_config.general['solver'] == 'opencl'): logger.exception(f'{self.params_str()} cannot currently be used ' +
logger.exception(self.params_str() + ' cannot currently be used ' +
'with the CUDA or OpenCL-based solver. Consider ' + 'with the CUDA or OpenCL-based solver. Consider ' +
'using a #voltage_source instead.') 'using a #voltage_source instead.')
raise ValueError raise ValueError
@@ -667,18 +641,15 @@ class TransmissionLine(UserObjectMulti):
logger.exception(self.params_str() + (' polarisation must be ' logger.exception(self.params_str() + (' polarisation must be '
'x, y, or z.')) 'x, y, or z.'))
raise ValueError raise ValueError
if ('2D TMx' in config.get_model_config().mode and if '2D TMx' in config.get_model_config().mode and polarisation in ['y','z',]:
(polarisation == 'y' or polarisation == 'z')):
logger.exception(self.params_str() + (' polarisation must be x in ' logger.exception(self.params_str() + (' polarisation must be x in '
'2D TMx mode.')) '2D TMx mode.'))
raise ValueError raise ValueError
elif ('2D TMy' in config.get_model_config().mode and elif '2D TMy' in config.get_model_config().mode and polarisation in ['x','z',]:
(polarisation == 'x' or polarisation == 'z')):
logger.exception(self.params_str() + (' polarisation must be y in ' logger.exception(self.params_str() + (' polarisation must be y in '
'2D TMy mode.')) '2D TMy mode.'))
raise ValueError raise ValueError
elif ('2D TMz' in config.get_model_config().mode and elif '2D TMz' in config.get_model_config().mode and polarisation in ['x','y',]:
(polarisation == 'x' or polarisation == 'y')):
logger.exception(self.params_str() + (' polarisation must be z in ' logger.exception(self.params_str() + (' polarisation must be z in '
'2D TMz mode.')) '2D TMz mode.'))
raise ValueError raise ValueError
@@ -688,14 +659,14 @@ class TransmissionLine(UserObjectMulti):
if resistance <= 0 or resistance >= config.sim_config.em_consts['z0']: if resistance <= 0 or resistance >= config.sim_config.em_consts['z0']:
logger.exception(self.params_str() + ' requires a resistance ' + logger.exception(f'{self.params_str()} requires a resistance ' +
'greater than zero and less than the impedance ' + 'greater than zero and less than the impedance ' +
'of free space, i.e. 376.73 Ohms.') 'of free space, i.e. 376.73 Ohms.')
raise ValueError raise ValueError
# Check if there is a waveformID in the waveforms list # Check if there is a waveformID in the waveforms list
if not any(x.ID == waveform_id for x in grid.waveforms): if not any(x.ID == waveform_id for x in grid.waveforms):
logger.exception(self.params_str() + ' there is no waveform ' + logger.exception(f'{self.params_str()} there is no waveform ' +
f'with the identifier {waveform_id}.') f'with the identifier {waveform_id}.')
raise ValueError raise ValueError
@@ -729,11 +700,8 @@ class TransmissionLine(UserObjectMulti):
'less.')) 'less.'))
raise ValueError raise ValueError
t.start = start t.start = start
if stop > grid.timewindow: t.stop = min(stop, grid.timewindow)
t.stop = grid.timewindow startstop = (f' start time {t.start:g} secs, finish time ' +
else:
t.stop = stop
startstop = (f' start time {t.start:g} secs, finish time '
f'{t.stop:g} secs ') f'{t.stop:g} secs ')
except KeyError: except KeyError:
t.start = 0 t.start = 0
@@ -743,7 +711,7 @@ class TransmissionLine(UserObjectMulti):
t.calculate_waveform_values(grid) t.calculate_waveform_values(grid)
t.calculate_incident_V_I(grid) t.calculate_incident_V_I(grid)
logger.info(self.grid_name(grid) + f'Transmission line with polarity ' + logger.info(f'{self.grid_name(grid)}Transmission line with polarity ' +
f'{t.polarisation} at {p2[0]:g}m, {p2[1]:g}m, ' + f'{t.polarisation} at {p2[0]:g}m, {p2[1]:g}m, ' +
f'{p2[2]:g}m, resistance {t.resistance:.1f} Ohms,' + f'{p2[2]:g}m, resistance {t.resistance:.1f} Ohms,' +
startstop + f'using waveform {t.waveformID} created.') startstop + f'using waveform {t.waveformID} created.')
@@ -817,16 +785,14 @@ class Rx(UserObjectMulti):
outputs = [self.kwargs['outputs']] outputs = [self.kwargs['outputs']]
except KeyError: except KeyError:
# If no ID or outputs are specified, use default # If no ID or outputs are specified, use default
r.ID = (r.__class__.__name__ + '(' + str(r.xcoord) + ',' + r.ID = f'{r.__class__.__name__}({str(r.xcoord)},{str(r.ycoord)},{str(r.zcoord)})'
str(r.ycoord) + ',' + str(r.zcoord) + ')')
for key in RxUser.defaultoutputs: for key in RxUser.defaultoutputs:
r.outputs[key] = np.zeros(grid.iterations, r.outputs[key] = np.zeros(grid.iterations,
dtype=config.sim_config.dtypes['float_or_double']) dtype=config.sim_config.dtypes['float_or_double'])
else: else:
outputs.sort() outputs.sort()
# Get allowable outputs # Get allowable outputs
if (config.sim_config.general['solver'] =='cuda' or if config.sim_config.general['solver'] in ['cuda', 'opencl']:
config.sim_config.general['solver'] =='opencl'):
allowableoutputs = RxUser.allowableoutputs_dev allowableoutputs = RxUser.allowableoutputs_dev
else: else:
allowableoutputs = RxUser.allowableoutputs allowableoutputs = RxUser.allowableoutputs
@@ -836,14 +802,14 @@ class Rx(UserObjectMulti):
r.outputs[field] = np.zeros(grid.iterations, r.outputs[field] = np.zeros(grid.iterations,
dtype=config.sim_config.dtypes['float_or_double']) dtype=config.sim_config.dtypes['float_or_double'])
else: else:
logger.exception(self.params_str() + ' contains an ' logger.exception(f'{self.params_str()} contains an output '
'output type that is not allowable. ' f'type that is not allowable. Allowable '
'Allowable outputs in current context are ' f'outputs in current context are '
f'{allowableoutputs}.') f'{allowableoutputs}.')
raise ValueError raise ValueError
logger.info(self.grid_name(grid) + f"Receiver at {p2[0]:g}m, " logger.info(f"{self.grid_name(grid)}Receiver at {p2[0]:g}m, {p2[1]:g}m, "
f"{p2[1]:g}m, {p2[2]:g}m with output component(s) " f"{p2[2]:g}m with output component(s) "
f"{', '.join(r.outputs)} created.") f"{', '.join(r.outputs)} created.")
grid.rxs.append(r) grid.rxs.append(r)
@@ -871,8 +837,7 @@ class RxArray(UserObjectMulti):
p2 = self.kwargs['p2'] p2 = self.kwargs['p2']
dl = self.kwargs['dl'] dl = self.kwargs['dl']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires exactly 9 ' logger.exception(f'{self.params_str()} requires exactly 9 parameters')
'parameters')
raise raise
xs, ys, zs = uip.check_src_rx_point(p1, self.params_str(), 'lower') xs, ys, zs = uip.check_src_rx_point(p1, self.params_str(), 'lower')
@@ -882,39 +847,39 @@ class RxArray(UserObjectMulti):
dx, dy, dz = uip.discretise_point(dl) dx, dy, dz = uip.discretise_point(dl)
if xs > xf or ys > yf or zs > zf: if xs > xf or ys > yf or zs > zf:
logger.exception(self.params_str() + ' the lower coordinates ' logger.exception(f'{self.params_str()} the lower coordinates ' +
'should be less than the upper coordinates.') 'should be less than the upper coordinates.')
raise ValueError raise ValueError
if dx < 0 or dy < 0 or dz < 0: if dx < 0 or dy < 0 or dz < 0:
logger.exception(self.params_str() + ' the step size should not ' logger.exception(f'{self.params_str()} the step size should not ' +
'be less than zero.') 'be less than zero.')
raise ValueError raise ValueError
if dx < 1: if dx < 1:
if dx == 0: if dx == 0:
dx = 1 dx = 1
else: else:
logger.exception(self.params_str() + ' the step size should ' logger.exception(f'{self.params_str()} the step size should ' +
'not be less than the spatial discretisation.') 'not be less than the spatial discretisation.')
raise ValueError raise ValueError
if dy < 1: if dy < 1:
if dy == 0: if dy == 0:
dy = 1 dy = 1
else: else:
logger.exception(self.params_str() + ' the step size should ' logger.exception(f'{self.params_str()} the step size should ' +
'not be less than the spatial discretisation.') 'not be less than the spatial discretisation.')
raise ValueError raise ValueError
if dz < 1: if dz < 1:
if dz == 0: if dz == 0:
dz = 1 dz = 1
else: else:
logger.exception(self.params_str() + ' the step size should ' logger.exception(f'{self.params_str()} the step size should ' +
'not be less than the spatial discretisation.') 'not be less than the spatial discretisation.')
raise ValueError raise ValueError
logger.info(self.grid_name(grid) + f'Receiver array {p3[0]:g}m, ' logger.info(f'{self.grid_name(grid)}Receiver array '
f'{p3[1]:g}m, {p3[2]:g}m, to {p4[0]:g}m, {p4[1]:g}m, ' f'{p3[0]:g}m, {p3[1]:g}m, {p3[2]:g}m, to '
f'{p4[2]:g}m with steps {dx * grid.dx:g}m, ' f'{p4[0]:g}m, {p4[1]:g}m, {p4[2]:g}m with steps '
f'{dy * grid.dy:g}m, {dz * grid.dz:g}m') f'{dx * grid.dx:g}m, {dy * grid.dy:g}m, {dz * grid.dz:g}m')
for x in range(xs, xf + 1, dx): for x in range(xs, xf + 1, dx):
for y in range(ys, yf + 1, dy): for y in range(ys, yf + 1, dy):
@@ -930,8 +895,7 @@ class RxArray(UserObjectMulti):
p5 = np.array([x, y, z]) p5 = np.array([x, y, z])
p5 = uip.descretised_to_continuous(p5) p5 = uip.descretised_to_continuous(p5)
p5 = uip.round_to_grid_static_point(p5) p5 = uip.round_to_grid_static_point(p5)
r.ID = (r.__class__.__name__ + '(' + str(x) + ',' + r.ID = f'{r.__class__.__name__}({str(x)},{str(y)},{str(z)})'
str(y) + ',' + str(z) + ')')
for key in RxUser.defaultoutputs: for key in RxUser.defaultoutputs:
r.outputs[key] = np.zeros(grid.iterations, dtype=config.sim_config.dtypes['float_or_double']) r.outputs[key] = np.zeros(grid.iterations, dtype=config.sim_config.dtypes['float_or_double'])
logger.info(f" Receiver at {p5[0]:g}m, {p5[1]:g}m, " logger.info(f" Receiver at {p5[0]:g}m, {p5[1]:g}m, "
@@ -968,8 +932,7 @@ class Snapshot(UserObjectMulti):
def create(self, grid, uip): def create(self, grid, uip):
if isinstance(grid, SubGridBaseGrid): if isinstance(grid, SubGridBaseGrid):
logger.exception(self.params_str() + ' do not add snapshots to ' logger.exception(f'{self.params_str()} do not add snapshots to subgrids.')
'subgrids.')
raise ValueError raise ValueError
try: try:
p1 = self.kwargs['p1'] p1 = self.kwargs['p1']
@@ -977,8 +940,7 @@ class Snapshot(UserObjectMulti):
dl = self.kwargs['dl'] dl = self.kwargs['dl']
filename = self.kwargs['filename'] filename = self.kwargs['filename']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires exactly 11 ' logger.exception(f'{self.params_str()} requires exactly 11 parameters.')
'parameters.')
raise raise
try: try:
@@ -986,7 +948,7 @@ class Snapshot(UserObjectMulti):
p4 = uip.round_to_grid_static_point(p2) p4 = uip.round_to_grid_static_point(p2)
p1, p2 = uip.check_box_points(p1, p2, self.params_str()) p1, p2 = uip.check_box_points(p1, p2, self.params_str())
except ValueError: except ValueError:
logger.exception(self.params_str() + ' point is outside the domain.') logger.exception(f'{self.params_str()} point is outside the domain.')
raise raise
xs, ys, zs = p1 xs, ys, zs = p1
xf, yf, zf = p2 xf, yf, zf = p2
@@ -1001,13 +963,12 @@ class Snapshot(UserObjectMulti):
try: try:
time = self.kwargs['time'] time = self.kwargs['time']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires exactly 5 ' logger.exception(f'{self.params_str()} requires exactly 5 parameters.')
'parameters.')
raise raise
if time > 0: if time > 0:
iterations = round_value((time / grid.dt)) + 1 iterations = round_value((time / grid.dt)) + 1
else: else:
logger.exception(self.params_str() + ' time value must be ' logger.exception(f'{self.params_str()} time value must be ' +
'greater than zero.') 'greater than zero.')
raise ValueError raise ValueError
@@ -1027,9 +988,9 @@ class Snapshot(UserObjectMulti):
# Check and set output names # Check and set output names
for output in tmp: for output in tmp:
if output not in SnapshotUser.allowableoutputs.keys(): if output not in SnapshotUser.allowableoutputs.keys():
logger.exception(self.params_str() + " contains an " logger.exception(f"{self.params_str()} contains an output "
"output type that is not allowable. " f"type that is not allowable. Allowable "
"Allowable outputs in current context are " f"outputs in current context are "
f"{', '.join(SnapshotUser.allowableoutputs.keys())}.") f"{', '.join(SnapshotUser.allowableoutputs.keys())}.")
raise ValueError raise ValueError
else: else:
@@ -1039,15 +1000,15 @@ class Snapshot(UserObjectMulti):
outputs = dict.fromkeys(SnapshotUser.allowableoutputs, True) outputs = dict.fromkeys(SnapshotUser.allowableoutputs, True)
if dx < 0 or dy < 0 or dz < 0: if dx < 0 or dy < 0 or dz < 0:
logger.exception(self.params_str() + ' the step size should not ' logger.exception(f'{self.params_str()} the step size should not ' +
'be less than zero.') 'be less than zero.')
raise ValueError raise ValueError
if dx < 1 or dy < 1 or dz < 1: if dx < 1 or dy < 1 or dz < 1:
logger.exception(self.params_str() + ' the step size should not ' logger.exception(f'{self.params_str()} the step size should not ' +
'be less than the spatial discretisation.') 'be less than the spatial discretisation.')
raise ValueError raise ValueError
if iterations <= 0 or iterations > grid.iterations: if iterations <= 0 or iterations > grid.iterations:
logger.exception(self.params_str() + ' time value is not valid.') logger.exception(f'{self.params_str()} time value is not valid.')
raise ValueError raise ValueError
s = SnapshotUser(xs, ys, zs, xf, yf, zf, dx, dy, dz, iterations, s = SnapshotUser(xs, ys, zs, xf, yf, zf, dx, dy, dz, iterations,
@@ -1088,39 +1049,37 @@ class Material(UserObjectMulti):
sm = self.kwargs['sm'] sm = self.kwargs['sm']
material_id = self.kwargs['id'] material_id = self.kwargs['id']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires exactly five ' logger.exception(f'{self.params_str()} requires exactly five parameters.')
'parameters.')
raise raise
if er < 1: if er < 1:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value of ' +
'of one or greater for static (DC) permittivity.') f'one or greater for static (DC) permittivity.')
raise ValueError raise ValueError
if se != 'inf': if se != 'inf':
se = float(se) se = float(se)
if se < 0: if se < 0:
logger.exception(self.params_str() + ' requires a positive ' logger.exception(f'{self.params_str()} requires a positive ' +
'value for electric conductivity.') f'value for electric conductivity.')
raise ValueError raise ValueError
else: else:
se = float('inf') se = float('inf')
if mr < 1: if mr < 1:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value of ' +
'of one or greater for magnetic permeability.') f'one or greater for magnetic permeability.')
raise ValueError raise ValueError
if sm < 0: if sm < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for magnetic loss.') f'for magnetic loss.')
raise ValueError raise ValueError
if any(x.ID == material_id for x in grid.materials): if any(x.ID == material_id for x in grid.materials):
logger.exception(self.params_str() + f' with ID {material_id} ' logger.exception(f'{self.params_str()} with ID {material_id} ' +
'already exists') f'already exists')
raise ValueError raise ValueError
# Create a new instance of the Material class material # Create a new instance of the Material class material
# (start index after pec & free_space) # (start index after pec & free_space)
m = MaterialUser(len(grid.materials), material_id) m = MaterialUser(len(grid.materials), material_id)
m.er = er
m.se = se m.se = se
m.mr = mr m.mr = mr
m.sm = sm m.sm = sm
@@ -1129,9 +1088,10 @@ class Material(UserObjectMulti):
if m.se == float('inf'): if m.se == float('inf'):
m.averagable = False m.averagable = False
logger.info(self.grid_name(grid) + f'Material {m.ID} with ' m.er = er
f'eps_r={m.er:g}, sigma={m.se:g} S/m; mu_r={m.mr:g}, ' logger.info(f'{self.grid_name(grid)}Material {m.ID} with eps_r={m.er:g}, '
f'sigma*={m.sm:g} Ohm/m created.') f'sigma={m.se:g} S/m; mu_r={m.mr:g}, sigma*={m.sm:g} Ohm/m '
f'created.')
grid.materials.append(m) grid.materials.append(m)
@@ -1162,12 +1122,11 @@ class AddDebyeDispersion(UserObjectMulti):
tau = self.kwargs['tau'] tau = self.kwargs['tau']
material_ids = self.kwargs['material_ids'] material_ids = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at least four ' logger.exception(f'{self.params_str()} requires at least four parameters.')
'parameters.')
raise raise
if poles < 0: if poles < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for number of poles.') 'for number of poles.')
raise ValueError raise ValueError
@@ -1176,8 +1135,7 @@ class AddDebyeDispersion(UserObjectMulti):
if len(materials) != len(material_ids): if len(materials) != len(material_ids):
notfound = [x for x in material_ids if x not in materials] notfound = [x for x in material_ids if x not in materials]
logger.exception(self.params_str() + f' material(s) {notfound} do ' logger.exception(f'{self.params_str()} material(s) {notfound} do not exist')
'not exist')
raise ValueError raise ValueError
for material in materials: for material in materials:
@@ -1189,14 +1147,14 @@ class AddDebyeDispersion(UserObjectMulti):
disp_material.type = 'debye' disp_material.type = 'debye'
disp_material.poles = poles disp_material.poles = poles
disp_material.averagable = False disp_material.averagable = False
for i in range(0, poles): for i in range(poles):
if tau[i] > 0: if tau[i] > 0:
logger.debug('Not checking if relaxation times are ' logger.debug('Not checking if relaxation times are '
'greater than time-step.') 'greater than time-step.')
disp_material.deltaer.append(er_delta[i]) disp_material.deltaer.append(er_delta[i])
disp_material.tau.append(tau[i]) disp_material.tau.append(tau[i])
else: else:
logger.exception(self.params_str() + ' requires positive ' logger.exception(f'{self.params_str()} requires positive ' +
'values for the permittivity difference.') 'values for the permittivity difference.')
raise ValueError raise ValueError
if disp_material.poles > config.get_model_config().materials['maxpoles']: if disp_material.poles > config.get_model_config().materials['maxpoles']:
@@ -1205,7 +1163,7 @@ class AddDebyeDispersion(UserObjectMulti):
# Replace original material with newly created DispersiveMaterial # Replace original material with newly created DispersiveMaterial
grid.materials = [disp_material if mat.numID==material.numID else mat for mat in grid.materials] grid.materials = [disp_material if mat.numID==material.numID else mat for mat in grid.materials]
logger.info(self.grid_name(grid) + f"Debye disperion added to {disp_material.ID} " logger.info(f"{self.grid_name(grid)}Debye disperion added to {disp_material.ID} "
f"with delta_eps_r={', '.join('%4.2f' % deltaer for deltaer in disp_material.deltaer)}, " f"with delta_eps_r={', '.join('%4.2f' % deltaer for deltaer in disp_material.deltaer)}, "
f"and tau={', '.join('%4.3e' % tau for tau in disp_material.tau)} secs created.") f"and tau={', '.join('%4.3e' % tau for tau in disp_material.tau)} secs created.")
@@ -1238,12 +1196,11 @@ class AddLorentzDispersion(UserObjectMulti):
delta = self.kwargs['delta'] delta = self.kwargs['delta']
material_ids = self.kwargs['material_ids'] material_ids = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at least five ' logger.exception(f'{self.params_str()} requires at least five parameters.')
'parameters.')
raise raise
if poles < 0: if poles < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for number of poles.') 'for number of poles.')
raise ValueError raise ValueError
@@ -1252,8 +1209,7 @@ class AddLorentzDispersion(UserObjectMulti):
if len(materials) != len(material_ids): if len(materials) != len(material_ids):
notfound = [x for x in material_ids if x not in materials] notfound = [x for x in material_ids if x not in materials]
logger.exception(self.params_str() + f' material(s) {notfound} do ' logger.exception(f'{self.params_str()} material(s) {notfound} do not exist')
'not exist')
raise ValueError raise ValueError
for material in materials: for material in materials:
@@ -1265,13 +1221,13 @@ class AddLorentzDispersion(UserObjectMulti):
disp_material.type = 'lorentz' disp_material.type = 'lorentz'
disp_material.poles = poles disp_material.poles = poles
disp_material.averagable = False disp_material.averagable = False
for i in range(0, poles): for i in range(poles):
if er_delta[i] > 0 and omega[i] > grid.dt and delta[i] > grid.dt: if er_delta[i] > 0 and omega[i] > grid.dt and delta[i] > grid.dt:
disp_material.deltaer.append(er_delta[i]) disp_material.deltaer.append(er_delta[i])
disp_material.tau.append(omega[i]) disp_material.tau.append(omega[i])
disp_material.alpha.append(delta[i]) disp_material.alpha.append(delta[i])
else: else:
logger.exception(self.params_str() + ' requires positive ' logger.exception(f'{self.params_str()} requires positive '
'values for the permittivity difference ' 'values for the permittivity difference '
'and frequencies, and associated times ' 'and frequencies, and associated times '
'that are greater than the time step for ' 'that are greater than the time step for '
@@ -1283,9 +1239,10 @@ class AddLorentzDispersion(UserObjectMulti):
# Replace original material with newly created DispersiveMaterial # Replace original material with newly created DispersiveMaterial
grid.materials = [disp_material if mat.numID==material.numID else mat for mat in grid.materials] grid.materials = [disp_material if mat.numID==material.numID else mat for mat in grid.materials]
logger.info(self.grid_name(grid) + f"Lorentz disperion added to {disp_material.ID} " logger.info(
f"with delta_eps_r={', '.join('%4.2f' % deltaer for deltaer in disp_material.deltaer)}, " f"{self.grid_name(grid)}Lorentz disperion added to {disp_material.ID} " +
f"omega={', '.join('%4.3e' % omega for omega in disp_material.tau)} secs, " f"with delta_eps_r={', '.join('%4.2f' % deltaer for deltaer in disp_material.deltaer)}, " +
f"omega={', '.join('%4.3e' % omega for omega in disp_material.tau)} secs, " +
f"and gamma={', '.join('%4.3e' % delta for delta in disp_material.alpha)} created.") f"and gamma={', '.join('%4.3e' % delta for delta in disp_material.alpha)} created.")
@@ -1313,12 +1270,12 @@ class AddDrudeDispersion(UserObjectMulti):
alpha = self.kwargs['alpha'] alpha = self.kwargs['alpha']
material_ids = self.kwargs['material_ids'] material_ids = self.kwargs['material_ids']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at least four ' logger.exception(f'{self.params_str()} requires at least four ' +
'parameters.') 'parameters.')
raise raise
if poles < 0: if poles < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for number of poles.') 'for number of poles.')
raise ValueError raise ValueError
@@ -1327,8 +1284,7 @@ class AddDrudeDispersion(UserObjectMulti):
if len(materials) != len(material_ids): if len(materials) != len(material_ids):
notfound = [x for x in material_ids if x not in materials] notfound = [x for x in material_ids if x not in materials]
logger.exception(self.params_str() + f' material(s) {notfound} do ' logger.exception(f'{self.params_str()} material(s) {notfound} do not exist.')
'not exist.')
raise ValueError raise ValueError
for material in materials: for material in materials:
@@ -1340,14 +1296,14 @@ class AddDrudeDispersion(UserObjectMulti):
disp_material.type = 'drude' disp_material.type = 'drude'
disp_material.poles = poles disp_material.poles = poles
disp_material.averagable = False disp_material.averagable = False
for i in range(0, poles): for i in range(poles):
if omega[i] > 0 and alpha[i] > grid.dt: if omega[i] > 0 and alpha[i] > grid.dt:
disp_material.tau.append(omega[i]) disp_material.tau.append(omega[i])
disp_material.alpha.append(alpha[i]) disp_material.alpha.append(alpha[i])
else: else:
logger.exception(self.params_str() + ' requires positive ' logger.exception(f'{self.params_str()} requires positive ' +
'values for the frequencies, and ' 'values for the frequencies, and ' +
'associated times that are greater than ' 'associated times that are greater than ' +
'the time step for the model.') 'the time step for the model.')
raise ValueError raise ValueError
if disp_material.poles > config.get_model_config().materials['maxpoles']: if disp_material.poles > config.get_model_config().materials['maxpoles']:
@@ -1356,7 +1312,8 @@ class AddDrudeDispersion(UserObjectMulti):
# Replace original material with newly created DispersiveMaterial # Replace original material with newly created DispersiveMaterial
grid.materials = [disp_material if mat.numID==material.numID else mat for mat in grid.materials] grid.materials = [disp_material if mat.numID==material.numID else mat for mat in grid.materials]
logger.info(self.grid_name(grid) + f"Drude disperion added to {disp_material.ID} " logger.info(
f"{self.grid_name(grid)}Drude disperion added to {disp_material.ID} "
f"with omega={', '.join('%4.3e' % omega for omega in disp_material.tau)} secs, " f"with omega={', '.join('%4.3e' % omega for omega in disp_material.tau)} secs, "
f"and gamma={', '.join('%4.3e' % alpha for alpha in disp_material.alpha)} secs created.") f"and gamma={', '.join('%4.3e' % alpha for alpha in disp_material.alpha)} secs created.")
@@ -1392,38 +1349,38 @@ class SoilPeplinski(UserObjectMulti):
water_fraction_upper = self.kwargs['water_fraction_upper'] water_fraction_upper = self.kwargs['water_fraction_upper']
ID = self.kwargs['id'] ID = self.kwargs['id']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at exactly seven ' logger.exception(f'{self.params_str()} requires at exactly seven ' +
'parameters.') 'parameters.')
raise raise
if sand_fraction < 0: if sand_fraction < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the sand fraction.') 'for the sand fraction.')
raise ValueError raise ValueError
if clay_fraction < 0: if clay_fraction < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the clay fraction.') 'for the clay fraction.')
raise ValueError raise ValueError
if bulk_density < 0: if bulk_density < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the bulk density.') 'for the bulk density.')
raise ValueError raise ValueError
if sand_density < 0: if sand_density < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the sand particle density.') 'for the sand particle density.')
raise ValueError raise ValueError
if water_fraction_lower < 0: if water_fraction_lower < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the lower limit of the water volumetric ' 'for the lower limit of the water volumetric '
'fraction.') 'fraction.')
raise ValueError raise ValueError
if water_fraction_upper < 0: if water_fraction_upper < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the upper limit of the water volumetric ' 'for the upper limit of the water volumetric '
'fraction.') 'fraction.')
raise ValueError raise ValueError
if any(x.ID == ID for x in grid.mixingmodels): if any(x.ID == ID for x in grid.mixingmodels):
logger.exception(self.params_str() + f' with ID {ID} already exists') logger.exception(f'{self.params_str()} with ID {ID} already exists')
raise ValueError raise ValueError
# Create a new instance of the Material class material # Create a new instance of the Material class material
@@ -1431,10 +1388,10 @@ class SoilPeplinski(UserObjectMulti):
s = PeplinskiSoilUser(ID, sand_fraction, clay_fraction, bulk_density, s = PeplinskiSoilUser(ID, sand_fraction, clay_fraction, bulk_density,
sand_density, (water_fraction_lower, water_fraction_upper)) sand_density, (water_fraction_lower, water_fraction_upper))
logger.info(self.grid_name(grid) + 'Mixing model (Peplinski) used to ' logger.info(f'{self.grid_name(grid)}Mixing model (Peplinski) used to ' +
f'create {s.ID} with sand fraction {s.S:g}, clay fraction ' f'create {s.ID} with sand fraction {s.S:g}, clay fraction ' +
f'{s.C:g}, bulk density {s.rb:g}g/cm3, sand particle ' f'{s.C:g}, bulk density {s.rb:g}g/cm3, sand particle ' +
f'density {s.rs:g}g/cm3, and water volumetric fraction ' f'density {s.rs:g}g/cm3, and water volumetric fraction ' +
f'{s.mu[0]:g} to {s.mu[1]:g} created.') f'{s.mu[0]:g} to {s.mu[1]:g} created.')
grid.mixingmodels.append(s) grid.mixingmodels.append(s)
@@ -1474,56 +1431,57 @@ class MaterialRange(UserObjectMulti):
ro_upper = self.kwargs['ro_upper'] ro_upper = self.kwargs['ro_upper']
ID = self.kwargs['id'] ID = self.kwargs['id']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at exactly nine ' logger.exception(f'{self.params_str()} requires at exactly nine ' +
'parameters.') 'parameters.')
raise raise
if er_lower < 1: if er_lower < 1:
logger.exception(self.params_str() + ' requires a value greater or equal to 1 ' logger.exception(f'{self.params_str()} requires a value greater or equal to 1 ' +
'for the lower range of relative permittivity.') 'for the lower range of relative permittivity.')
raise ValueError raise ValueError
if mr_lower < 1: if mr_lower < 1:
logger.exception(self.params_str() + ' requires a value greater or equal to 1 ' logger.exception(f'{self.params_str()} requires a value greater or equal to 1 ' +
'for the lower range of relative magnetic permeability.') 'for the lower range of relative magnetic permeability.')
raise ValueError raise ValueError
if sigma_lower < 0: if sigma_lower < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the lower limit of conductivity.') 'for the lower limit of conductivity.')
raise ValueError raise ValueError
if ro_lower < 0: if ro_lower < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the lower range magnetic loss.') 'for the lower range magnetic loss.')
raise ValueError raise ValueError
if er_upper < 1: if er_upper < 1:
logger.exception(self.params_str() + ' requires a value greater or equal to 1' logger.exception(f'{self.params_str()} requires a value greater or equal to 1' +
'for the upper range of relative permittivity.') 'for the upper range of relative permittivity.')
raise ValueError raise ValueError
if mr_upper < 1: if mr_upper < 1:
logger.exception(self.params_str() + ' requires a value greater or equal to 1' logger.exception(f'{self.params_str()} requires a value greater or equal to 1' +
'for the upper range of relative magnetic permeability') 'for the upper range of relative magnetic permeability')
raise ValueError raise ValueError
if sigma_upper < 0: if sigma_upper < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the upper range of conductivity.') 'for the upper range of conductivity.')
raise ValueError raise ValueError
if ro_upper < 0: if ro_upper < 0:
logger.exception(self.params_str() + ' requires a positive value ' logger.exception(f'{self.params_str()} requires a positive value ' +
'for the upper range of magnetic loss.') 'for the upper range of magnetic loss.')
if any(x.ID == ID for x in grid.mixingmodels): if any(x.ID == ID for x in grid.mixingmodels):
logger.exception(self.params_str() + f' with ID {ID} already exists') logger.exception(f'{self.params_str()} with ID {ID} already exists')
raise ValueError raise ValueError
# Create a new instance of the Material class material # Create a new instance of the Material class material
# (start index after pec & free_space) # (start index after pec & free_space)
s = RangeMaterialUser(ID, (er_lower, er_upper), (sigma_lower, sigma_upper), (mr_lower, mr_upper), (ro_lower, ro_upper)) s = RangeMaterialUser(ID, (er_lower, er_upper), (sigma_lower, sigma_upper),
(mr_lower, mr_upper), (ro_lower, ro_upper))
logger.info(self.grid_name(grid) + 'Material properties used to ' logger.info(f'{self.grid_name(grid)}Material properties used to ' +
f'create {s.ID} with range(s) {s.er[0]:g} to {s.er[1]:g}, relative permittivity ' f'create {s.ID} with range(s) {s.er[0]:g} to {s.er[1]:g}, relative permittivity ' +
f'{s.sig[0]:g} to {s.sig[1]:g}, S/m conductivity, {s.mu[0]:g} to {s.mu[1]:g} relative magnetic permeability ' f'{s.sig[0]:g} to {s.sig[1]:g}, S/m conductivity, {s.mu[0]:g} to {s.mu[1]:g} relative magnetic permeability ' +
f'{s.ro[0]:g} to {s.ro[1]:g} Ohm/m magnetic loss, created') f'{s.ro[0]:g} to {s.ro[1]:g} Ohm/m magnetic loss, created')
@@ -1550,27 +1508,25 @@ class MaterialList(UserObjectMulti):
list_of_materials = self.kwargs['list_of_materials'] list_of_materials = self.kwargs['list_of_materials']
ID = self.kwargs['id'] ID = self.kwargs['id']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires at at least 2 ' logger.exception(f'{self.params_str()} requires at at least 2 ' +
'parameters.') 'parameters.')
raise raise
if any(x.ID == ID for x in grid.mixingmodels): if any(x.ID == ID for x in grid.mixingmodels):
logger.exception(self.params_str() + f' with ID {ID} already exists') logger.exception(f'{self.params_str()} with ID {ID} already exists')
raise ValueError raise ValueError
# Create a new instance of the Material class material # Create a new instance of the Material class material
# (start index after pec & free_space) # (start index after pec & free_space)
s = ListMaterialUser(ID, list_of_materials) s = ListMaterialUser(ID, list_of_materials)
logger.info(self.grid_name(grid) + 'A list of materials used to ' logger.info(f'{self.grid_name(grid)}A list of materials used to ' +
f'create {s.ID} that includes {s.mat}, created') f'create {s.ID} that includes {s.mat}, created')
grid.mixingmodels.append(s) grid.mixingmodels.append(s)
class GeometryView(UserObjectMulti): class GeometryView(UserObjectMulti):
"""Outputs to file(s) information about the geometry (mesh) of model. """Outputs to file(s) information about the geometry (mesh) of model.
@@ -1613,7 +1569,7 @@ class GeometryView(UserObjectMulti):
output_type = self.kwargs['output_type'].lower() output_type = self.kwargs['output_type'].lower()
filename = self.kwargs['filename'] filename = self.kwargs['filename']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires exactly eleven ' logger.exception(f'{self.params_str()} requires exactly eleven ' +
'parameters.') 'parameters.')
raise raise
@@ -1624,45 +1580,44 @@ class GeometryView(UserObjectMulti):
p4 = uip.round_to_grid_static_point(p2) p4 = uip.round_to_grid_static_point(p2)
p1, p2 = uip.check_box_points(p1, p2, self.params_str()) p1, p2 = uip.check_box_points(p1, p2, self.params_str())
except ValueError: except ValueError:
logger.exception(self.params_str() + ' point is outside the domain.') logger.exception(f'{self.params_str()} point is outside the domain.')
raise raise
xs, ys, zs = p1 xs, ys, zs = p1
xf, yf, zf = p2 xf, yf, zf = p2
dx, dy, dz = uip.discretise_static_point(dl) dx, dy, dz = uip.discretise_static_point(dl)
if dx < 0 or dy < 0 or dz < 0: if dx < 0 or dy < 0 or dz < 0:
logger.exception(self.params_str() + ' the step size should not be ' logger.exception(f'{self.params_str()} the step size should not be ' +
'less than zero.') 'less than zero.')
raise ValueError raise ValueError
if dx > grid.nx or dy > grid.ny or dz > grid.nz: if dx > grid.nx or dy > grid.ny or dz > grid.nz:
logger.exception(self.params_str() + ' the step size should be ' logger.exception(f'{self.params_str()} the step size should be ' +
'less than the domain size.') 'less than the domain size.')
raise ValueError raise ValueError
if dx < 1 or dy < 1 or dz < 1: if dx < 1 or dy < 1 or dz < 1:
logger.exception(self.params_str() + ' the step size should not ' logger.exception(f'{self.params_str()} the step size should not ' +
'be less than the spatial discretisation.') 'be less than the spatial discretisation.')
raise ValueError raise ValueError
if output_type != 'n' and output_type != 'f': if output_type not in ['n', 'f']:
logger.exception(self.params_str() + ' requires type to be either ' logger.exception(f'{self.params_str()} requires type to be either ' +
'n (normal) or f (fine).') 'n (normal) or f (fine).')
raise ValueError raise ValueError
if output_type == 'f' and (dx * grid.dx != grid.dx or if output_type == 'f' and (dx * grid.dx != grid.dx or
dy * grid.dy != grid.dy or dy * grid.dy != grid.dy or
dz * grid.dz != grid.dz): dz * grid.dz != grid.dz):
logger.exception(self.params_str() + ' requires the spatial ' logger.exception(f'{self.params_str()} requires the spatial ' +
'discretisation for the geometry view to be the ' 'discretisation for the geometry view to be the ' +
'same as the model for geometry view of ' 'same as the model for geometry view of ' +
'type f (fine)') 'type f (fine)')
raise ValueError raise ValueError
g = GeometryViewUser(xs, ys, zs, xf, yf, zf, dx, dy, dz, filename, grid) g = GeometryViewUser(xs, ys, zs, xf, yf, zf, dx, dy, dz, filename, grid)
logger.info(self.grid_name(grid) + f'Geometry view from {p3[0]:g}m, ' logger.info(f'{self.grid_name(grid)}Geometry view from {p3[0]:g}m, ' +
f'{p3[1]:g}m, {p3[2]:g}m, to {p4[0]:g}m, {p4[1]:g}m, ' f'{p3[1]:g}m, {p3[2]:g}m, to {p4[0]:g}m, {p4[1]:g}m, ' +
f'{p4[2]:g}m, discretisation {dx * grid.dx:g}m, ' f'{p4[2]:g}m, discretisation {dx * grid.dx:g}m, ' +
f'{dy * grid.dy:g}m, {dz * grid.dz:g}m, with filename ' f'{dy * grid.dy:g}m, {dz * grid.dz:g}m, with filename ' +
f'base {g.filename} created.') f'base {g.filename} created.')
grid.geometryviews.append(g) grid.geometryviews.append(g)
@@ -1692,7 +1647,7 @@ class GeometryObjectsWrite(UserObjectMulti):
p2 = self.kwargs['p2'] p2 = self.kwargs['p2']
basefilename = self.kwargs['filename'] basefilename = self.kwargs['filename']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires exactly seven ' logger.exception(f'{self.params_str()} requires exactly seven ' +
'parameters.') 'parameters.')
raise raise
@@ -1756,30 +1711,30 @@ class PMLCFS(UserObjectMulti):
sigmamin = self.kwargs['sigmamin'] sigmamin = self.kwargs['sigmamin']
sigmamax = self.kwargs['sigmamax'] sigmamax = self.kwargs['sigmamax']
except KeyError: except KeyError:
logger.exception(self.params_str() + ' requires exactly twelve ' logger.exception(f'{self.params_str()} requires exactly twelve ' +
'parameters.') 'parameters.')
raise raise
if (alphascalingprofile not in CFSParameter.scalingprofiles.keys() or if (alphascalingprofile not in CFSParameter.scalingprofiles.keys() or
kappascalingprofile not in CFSParameter.scalingprofiles.keys() or kappascalingprofile not in CFSParameter.scalingprofiles.keys() or
sigmascalingprofile not in CFSParameter.scalingprofiles.keys()): sigmascalingprofile not in CFSParameter.scalingprofiles.keys()):
logger.exception(self.params_str() + ' must have scaling type ' logger.exception(f'{self.params_str()} must have scaling type ' +
f"{','.join(CFSParameter.scalingprofiles.keys())}") f"{','.join(CFSParameter.scalingprofiles.keys())}")
raise ValueError raise ValueError
if (alphascalingdirection not in CFSParameter.scalingdirections or if (alphascalingdirection not in CFSParameter.scalingdirections or
kappascalingdirection not in CFSParameter.scalingdirections or kappascalingdirection not in CFSParameter.scalingdirections or
sigmascalingdirection not in CFSParameter.scalingdirections): sigmascalingdirection not in CFSParameter.scalingdirections):
logger.exception(self.params_str() + ' must have scaling type ' logger.exception(f'{self.params_str()} must have scaling type ' +
f"{','.join(CFSParameter.scalingdirections)}") f"{','.join(CFSParameter.scalingdirections)}")
raise ValueError raise ValueError
if (float(alphamin) < 0 or float(alphamax) < 0 or if (float(alphamin) < 0 or float(alphamax) < 0 or
float(kappamin) < 0 or float(kappamax) < 0 or float(sigmamin) < 0): float(kappamin) < 0 or float(kappamax) < 0 or float(sigmamin) < 0):
logger.exception(self.params_str() + ' minimum and maximum scaling ' logger.exception(f'{self.params_str()} minimum and maximum scaling ' +
'values must be greater than zero.') 'values must be greater than zero.')
raise ValueError raise ValueError
# TODO: Fix handling of kappa for 2nd order PMLs # TODO: Fix handling of kappa for 2nd order PMLs
# if float(kappamin) < 1: # if float(kappamin) < 1:
# logger.exception(self.params_str() + ' minimum scaling value for ' # logger.exception(f'{self.params_str()} minimum scaling value for '
# 'kappa must be greater than or equal to one.') # 'kappa must be greater than or equal to one.')
# raise ValueError # raise ValueError
@@ -1819,7 +1774,7 @@ class PMLCFS(UserObjectMulti):
grid.pmls['cfs'].append(cfs) grid.pmls['cfs'].append(cfs)
if len(grid.pmls['cfs']) > 2: if len(grid.pmls['cfs']) > 2:
logger.exception(self.params_str() + ' can only be used up to two ' logger.exception(f'{self.params_str()} can only be used up to two ' +
'times, for up to a 2nd order PML.') 'times, for up to a 2nd order PML.')
raise ValueError raise ValueError

49
gprMax/cmds_singleuse.py
查看文件

@@ -93,20 +93,20 @@ class Discretisation(UserObjectSingle):
G.dl = np.array(self.kwargs['p1']) G.dl = np.array(self.kwargs['p1'])
G.dx, G.dy, G.dz = self.kwargs['p1'] G.dx, G.dy, G.dz = self.kwargs['p1']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' discretisation requires a point') logger.exception(f'{self.__str__()} discretisation requires a point')
raise raise
if G.dl[0] <= 0: if G.dl[0] <= 0:
logger.exception(self.__str__() + ' discretisation requires the ' + logger.exception(f'{self.__str__()} discretisation requires the '
'x-direction spatial step to be greater than zero') f'x-direction spatial step to be greater than zero')
raise ValueError raise ValueError
if G.dl[1] <= 0: if G.dl[1] <= 0:
logger.exception(self.__str__() + ' discretisation requires the ' + logger.exception(f'{self.__str__()} discretisation requires the '
'y-direction spatial step to be greater than zero') f'y-direction spatial step to be greater than zero')
raise ValueError raise ValueError
if G.dl[2] <= 0: if G.dl[2] <= 0:
logger.exception(self.__str__() + ' discretisation requires the ' + logger.exception(f'{self.__str__()} discretisation requires the '
'z-direction spatial step to be greater than zero') f'z-direction spatial step to be greater than zero')
raise ValueError raise ValueError
logger.info(f'Spatial discretisation: {G.dl[0]:g} x {G.dl[1]:g} x {G.dl[2]:g}m') logger.info(f'Spatial discretisation: {G.dl[0]:g} x {G.dl[1]:g} x {G.dl[2]:g}m')
@@ -127,11 +127,12 @@ class Domain(UserObjectSingle):
try: try:
G.nx, G.ny, G.nz = uip.discretise_point(self.kwargs['p1']) G.nx, G.ny, G.nz = uip.discretise_point(self.kwargs['p1'])
except KeyError: except KeyError:
logger.exception(self.__str__() + ' please specify a point') logger.exception(f'{self.__str__()} please specify a point')
raise raise
if G.nx == 0 or G.ny == 0 or G.nz == 0: if G.nx == 0 or G.ny == 0 or G.nz == 0:
logger.exception(self.__str__() + ' requires at least one cell in every dimension') logger.exception(f'{self.__str__()} requires at least one cell in '
f'every dimension')
raise ValueError raise ValueError
logger.info(f"Domain size: {self.kwargs['p1'][0]:g} x {self.kwargs['p1'][1]:g} x " + logger.info(f"Domain size: {self.kwargs['p1'][0]:g} x {self.kwargs['p1'][1]:g} x " +
@@ -181,12 +182,12 @@ class TimeStepStabilityFactor(UserObjectSingle):
try: try:
f = self.kwargs['f'] f = self.kwargs['f']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires exactly one parameter') logger.exception(f'{self.__str__()} requires exactly one parameter')
raise raise
if f <= 0 or f > 1: if f <= 0 or f > 1:
logger.exception(self.__str__() + ' requires the value of the time ' + logger.exception(f'{self.__str__()} requires the value of the time '
'step stability factor to be between zero and one') f'step stability factor to be between zero and one')
raise ValueError raise ValueError
G.dt = G.dt * f G.dt = G.dt * f
@@ -250,12 +251,12 @@ class OMPThreads(UserObjectSingle):
try: try:
n = self.kwargs['n'] n = self.kwargs['n']
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires exactly one parameter ' + logger.exception(f'{self.__str__()} requires exactly one parameter '
'to specify the number of CPU OpenMP threads to use') f'to specify the number of CPU OpenMP threads to use')
raise raise
if n < 1: if n < 1:
logger.exception(self.__str__() + ' requires the value to be an ' + logger.exception(f'{self.__str__()} requires the value to be an '
'integer not less than one') f'integer not less than one')
raise ValueError raise ValueError
config.get_model_config().ompthreads = set_omp_threads(n) config.get_model_config().ompthreads = set_omp_threads(n)
@@ -302,7 +303,7 @@ class PMLProps(UserObjectSingle):
G.pmls['thickness']['ymax'] = int(self.kwargs['ymax']) G.pmls['thickness']['ymax'] = int(self.kwargs['ymax'])
G.pmls['thickness']['zmax'] = int(self.kwargs['zmax']) G.pmls['thickness']['zmax'] = int(self.kwargs['zmax'])
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires either one or six parameter(s)') logger.exception(f'{self.__str__()} requires either one or six parameter(s)')
raise raise
if (2 * G.pmls['thickness']['x0'] >= G.nx or if (2 * G.pmls['thickness']['x0'] >= G.nx or
@@ -311,7 +312,7 @@ class PMLProps(UserObjectSingle):
2 * G.pmls['thickness']['xmax'] >= G.nx or 2 * G.pmls['thickness']['xmax'] >= G.nx or
2 * G.pmls['thickness']['ymax'] >= G.ny or 2 * G.pmls['thickness']['ymax'] >= G.ny or
2 * G.pmls['thickness']['zmax'] >= G.nz): 2 * G.pmls['thickness']['zmax'] >= G.nz):
logger.exception(self.__str__() + ' has too many cells for the domain size') logger.exception(f'{self.__str__()} has too many cells for the domain size')
raise ValueError raise ValueError
@@ -330,12 +331,12 @@ class SrcSteps(UserObjectSingle):
try: try:
G.srcsteps = uip.discretise_point(self.kwargs['p1']) G.srcsteps = uip.discretise_point(self.kwargs['p1'])
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires exactly three parameters') logger.exception(f'{self.__str__()} requires exactly three parameters')
raise raise
logger.info(f'Simple sources will step {G.srcsteps[0] * G.dx:g}m, ' + logger.info(f'Simple sources will step {G.srcsteps[0] * G.dx:g}m, ' +
f'{G.srcsteps[1] * G.dy:g}m, {G.srcsteps[2] * G.dz:g}m ' + f'{G.srcsteps[1] * G.dy:g}m, {G.srcsteps[2] * G.dz:g}m ' +
f'for each model run.') 'for each model run.')
class RxSteps(UserObjectSingle): class RxSteps(UserObjectSingle):
@@ -353,12 +354,12 @@ class RxSteps(UserObjectSingle):
try: try:
G.rxsteps = uip.discretise_point(self.kwargs['p1']) G.rxsteps = uip.discretise_point(self.kwargs['p1'])
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires exactly three parameters') logger.exception(f'{self.__str__()} requires exactly three parameters')
raise raise
logger.info(f'All receivers will step {G.rxsteps[0] * G.dx:g}m, ' + logger.info(f'All receivers will step {G.rxsteps[0] * G.dx:g}m, ' +
f'{G.rxsteps[1] * G.dy:g}m, {G.rxsteps[2] * G.dz:g}m ' + f'{G.rxsteps[1] * G.dy:g}m, {G.rxsteps[2] * G.dz:g}m ' +
f'for each model run.') 'for each model run.')
class ExcitationFile(UserObjectSingle): class ExcitationFile(UserObjectSingle):
@@ -378,7 +379,7 @@ class ExcitationFile(UserObjectSingle):
def create(self, G, uip): def create(self, G, uip):
try: try:
kwargs = dict() kwargs = {}
excitationfile = self.kwargs['filepath'] excitationfile = self.kwargs['filepath']
kwargs['kind'] = self.kwargs['kind'] kwargs['kind'] = self.kwargs['kind']
kwargs['fill_value'] = self.kwargs['fill_value'] kwargs['fill_value'] = self.kwargs['fill_value']
@@ -389,7 +390,7 @@ class ExcitationFile(UserObjectSingle):
args, varargs, keywords, defaults = inspect.getargspec(interpolate.interp1d) args, varargs, keywords, defaults = inspect.getargspec(interpolate.interp1d)
kwargs = dict(zip(reversed(args), reversed(defaults))) kwargs = dict(zip(reversed(args), reversed(defaults)))
except KeyError: except KeyError:
logger.exception(self.__str__() + ' requires either one or three parameter(s)') logger.exception(f'{self.__str__()} requires either one or three parameter(s)')
raise raise
# See if file exists at specified path and if not try input file directory # See if file exists at specified path and if not try input file directory

17
gprMax/config.py
查看文件

@@ -224,14 +224,14 @@ class SimulationConfig:
# solver: cpu, cuda, opencl. # solver: cpu, cuda, opencl.
# subgrid: whether the simulation uses sub-grids. # subgrid: whether the simulation uses sub-grids.
# precision: data type for electromagnetic field output (single/double). # precision: data type for electromagnetic field output (single/double).
# progressbars: progress bars on stdoout or not - switch off
# progressbars when logging level is greater than
# info (20)
self.general = {'solver': 'cpu', self.general = {'solver': 'cpu',
'subgrid': False, 'subgrid': False,
'precision': 'single'} 'precision': 'single',
'progressbars': args.log_level <= 20}
# Progress bars on stdoout or not - switch off progressbars when
# logging level is greater than info (20)
self.general['progressbars'] = False if args.log_level > 20 else True
self.em_consts = {'c': c, # Speed of light in free space (m/s) self.em_consts = {'c': c, # Speed of light in free space (m/s)
'e0': e0, # Permittivity of free space (F/m) 'e0': e0, # Permittivity of free space (F/m)
@@ -289,7 +289,6 @@ class SimulationConfig:
# Set more complex parameters # Set more complex parameters
self._set_precision() self._set_precision()
self._get_byteorder()
self._set_input_file_path() self._set_input_file_path()
self._set_model_start_end() self._set_model_start_end()
@@ -348,12 +347,6 @@ class SimulationConfig:
elif self.general['solver'] == 'opencl': elif self.general['solver'] == 'opencl':
self.dtypes['C_complex'] = 'cdouble' self.dtypes['C_complex'] = 'cdouble'
def _get_byteorder(self):
"""Checks the byte order of system to use for VTK files, i.e. geometry
views and snapshots.
"""
self.vtk_byteorder = 'LittleEndian' if sys.byteorder == 'little' else 'BigEndian'
def _set_model_start_end(self): def _set_model_start_end(self):
"""Sets range for number of models to run (internally 0 index).""" """Sets range for number of models to run (internally 0 index)."""
if self.args.i: if self.args.i:

14
gprMax/contexts.py
查看文件

@@ -92,7 +92,7 @@ class Context:
def print_logo_copyright(self): def print_logo_copyright(self):
"""Prints gprMax logo, version, and copyright/licencing information.""" """Prints gprMax logo, version, and copyright/licencing information."""
logo_copyright = logo(__version__ + ' (' + codename + ')') logo_copyright = logo(f'{__version__} ({codename})')
logger.basic(logo_copyright) logger.basic(logo_copyright)
def print_sim_time_taken(self): def print_sim_time_taken(self):
@@ -164,19 +164,15 @@ class MPIContext(Context):
executor = self.MPIExecutor(self._run_model, comm=self.comm) executor = self.MPIExecutor(self._run_model, comm=self.comm)
# Check GPU resources versus number of MPI tasks # Check GPU resources versus number of MPI tasks
if executor.is_master(): if (executor.is_master() and
if config.sim_config.general['solver'] == 'cuda': config.sim_config.general['solver'] == 'cuda' and
if executor.size - 1 > len(config.sim_config.devices['devs']): executor.size - 1 > len(config.sim_config.devices['devs'])):
logger.exception('Not enough GPU resources for number of ' logger.exception('Not enough GPU resources for number of '
'MPI tasks requested. Number of MPI tasks ' 'MPI tasks requested. Number of MPI tasks '
'should be equal to number of GPUs + 1.') 'should be equal to number of GPUs + 1.')
raise ValueError raise ValueError
# Create job list jobs = [{'i': i} for i in self.model_range]
jobs = []
for i in self.model_range:
jobs.append({'i': i})
# Send the workers to their work loop # Send the workers to their work loop
executor.start() executor.start()
if executor.is_master(): if executor.is_master():

70
gprMax/cython/pml_build.pyx 普通文件
查看文件

@@ -0,0 +1,70 @@
# Copyright (C) 2015-2023: The University of Edinburgh, United Kingdom
# Authors: Craig Warren, Antonis Giannopoulos, and John Hartley
#
# This file is part of gprMax.
#
# gprMax is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# gprMax is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with gprMax. If not, see <http://www.gnu.org/licenses/>.
cpdef pml_average_er_mr(
str dir,
int s,
G
):
"""Calculates average permittivity and permeability for building PML
(based on underlying material er and mr from solid array).
Args:
dir: string identifier for direction of PML.
s: int for starting cell of PML.
G: FDTDGrid class describing a grid in a model.
"""
sumer = 0 # Sum of relative permittivities in PML slab
summr = 0 # Sum of relative permeabilities in PML slab
if dir == 'x':
for j in range(G.ny):
for k in range(G.nz):
numID = G.solid[s, j, k]
material = [x for x in G.materials if x.numID == numID]
material = material[0]
sumer += material.er
summr += material.mr
averageer = sumer / (G.ny * G.nz)
averagemr = summr / (G.ny * G.nz)
elif dir == 'y':
for i in range(G.nx):
for k in range(G.nz):
numID = G.solid[i, s, k]
material = [x for x in G.materials if x.numID == numID]
material = material[0]
sumer += material.er
summr += material.mr
averageer = sumer / (G.nx * G.nz)
averagemr = summr / (G.nx * G.nz)
elif dir == 'z':
for i in range(G.nx):
for j in range(G.ny):
numID = G.solid[i, j, s]
material = [x for x in G.materials if x.numID == numID]
material = material[0]
sumer += material.er
summr += material.mr
averageer = sumer / (G.nx * G.ny)
averagemr = summr / (G.nx * G.ny)
return averageer, averagemr

4
gprMax/fields_outputs.py
查看文件

@@ -76,7 +76,7 @@ def write_hdf5_outputfile(outputfile, G):
# Write meta data and data for any subgrids # Write meta data and data for any subgrids
if sg_rxs: if sg_rxs:
for sg in G.subgrids: for sg in G.subgrids:
grp = f.create_group('/subgrids/' + sg.name) grp = f.create_group(f'/subgrids/{sg.name}')
write_hd5_data(grp, sg, is_subgrid=True) write_hd5_data(grp, sg, is_subgrid=True)
if G.rxs or sg_rxs: if G.rxs or sg_rxs:
@@ -115,7 +115,7 @@ def write_hd5_data(basegrp, G, is_subgrid=False):
# Create group for sources (except transmission lines); add type and positional data attributes # Create group for sources (except transmission lines); add type and positional data attributes
srclist = G.voltagesources + G.hertziandipoles + G.magneticdipoles srclist = G.voltagesources + G.hertziandipoles + G.magneticdipoles
for srcindex, src in enumerate(srclist): for srcindex, src in enumerate(srclist):
grp = basegrp.create_group('srcs/src' + str(srcindex + 1)) grp = basegrp.create_group(f'srcs/src{str(srcindex + 1)}')
grp.attrs['Type'] = type(src).__name__ grp.attrs['Type'] = type(src).__name__
grp.attrs['Position'] = (src.xcoord * G.dx, src.ycoord * G.dy, src.zcoord * G.dz) grp.attrs['Position'] = (src.xcoord * G.dx, src.ycoord * G.dy, src.zcoord * G.dz)

17
gprMax/geometry_outputs.py
查看文件

@@ -277,14 +277,10 @@ class Comments():
""" """
# Comments for Paraview macro # Comments for Paraview macro
comments = {} comments = {'gprMax_version': __version__,
'dx_dy_dz': self.dx_dy_dz_comment(),
comments['gprMax_version'] = __version__ 'nx_ny_nz': self.nx_ny_nz_comment(),
comments['dx_dy_dz'] = self.dx_dy_dz_comment() 'Materials': self.materials_comment()} # Write the name and numeric ID for each material
comments['nx_ny_nz'] = self.nx_ny_nz_comment()
# Write the name and numeric ID for each material
comments['Materials'] = self.materials_comment()
# Write information on PMLs, sources, and receivers # Write information on PMLs, sources, and receivers
if not self.materials_only: if not self.materials_only:
@@ -345,7 +341,7 @@ class Comments():
def materials_comment(self): def materials_comment(self):
if not self.averaged_materials: if not self.averaged_materials:
return [m.ID for m in self.grid.materials if '+' not in m.ID] return [m.ID for m in self.grid.materials if m.type is not 'dielectric-smoothed']
else: else:
return [m.ID for m in self.grid.materials] return [m.ID for m in self.grid.materials]
@@ -375,8 +371,7 @@ class GeometryObjects:
parts = config.sim_config.input_file_path.with_suffix('').parts parts = config.sim_config.input_file_path.with_suffix('').parts
self.filename_hdf5 = Path(*parts[:-1], self.basefilename) self.filename_hdf5 = Path(*parts[:-1], self.basefilename)
self.filename_hdf5 = self.filename_hdf5.with_suffix('.h5') self.filename_hdf5 = self.filename_hdf5.with_suffix('.h5')
self.filename_materials = Path( self.filename_materials = Path(*parts[:-1], f'{self.basefilename}_materials')
*parts[:-1], self.basefilename + '_materials')
self.filename_materials = self.filename_materials.with_suffix('.txt') self.filename_materials = self.filename_materials.with_suffix('.txt')
# Sizes of arrays to write necessary to update progress bar # Sizes of arrays to write necessary to update progress bar

13
gprMax/hash_cmds_file.py
查看文件

@@ -178,8 +178,8 @@ def write_processed_file(processedlines):
for item in processedlines: for item in processedlines:
f.write(f'{item}') f.write(f'{item}')
logger.info(f'Written input commands, after processing any Python code and ' + logger.info(f'Written input commands, after processing any Python ' +
f'include commands, to file: {processedfile}\n') f'code and include commands, to file: {processedfile}\n')
def check_cmd_names(processedlines, checkessential=True): def check_cmd_names(processedlines, checkessential=True):
@@ -212,6 +212,7 @@ def check_cmd_names(processedlines, checkessential=True):
# - these will be lists within the dictionary # - these will be lists within the dictionary
multiplecmds = {key: [] for key in ['#geometry_view', multiplecmds = {key: [] for key in ['#geometry_view',
'#geometry_objects_write', '#material', '#geometry_objects_write', '#material',
'#material_range', '#material_list',
'#soil_peplinski', '#soil_peplinski',
'#add_dispersion_debye', '#add_dispersion_debye',
'#add_dispersion_lorentz', '#add_dispersion_lorentz',
@@ -219,14 +220,14 @@ def check_cmd_names(processedlines, checkessential=True):
'#waveform', '#voltage_source', '#waveform', '#voltage_source',
'#hertzian_dipole', '#magnetic_dipole', '#hertzian_dipole', '#magnetic_dipole',
'#transmission_line', '#rx', '#rx_array', '#transmission_line', '#rx', '#rx_array',
'#snapshot', '#include_file', '#material_range', '#material_list']} '#snapshot', '#include_file']}
# Geometry object building commands that there can be multiple instances # Geometry object building commands that there can be multiple instances
# of in a model - these will be lists within the dictionary # of in a model - these will be lists within the dictionary
geometrycmds = ['#geometry_objects_read', '#edge', '#plate', '#triangle', geometrycmds = ['#geometry_objects_read', '#edge', '#plate', '#triangle',
'#box', '#sphere', '#ellipsoid', '#cone', '#cylinder', '#cylindrical_sector', '#box', '#sphere', '#ellipsoid', '#cone', '#cylinder',
'#fractal_box', '#add_surface_roughness', '#cylindrical_sector', '#fractal_box',
'#add_surface_water', '#add_grass'] '#add_surface_roughness', '#add_surface_water', '#add_grass']
# List to store all geometry object commands in order from input file # List to store all geometry object commands in order from input file
geometry = [] geometry = []

23
gprMax/hash_cmds_geometry.py
查看文件

@@ -17,21 +17,22 @@
# along with gprMax. If not, see <http://www.gnu.org/licenses/>. # along with gprMax. If not, see <http://www.gnu.org/licenses/>.
import logging import logging
import numpy as np import numpy as np
from .cmds_geometry.add_grass import AddGrass from .cmds_geometry.add_grass import AddGrass
from .cmds_geometry.add_surface_roughness import AddSurfaceRoughness from .cmds_geometry.add_surface_roughness import AddSurfaceRoughness
from .cmds_geometry.add_surface_water import AddSurfaceWater from .cmds_geometry.add_surface_water import AddSurfaceWater
from .cmds_geometry.box import Box from .cmds_geometry.box import Box
from .cmds_geometry.cylinder import Cylinder
from .cmds_geometry.cone import Cone from .cmds_geometry.cone import Cone
from .cmds_geometry.cylinder import Cylinder
from .cmds_geometry.cylindrical_sector import CylindricalSector from .cmds_geometry.cylindrical_sector import CylindricalSector
from .cmds_geometry.edge import Edge from .cmds_geometry.edge import Edge
from .cmds_geometry.ellipsoid import Ellipsoid
from .cmds_geometry.fractal_box import FractalBox from .cmds_geometry.fractal_box import FractalBox
from .cmds_geometry.plate import Plate from .cmds_geometry.plate import Plate
from .cmds_geometry.sphere import Sphere from .cmds_geometry.sphere import Sphere
from .cmds_geometry.triangle import Triangle from .cmds_geometry.triangle import Triangle
from .cmds_geometry.ellipsoid import Ellipsoid
from .utilities.utilities import round_value from .utilities.utilities import round_value
logger = logging.getLogger(__name__) logger = logging.getLogger(__name__)
@@ -76,7 +77,8 @@ def process_geometrycmds(geometry):
tmp = object.split() tmp = object.split()
if tmp[0] == '#geometry_objects_read:': if tmp[0] == '#geometry_objects_read:':
from .cmds_geometry.geometry_objects_read import GeometryObjectsRead from .cmds_geometry.geometry_objects_read import \
GeometryObjectsRead
if len(tmp) != 6: if len(tmp) != 6:
logger.exception("'" + ' '.join(tmp) + "'" + logger.exception("'" + ' '.join(tmp) + "'" +
@@ -219,7 +221,6 @@ def process_geometrycmds(geometry):
scene_objects.append(cylinder) scene_objects.append(cylinder)
elif tmp[0] == '#cone:': elif tmp[0] == '#cone:':
if len(tmp) < 10: if len(tmp) < 10:
logger.exception("'" + ' '.join(tmp) + "'" + logger.exception("'" + ' '.join(tmp) + "'" +
@@ -252,7 +253,6 @@ def process_geometrycmds(geometry):
scene_objects.append(cone) scene_objects.append(cone)
elif tmp[0] == '#cylindrical_sector:': elif tmp[0] == '#cylindrical_sector:':
if len(tmp) < 10: if len(tmp) < 10:
logger.exception("'" + ' '.join(tmp) + "'" + logger.exception("'" + ' '.join(tmp) + "'" +
@@ -331,7 +331,6 @@ def process_geometrycmds(geometry):
scene_objects.append(sphere) scene_objects.append(sphere)
elif tmp[0] == '#ellipsoid:': elif tmp[0] == '#ellipsoid:':
if len(tmp) < 8: if len(tmp) < 8:
logger.exception("'" + ' '.join(tmp) + "'" + logger.exception("'" + ' '.join(tmp) + "'" +
@@ -345,17 +344,19 @@ def process_geometrycmds(geometry):
# Isotropic case with no user specified averaging # Isotropic case with no user specified averaging
if len(tmp) == 8: if len(tmp) == 8:
ellipsoid = Ellipsoid(p1=p1, xr=xr, yr=yr, zr=zr, material_id=tmp[7]) ellipsoid = Ellipsoid(p1=p1, xr=xr, yr=yr, zr=zr,
material_id=tmp[7])
# Isotropic case with user specified averaging # Isotropic case with user specified averaging
elif len(tmp) == 9: elif len(tmp) == 9:
averaging = check_averaging(tmp[8].lower()) averaging = check_averaging(tmp[8].lower())
ellipsoid = Ellipsoid(p1=p1, xr=xr, yr=yr, zr=zr, material_id=tmp[7], ellipsoid = Ellipsoid(p1=p1, xr=xr, yr=yr, zr=zr,
averaging=averaging) material_id=tmp[7], averaging=averaging)
# Uniaxial anisotropic case # Uniaxial anisotropic case
elif len(tmp) == 8: elif len(tmp) == 8:
ellipsoid = Ellipsoid(p1=p1, xr=xr, yr=yr, zr=zr, material_id=tmp[7:]) ellipsoid = Ellipsoid(p1=p1, xr=xr, yr=yr, zr=zr,
material_id=tmp[7:])
else: else:
logger.exception("'" + ' '.join(tmp) + "'" + logger.exception("'" + ' '.join(tmp) + "'" +
@@ -364,8 +365,6 @@ def process_geometrycmds(geometry):
scene_objects.append(ellipsoid) scene_objects.append(ellipsoid)
elif tmp[0] == '#fractal_box:': elif tmp[0] == '#fractal_box:':
# Default is no dielectric smoothing for a fractal box # Default is no dielectric smoothing for a fractal box

11
gprMax/hash_cmds_multiuse.py
查看文件

@@ -21,8 +21,9 @@ import logging
from .cmds_multiuse import (AddDebyeDispersion, AddDrudeDispersion, from .cmds_multiuse import (AddDebyeDispersion, AddDrudeDispersion,
AddLorentzDispersion, GeometryObjectsWrite, AddLorentzDispersion, GeometryObjectsWrite,
GeometryView, HertzianDipole, MagneticDipole, GeometryView, HertzianDipole, MagneticDipole,
Material, Rx, RxArray, Snapshot, SoilPeplinski, Material, MaterialList, MaterialRange, Rx, RxArray,
TransmissionLine, VoltageSource, Waveform, MaterialRange, MaterialList) Snapshot, SoilPeplinski, TransmissionLine,
VoltageSource, Waveform)
logger = logging.getLogger(__name__) logger = logging.getLogger(__name__)
@@ -362,7 +363,6 @@ def process_multicmds(multicmds):
id=tmp[8]) id=tmp[8])
scene_objects.append(material_range) scene_objects.append(material_range)
cmdname = '#material_list' cmdname = '#material_list'
if multicmds[cmdname] is not None: if multicmds[cmdname] is not None:
for cmdinstance in multicmds[cmdname]: for cmdinstance in multicmds[cmdname]:
@@ -375,14 +375,11 @@ def process_multicmds(multicmds):
tokens = len(tmp) tokens = len(tmp)
lmats = [] lmats = []
for iter in range(0,tokens-1): for iter in range(tokens-1):
lmats.append(tmp[iter]) lmats.append(tmp[iter])
material_list = MaterialList(list_of_materials=lmats, material_list = MaterialList(list_of_materials=lmats,
id=tmp[tokens-1]) id=tmp[tokens-1])
scene_objects.append(material_list) scene_objects.append(material_list)
return scene_objects return scene_objects

24
gprMax/hash_cmds_singleuse.py
查看文件

@@ -55,8 +55,8 @@ def process_singlecmds(singlecmds):
if singlecmds[cmd] is not None: if singlecmds[cmd] is not None:
tmp = tuple(int(x) for x in singlecmds[cmd].split()) tmp = tuple(int(x) for x in singlecmds[cmd].split())
if len(tmp) != 1: if len(tmp) != 1:
logger.exception(cmd + ' requires exactly one parameter to specify ' + logger.exception(f'{cmd} requires exactly one parameter to specify ' +
'the number of CPU OpenMP threads to use') f'the number of CPU OpenMP threads to use')
raise ValueError raise ValueError
omp_threads = OMPThreads(n=tmp[0]) omp_threads = OMPThreads(n=tmp[0])
@@ -66,7 +66,7 @@ def process_singlecmds(singlecmds):
if singlecmds[cmd] is not None: if singlecmds[cmd] is not None:
tmp = [float(x) for x in singlecmds[cmd].split()] tmp = [float(x) for x in singlecmds[cmd].split()]
if len(tmp) != 3: if len(tmp) != 3:
logger.exception(cmd + ' requires exactly three parameters') logger.exception(f'{cmd} requires exactly three parameters')
raise ValueError raise ValueError
dl = (tmp[0], tmp[1], tmp[2]) dl = (tmp[0], tmp[1], tmp[2])
@@ -77,7 +77,7 @@ def process_singlecmds(singlecmds):
if singlecmds[cmd] is not None: if singlecmds[cmd] is not None:
tmp = [float(x) for x in singlecmds[cmd].split()] tmp = [float(x) for x in singlecmds[cmd].split()]
if len(tmp) != 3: if len(tmp) != 3:
logger.exception(cmd + ' requires exactly three parameters') logger.exception(f'{cmd} requires exactly three parameters')
raise ValueError raise ValueError
p1 = (tmp[0], tmp[1], tmp[2]) p1 = (tmp[0], tmp[1], tmp[2])
@@ -94,8 +94,8 @@ def process_singlecmds(singlecmds):
if singlecmds[cmd] is not None: if singlecmds[cmd] is not None:
tmp = singlecmds[cmd].split() tmp = singlecmds[cmd].split()
if len(tmp) != 1: if len(tmp) != 1:
logger.exception(cmd + ' requires exactly one parameter to specify ' + logger.exception(f'{cmd} requires exactly one parameter to specify the ' +
'the time window. Either in seconds or number of iterations.') f'time window. Either in seconds or number of iterations.')
raise ValueError raise ValueError
tmp = tmp[0].lower() tmp = tmp[0].lower()
@@ -113,8 +113,8 @@ def process_singlecmds(singlecmds):
cmd = '#pml_cells' cmd = '#pml_cells'
if singlecmds[cmd] is not None: if singlecmds[cmd] is not None:
tmp = singlecmds[cmd].split() tmp = singlecmds[cmd].split()
if len(tmp) != 1 and len(tmp) != 6: if len(tmp) not in [1, 6]:
logger.exception(cmd + ' requires either one or six parameter(s)') logger.exception(f'{cmd} requires either one or six parameter(s)')
raise ValueError raise ValueError
if len(tmp) == 1: if len(tmp) == 1:
pml_cells = PMLProps(thickness=int(tmp[0])) pml_cells = PMLProps(thickness=int(tmp[0]))
@@ -132,7 +132,7 @@ def process_singlecmds(singlecmds):
if singlecmds[cmd] is not None: if singlecmds[cmd] is not None:
tmp = singlecmds[cmd].split() tmp = singlecmds[cmd].split()
if len(tmp) != 3: if len(tmp) != 3:
logger.exception(cmd + ' requires exactly three parameters') logger.exception(f'{cmd} requires exactly three parameters')
raise ValueError raise ValueError
p1 = (float(tmp[0]), float(tmp[1]), float(tmp[2])) p1 = (float(tmp[0]), float(tmp[1]), float(tmp[2]))
@@ -143,7 +143,7 @@ def process_singlecmds(singlecmds):
if singlecmds[cmd] is not None: if singlecmds[cmd] is not None:
tmp = singlecmds[cmd].split() tmp = singlecmds[cmd].split()
if len(tmp) != 3: if len(tmp) != 3:
logger.exception(cmd + ' requires exactly three parameters') logger.exception(f'{cmd} requires exactly three parameters')
raise ValueError raise ValueError
p1 = (float(tmp[0]), float(tmp[1]), float(tmp[2])) p1 = (float(tmp[0]), float(tmp[1]), float(tmp[2]))
@@ -154,8 +154,8 @@ def process_singlecmds(singlecmds):
cmd = '#excitation_file' cmd = '#excitation_file'
if singlecmds[cmd] is not None: if singlecmds[cmd] is not None:
tmp = singlecmds[cmd].split() tmp = singlecmds[cmd].split()
if len(tmp) != 1 and len(tmp) != 3: if len(tmp) not in [1, 3]:
logger.exception(cmd + ' requires either one or three parameter(s)') logger.exception(f'{cmd} requires either one or three parameter(s)')
raise ValueError raise ValueError
if len(tmp) > 1: if len(tmp) > 1:

192
gprMax/materials.py
查看文件

@@ -215,77 +215,6 @@ class DispersiveMaterial(Material):
return er return er
def process_materials(G):
"""Processes complete list of materials - calculates update coefficients,
stores in arrays, and builds text list of materials/properties
Args:
G: FDTDGrid class describing a grid in a model.
Returns:
materialsdata: list of material IDs, names, and properties to
print a table.
"""
if config.get_model_config().materials['maxpoles'] == 0:
materialsdata = [['\nID', '\nName', '\nType', '\neps_r', 'sigma\n[S/m]',
'\nmu_r', 'sigma*\n[Ohm/m]', 'Dielectric\nsmoothable']]
else:
materialsdata = [['\nID', '\nName', '\nType', '\neps_r', 'sigma\n[S/m]',
'Delta\neps_r', 'tau\n[s]', 'omega\n[Hz]', 'delta\n[Hz]',
'gamma\n[Hz]', '\nmu_r', 'sigma*\n[Ohm/m]', 'Dielectric\nsmoothable']]
for material in G.materials:
# Calculate update coefficients for specific material
material.calculate_update_coeffsE(G)
material.calculate_update_coeffsH(G)
# Add update coefficients to overall storage for all materials
G.updatecoeffsE[material.numID, :] = material.CA, material.CBx, material.CBy, material.CBz, material.srce
G.updatecoeffsH[material.numID, :] = material.DA, material.DBx, material.DBy, material.DBz, material.srcm
# Add update coefficients to overall storage for dispersive materials
if hasattr(material, 'poles'):
z = 0
for pole in range(config.get_model_config().materials['maxpoles']):
G.updatecoeffsdispersive[material.numID, z:z + 3] = (config.sim_config.em_consts['e0'] *
material.eqt2[pole], material.eqt[pole], material.zt[pole])
z += 3
# Construct information on material properties for printing table
materialtext = []
materialtext.append(str(material.numID))
materialtext.append(material.ID[:50] if len(material.ID) > 50 else material.ID)
materialtext.append(material.type)
materialtext.append(f'{material.er:g}')
materialtext.append(f'{material.se:g}')
if config.get_model_config().materials['maxpoles'] > 0:
if 'debye' in material.type:
materialtext.append('\n'.join('{:g}'.format(deltaer) for deltaer in material.deltaer))
materialtext.append('\n'.join('{:g}'.format(tau) for tau in material.tau))
materialtext.extend(['', '', ''])
elif 'lorentz' in material.type:
materialtext.append(', '.join('{:g}'.format(deltaer) for deltaer in material.deltaer))
materialtext.append('')
materialtext.append(', '.join('{:g}'.format(tau) for tau in material.tau))
materialtext.append(', '.join('{:g}'.format(alpha) for alpha in material.alpha))
materialtext.append('')
elif 'drude' in material.type:
materialtext.extend(['', ''])
materialtext.append(', '.join('{:g}'.format(tau) for tau in material.tau))
materialtext.append('')
materialtext.append(', '.join('{:g}'.format(alpha) for alpha in material.alpha))
else:
materialtext.extend(['', '', '', '', ''])
materialtext.append(f'{material.mr:g}')
materialtext.append(f'{material.sm:g}')
materialtext.append(material.averagable)
materialsdata.append(materialtext)
return materialsdata
class PeplinskiSoil: class PeplinskiSoil:
"""Soil objects that are characterised according to a mixing model """Soil objects that are characterised according to a mixing model
by Peplinski (http://dx.doi.org/10.1109/36.387598). by Peplinski (http://dx.doi.org/10.1109/36.387598).
@@ -356,7 +285,6 @@ class PeplinskiSoil:
#mumaterials = mubins + (mubins[1] - mubins[0]) / 2 #mumaterials = mubins + (mubins[1] - mubins[0]) / 2
mumaterials = 0.5 * (mubins[1:nbins+1] + mubins[0:nbins]) mumaterials = 0.5 * (mubins[1:nbins+1] + mubins[0:nbins])
# Create an iterator # Create an iterator
muiter = np.nditer(mumaterials, flags=['c_index']) muiter = np.nditer(mumaterials, flags=['c_index'])
while not muiter.finished: while not muiter.finished:
@@ -399,7 +327,6 @@ class PeplinskiSoil:
muiter.iternext() muiter.iternext()
class RangeMaterial: class RangeMaterial:
"""Material objects defined by a given range of their parameters to be used for """Material objects defined by a given range of their parameters to be used for
factal spatial disttibutions. factal spatial disttibutions.
@@ -425,7 +352,6 @@ class RangeMaterial:
# and assume that all must be sequentially numbered. This allows for more general mixing models # and assume that all must be sequentially numbered. This allows for more general mixing models
self.matID = [] self.matID = []
def calculate_properties(self, nbins, G): def calculate_properties(self, nbins, G):
"""Calculates the specific properties of each of the materials. """Calculates the specific properties of each of the materials.
@@ -436,6 +362,7 @@ class RangeMaterial:
# Generate a set of relative permittivity bins based on the given range # Generate a set of relative permittivity bins based on the given range
erbins = np.linspace(self.er[0], self.er[1], nbins+1) erbins = np.linspace(self.er[0], self.er[1], nbins+1)
# Generate a range of relative permittivity values the mid-point of # Generate a range of relative permittivity values the mid-point of
# each bin to make materials from # each bin to make materials from
#ermaterials = erbins + np.abs((erbins[1] - erbins[0])) / 2 #ermaterials = erbins + np.abs((erbins[1] - erbins[0])) / 2
@@ -443,6 +370,7 @@ class RangeMaterial:
# Generate a set of conductivity bins based on the given range # Generate a set of conductivity bins based on the given range
sigmabins = np.linspace(self.sig[0], self.sig[1], nbins + 1) sigmabins = np.linspace(self.sig[0], self.sig[1], nbins + 1)
# Generate a range of conductivity values the mid-point of # Generate a range of conductivity values the mid-point of
# each bin to make materials from # each bin to make materials from
#sigmamaterials = sigmabins + (sigmabins[1] - sigmabins[0]) / 2 #sigmamaterials = sigmabins + (sigmabins[1] - sigmabins[0]) / 2
@@ -450,6 +378,7 @@ class RangeMaterial:
# Generate a set of magnetic permeability bins based on the given range # Generate a set of magnetic permeability bins based on the given range
mubins = np.linspace(self.mu[0], self.mu[1], nbins + 1) mubins = np.linspace(self.mu[0], self.mu[1], nbins + 1)
# Generate a range of magnetic permeability values the mid-point of # Generate a range of magnetic permeability values the mid-point of
# each bin to make materials from # each bin to make materials from
#mumaterials = mubins + np.abs((mubins[1] - mubins[0])) / 2 #mumaterials = mubins + np.abs((mubins[1] - mubins[0])) / 2
@@ -457,29 +386,26 @@ class RangeMaterial:
# Generate a set of magnetic loss bins based on the given range # Generate a set of magnetic loss bins based on the given range
robins = np.linspace(self.ro[0], self.ro[1], nbins + 1) robins = np.linspace(self.ro[0], self.ro[1], nbins + 1)
# Generate a range of magnetic loss values the mid-point of
# each bin to make materials from # Generate a range of magnetic loss values the mid-point of each bin to
# make materials from
#romaterials = robins + np.abs((robins[1] - robins[0])) / 2 #romaterials = robins + np.abs((robins[1] - robins[0])) / 2
romaterials = 0.5 * (robins[1:nbins+1] + robins[0:nbins]) romaterials = 0.5 * (robins[1:nbins+1] + robins[0:nbins])
# Iterate over the bins # Iterate over the bins
for iter in np.arange(0,nbins): for iter in np.arange(nbins):
# Relative permittivity # Relative permittivity
er = ermaterials[iter] er = ermaterials[iter]
# Effective conductivity # Effective conductivity
se = sigmamaterials[iter] se = sigmamaterials[iter]
# Magnetic permeability
# magnetic permeability
mr = mumaterials[iter] mr = mumaterials[iter]
# Magnetic loss
# magnetic loss
sm = romaterials[iter] sm = romaterials[iter]
# Check to see if the material already exists before creating a new one # Check to see if the material already exists before creating a new one
requiredID = '|{:.4f}+{:.4f}+{:.4f}+{:.4f}|'.format(float(er),float(se),float(mr),float(sm)) requiredID = f'|{float(er):.4f}+{float(se):.4f}+{float(mr):.4f}+{float(sm):.4f}|'
material = next((x for x in G.materials if x.ID == requiredID), None) material = next((x for x in G.materials if x.ID == requiredID), None)
if iter == 0: if iter == 0:
if material: if material:
@@ -499,11 +425,10 @@ class RangeMaterial:
self.matID.append(m.numID) self.matID.append(m.numID)
class ListMaterial: class ListMaterial:
"""A list of predefined materials that are to be used for """A list of predefined materials to be used for
factal spatial disttibutions. No new materials are created but the ones specified are factal spatial disttibutions. This command does not create new materials but collects them to be used in a
grouped together to be used in fractal spatial distributions. stochastic distribution by a fractal box.
""" """
def __init__(self, ID, listofmaterials): def __init__(self, ID, listofmaterials):
@@ -524,21 +449,18 @@ class ListMaterial:
def calculate_properties(self, nbins, G): def calculate_properties(self, nbins, G):
"""Calculates the properties of the materials. No Debye is used but name kept the same as used in other """Calculates the properties of the materials.
class that needs Debye
Args: Args:
nbins: int for number of bins to use to create the different materials. nbins: int for number of bins to use to create the different materials.
G: FDTDGrid class describing a grid in a model. G: FDTDGrid class describing a grid in a model.
""" """
# Iterate over the bins # Iterate over the bins
for iter in np.arange(0,nbins): for iter in np.arange(nbins):
# Check to see if the material already exists before creating a new one
#requiredID = '|{:}_in_{:}|'.format((self.mat[iter]),(self.ID)) #requiredID = '|{:}_in_{:}|'.format((self.mat[iter]),(self.ID))
requiredID = self.mat[iter] requiredID = self.mat[iter]
# Check if the material already exists before creating a new one
material = next((x for x in G.materials if x.ID == requiredID), None) material = next((x for x in G.materials if x.ID == requiredID), None)
self.matID.append(material.numID) self.matID.append(material.numID)
@@ -555,15 +477,11 @@ class ListMaterial:
# m.numID = len(G.materials) # m.numID = len(G.materials)
# G.materials.append(m) # G.materials.append(m)
if not material: if not material:
logger.exception(self.__str__() + f' material(s) {material} do not exist') logger.exception(self.__str__() + f' material(s) {material} do not exist')
raise ValueError raise ValueError
def create_built_in_materials(G): def create_built_in_materials(G):
"""Creates pre-defined (built-in) materials. """Creates pre-defined (built-in) materials.
@@ -571,8 +489,6 @@ def create_built_in_materials(G):
G: FDTDGrid class describing a grid in a model. G: FDTDGrid class describing a grid in a model.
""" """
G.n_built_in_materials = len(G.materials)
m = Material(0, 'pec') m = Material(0, 'pec')
m.se = float('inf') m.se = float('inf')
m.type = 'builtin' m.type = 'builtin'
@@ -583,8 +499,6 @@ def create_built_in_materials(G):
m.type = 'builtin' m.type = 'builtin'
G.materials.append(m) G.materials.append(m)
G.n_built_in_materials = len(G.materials)
def calculate_water_properties(T=25, S=0): def calculate_water_properties(T=25, S=0):
"""Get extended Debye model properties for water. """Get extended Debye model properties for water.
@@ -627,8 +541,6 @@ def create_water(G, T=25, S=0):
eri, er, tau, sig = calculate_water_properties(T, S) eri, er, tau, sig = calculate_water_properties(T, S)
G.n_built_in_materials = len(G.materials)
m = DispersiveMaterial(len(G.materials), 'water') m = DispersiveMaterial(len(G.materials), 'water')
m.averagable = False m.averagable = False
m.type = 'builtin, debye' m.type = 'builtin, debye'
@@ -641,8 +553,6 @@ def create_water(G, T=25, S=0):
if config.get_model_config().materials['maxpoles'] == 0: if config.get_model_config().materials['maxpoles'] == 0:
config.get_model_config().materials['maxpoles'] = 1 config.get_model_config().materials['maxpoles'] = 1
G.n_built_in_materials = len(G.materials)
def create_grass(G): def create_grass(G):
"""Creates single-pole Debye model for grass """Creates single-pole Debye model for grass
@@ -657,8 +567,6 @@ def create_grass(G):
tau = 1.0793e-11 tau = 1.0793e-11
sig = 0 sig = 0
G.n_built_in_materials = len(G.materials)
m = DispersiveMaterial(len(G.materials), 'grass') m = DispersiveMaterial(len(G.materials), 'grass')
m.averagable = False m.averagable = False
m.type = 'builtin, debye' m.type = 'builtin, debye'
@@ -671,4 +579,70 @@ def create_grass(G):
if config.get_model_config().materials['maxpoles'] == 0: if config.get_model_config().materials['maxpoles'] == 0:
config.get_model_config().materials['maxpoles'] = 1 config.get_model_config().materials['maxpoles'] = 1
G.n_built_in_materials = len(G.materials)
def process_materials(G):
"""Processes complete list of materials - calculates update coefficients,
stores in arrays, and builds text list of materials/properties
Args:
G: FDTDGrid class describing a grid in a model.
Returns:
materialsdata: list of material IDs, names, and properties to
print a table.
"""
if config.get_model_config().materials['maxpoles'] == 0:
materialsdata = [['\nID', '\nName', '\nType', '\neps_r', 'sigma\n[S/m]',
'\nmu_r', 'sigma*\n[Ohm/m]', 'Dielectric\nsmoothable']]
else:
materialsdata = [['\nID', '\nName', '\nType', '\neps_r', 'sigma\n[S/m]',
'Delta\neps_r', 'tau\n[s]', 'omega\n[Hz]', 'delta\n[Hz]',
'gamma\n[Hz]', '\nmu_r', 'sigma*\n[Ohm/m]', 'Dielectric\nsmoothable']]
for material in G.materials:
# Calculate update coefficients for specific material
material.calculate_update_coeffsE(G)
material.calculate_update_coeffsH(G)
# Add update coefficients to overall storage for all materials
G.updatecoeffsE[material.numID, :] = material.CA, material.CBx, material.CBy, material.CBz, material.srce
G.updatecoeffsH[material.numID, :] = material.DA, material.DBx, material.DBy, material.DBz, material.srcm
# Add update coefficients to overall storage for dispersive materials
if hasattr(material, 'poles'):
z = 0
for pole in range(config.get_model_config().materials['maxpoles']):
G.updatecoeffsdispersive[material.numID, z:z + 3] = (config.sim_config.em_consts['e0'] *
material.eqt2[pole], material.eqt[pole], material.zt[pole])
z += 3
# Construct information on material properties for printing table
materialtext = [str(material.numID),
material.ID[:50] if len(material.ID) > 50 else material.ID,
material.type,
f'{material.er:g}',
f'{material.se:g}']
if config.get_model_config().materials['maxpoles'] > 0:
if 'debye' in material.type:
materialtext.append('\n'.join('{:g}'.format(deltaer) for deltaer in material.deltaer))
materialtext.append('\n'.join('{:g}'.format(tau) for tau in material.tau))
materialtext.extend(['', '', ''])
elif 'lorentz' in material.type:
materialtext.append(', '.join('{:g}'.format(deltaer) for deltaer in material.deltaer))
materialtext.append('')
materialtext.append(', '.join('{:g}'.format(tau) for tau in material.tau))
materialtext.append(', '.join('{:g}'.format(alpha) for alpha in material.alpha))
materialtext.append('')
elif 'drude' in material.type:
materialtext.extend(['', ''])
materialtext.append(', '.join('{:g}'.format(tau) for tau in material.tau))
materialtext.append('')
materialtext.append(', '.join('{:g}'.format(alpha) for alpha in material.alpha))
else:
materialtext.extend(['', '', '', '', ''])
materialtext.extend((f'{material.mr:g}', f'{material.sm:g}', material.averagable))
materialsdata.append(materialtext)
return materialsdata

13
gprMax/model_build_run.py
查看文件

@@ -73,7 +73,7 @@ class ModelBuildRun:
self.p = psutil.Process() self.p = psutil.Process()
# Normal model reading/building process; bypassed if geometry information to be reused # Normal model reading/building process; bypassed if geometry information to be reused
self.build_geometry() if not config.get_model_config().reuse_geometry else self.reuse_geometry() self.reuse_geometry() if config.get_model_config().reuse_geometry else self.build_geometry()
logger.info(f'\nOutput directory: {config.get_model_config().output_file_path.parent.resolve()}') logger.info(f'\nOutput directory: {config.get_model_config().output_file_path.parent.resolve()}')
@@ -109,7 +109,7 @@ class ModelBuildRun:
# Write files for any geometry views and geometry object outputs # Write files for any geometry views and geometry object outputs
gvs = G.geometryviews + [gv for sg in G.subgrids for gv in sg.geometryviews] gvs = G.geometryviews + [gv for sg in G.subgrids for gv in sg.geometryviews]
if not (gvs or G.geometryobjectswrite) and config.sim_config.args.geometry_only: if (not gvs and not G.geometryobjectswrite and config.sim_config.args.geometry_only):
logger.exception('\nNo geometry views or geometry objects found.') logger.exception('\nNo geometry views or geometry objects found.')
raise ValueError raise ValueError
save_geometry_views(gvs) save_geometry_views(gvs)
@@ -271,8 +271,7 @@ class ModelBuildRun:
logger.warning(f"You have specified more threads ({config.get_model_config().ompthreads}) " logger.warning(f"You have specified more threads ({config.get_model_config().ompthreads}) "
f"than available physical CPU cores ({config.sim_config.hostinfo['physicalcores']}). " f"than available physical CPU cores ({config.sim_config.hostinfo['physicalcores']}). "
f"This may lead to degraded performance.") f"This may lead to degraded performance.")
# Print information about any compute device, e.g. GPU, in use elif config.sim_config.general['solver'] in ['cuda', 'opencl']:
elif config.sim_config.general['solver'] == 'cuda' or config.sim_config.general['solver'] == 'opencl':
if config.sim_config.general['solver'] == 'opencl': if config.sim_config.general['solver'] == 'opencl':
solvername = 'OpenCL' solvername = 'OpenCL'
platformname = ' on ' + ' '.join(config.get_model_config().device['dev'].platform.name.split()) platformname = ' on ' + ' '.join(config.get_model_config().device['dev'].platform.name.split())
@@ -338,11 +337,11 @@ class GridBuilder:
pbar.close() pbar.close()
def tm_grid_update(self): def tm_grid_update(self):
if '2D TMx' == config.get_model_config().mode: if config.get_model_config().mode == '2D TMx':
self.grid.tmx() self.grid.tmx()
elif '2D TMy' == config.get_model_config().mode: elif config.get_model_config().mode == '2D TMy':
self.grid.tmy() self.grid.tmy()
elif '2D TMz' == config.get_model_config().mode: elif config.get_model_config().mode == '2D TMz':
self.grid.tmz() self.grid.tmz()
def update_voltage_source_materials(self): def update_voltage_source_materials(self):

4
gprMax/mpi.py
查看文件

@@ -220,8 +220,8 @@ class MPIExecutor(object):
def join(self): def join(self):
"""Joins the workers.""" """Joins the workers."""
if self.is_master(): if not self.is_master():
return
logger.debug(f'({self.comm.name}) - Terminating. Sending sentinel to all workers.') logger.debug(f'({self.comm.name}) - Terminating. Sending sentinel to all workers.')
# Send sentinel to all workers # Send sentinel to all workers
for worker in self.workers: for worker in self.workers:

32
gprMax/pml.py
查看文件

@@ -22,6 +22,8 @@ import numpy as np
import gprMax.config as config import gprMax.config as config
from .cython.pml_build import pml_average_er_mr
class CFSParameter: class CFSParameter:
"""Individual CFS parameter (e.g. alpha, kappa, or sigma).""" """Individual CFS parameter (e.g. alpha, kappa, or sigma)."""
@@ -609,9 +611,6 @@ def build_pml(G, key, value):
elif config.sim_config.general['solver'] == 'opencl': elif config.sim_config.general['solver'] == 'opencl':
pml_type = OpenCLPML pml_type = OpenCLPML
sumer = 0 # Sum of relative permittivities in PML slab
summr = 0 # Sum of relative permeabilities in PML slab
if key[0] == 'x': if key[0] == 'x':
if key == 'x0': if key == 'x0':
pml = pml_type(G, ID=key, direction='xminus', pml = pml_type(G, ID=key, direction='xminus',
@@ -621,14 +620,7 @@ def build_pml(G, key, value):
xs=G.nx - value, xf=G.nx, yf=G.ny, zf=G.nz) xs=G.nx - value, xf=G.nx, yf=G.ny, zf=G.nz)
pml.CFS = G.pmls['cfs'] pml.CFS = G.pmls['cfs']
G.pmls['slabs'].append(pml) G.pmls['slabs'].append(pml)
for j in range(G.ny): averageer, averagemr = pml_average_er_mr('x', pml.xs, G)
for k in range(G.nz):
numID = G.solid[pml.xs, j, k]
material = next(x for x in G.materials if x.numID == numID)
sumer += material.er
summr += material.mr
averageer = sumer / (G.ny * G.nz)
averagemr = summr / (G.ny * G.nz)
elif key[0] == 'y': elif key[0] == 'y':
if key == 'y0': if key == 'y0':
@@ -639,14 +631,7 @@ def build_pml(G, key, value):
ys=G.ny - value, xf=G.nx, yf=G.ny, zf=G.nz) ys=G.ny - value, xf=G.nx, yf=G.ny, zf=G.nz)
pml.CFS = G.pmls['cfs'] pml.CFS = G.pmls['cfs']
G.pmls['slabs'].append(pml) G.pmls['slabs'].append(pml)
for i in range(G.nx): averageer, averagemr = pml_average_er_mr('y', pml.ys, G)
for k in range(G.nz):
numID = G.solid[i, pml.ys, k]
material = next(x for x in G.materials if x.numID == numID)
sumer += material.er
summr += material.mr
averageer = sumer / (G.nx * G.nz)
averagemr = summr / (G.nx * G.nz)
elif key[0] == 'z': elif key[0] == 'z':
if key == 'z0': if key == 'z0':
@@ -657,13 +642,6 @@ def build_pml(G, key, value):
zs=G.nz - value, xf=G.nx, yf=G.ny, zf=G.nz) zs=G.nz - value, xf=G.nx, yf=G.ny, zf=G.nz)
pml.CFS = G.pmls['cfs'] pml.CFS = G.pmls['cfs']
G.pmls['slabs'].append(pml) G.pmls['slabs'].append(pml)
for i in range(G.nx): averageer, averagemr = pml_average_er_mr('z', pml.zs, G)
for j in range(G.ny):
numID = G.solid[i, j, pml.zs]
material = next(x for x in G.materials if x.numID == numID)
sumer += material.er
summr += material.mr
averageer = sumer / (G.nx * G.ny)
averagemr = summr / (G.nx * G.ny)
pml.calculate_update_coeffs(averageer, averagemr) pml.calculate_update_coeffs(averageer, averagemr)

2
gprMax/scene.py
查看文件

@@ -105,7 +105,7 @@ class Scene:
# Check essential commands and warn user if missing # Check essential commands and warn user if missing
for cmd_type in self.essential_cmds: for cmd_type in self.essential_cmds:
d = any([isinstance(cmd, cmd_type) for cmd in cmds_unique]) d = any(isinstance(cmd, cmd_type) for cmd in cmds_unique)
if not d: if not d:
logger.exception('Your input file is missing essential commands ' + logger.exception('Your input file is missing essential commands ' +
'required to run a model. Essential commands ' + 'required to run a model. Essential commands ' +

13
gprMax/sources.py
查看文件

@@ -96,7 +96,7 @@ class VoltageSource(Source):
i = self.xcoord i = self.xcoord
j = self.ycoord j = self.ycoord
k = self.zcoord k = self.zcoord
componentID = 'E' + self.polarisation componentID = f'E{self.polarisation}'
if self.polarisation == 'x': if self.polarisation == 'x':
if self.resistance != 0: if self.resistance != 0:
@@ -130,16 +130,17 @@ class VoltageSource(Source):
G: FDTDGrid class describing a grid in a model. G: FDTDGrid class describing a grid in a model.
""" """
if self.resistance != 0: if self.resistance == 0:
return
i = self.xcoord i = self.xcoord
j = self.ycoord j = self.ycoord
k = self.zcoord k = self.zcoord
componentID = 'E' + self.polarisation componentID = f'E{self.polarisation}'
requirednumID = G.ID[G.IDlookup[componentID], i, j, k] requirednumID = G.ID[G.IDlookup[componentID], i, j, k]
material = next(x for x in G.materials if x.numID == requirednumID) material = next(x for x in G.materials if x.numID == requirednumID)
newmaterial = deepcopy(material) newmaterial = deepcopy(material)
newmaterial.ID = material.ID + '+' + self.ID newmaterial.ID = f'{material.ID}+{self.ID}'
newmaterial.numID = len(G.materials) newmaterial.numID = len(G.materials)
newmaterial.averagable = False newmaterial.averagable = False
newmaterial.type += ',\nvoltage-source' if newmaterial.type else 'voltage-source' newmaterial.type += ',\nvoltage-source' if newmaterial.type else 'voltage-source'
@@ -180,7 +181,7 @@ class HertzianDipole(Source):
i = self.xcoord i = self.xcoord
j = self.ycoord j = self.ycoord
k = self.zcoord k = self.zcoord
componentID = 'E' + self.polarisation componentID = f'E{self.polarisation}'
if self.polarisation == 'x': if self.polarisation == 'x':
Ex[i, j, k] -= (updatecoeffsE[ID[G.IDlookup[componentID], i, j, k], 4] * Ex[i, j, k] -= (updatecoeffsE[ID[G.IDlookup[componentID], i, j, k], 4] *
self.waveformvalues_halfdt[iteration] * self.dl * self.waveformvalues_halfdt[iteration] * self.dl *
@@ -217,7 +218,7 @@ class MagneticDipole(Source):
i = self.xcoord i = self.xcoord
j = self.ycoord j = self.ycoord
k = self.zcoord k = self.zcoord
componentID = 'H' + self.polarisation componentID = f'H{self.polarisation}'
if self.polarisation == 'x': if self.polarisation == 'x':
Hx[i, j, k] -= (updatecoeffsH[ID[G.IDlookup[componentID], i, j, k], 4] * Hx[i, j, k] -= (updatecoeffsH[ID[G.IDlookup[componentID], i, j, k], 4] *

8
gprMax/subgrids/precursor_nodes.py
查看文件

@@ -174,7 +174,7 @@ class PrecursorNodesBase:
for f in field_names: for f in field_names:
try: try:
val = c1 * getattr(self, f + '_0') + c2 * getattr(self, f + '_1') val = c1 * getattr(self, f'{f}_0') + c2 * getattr(self, f'{f}_1')
except ValueError: except ValueError:
raise raise
setattr(self, f, val) setattr(self, f, val)
@@ -184,7 +184,7 @@ class PrecursorNodesBase:
current main time step, i.e. ey_left = copy.ey_left_1 current main time step, i.e. ey_left = copy.ey_left_1
""" """
for f in field_names: for f in field_names:
val = np.copy(getattr(self, f + '_1')) val = np.copy(getattr(self, f'{f}_1'))
setattr(self, f, val) setattr(self, f, val)
def calc_exact_magnetic_in_time(self): def calc_exact_magnetic_in_time(self):
@@ -218,9 +218,9 @@ class PrecursorNodesBase:
def update_previous_timestep_fields(self, field_names): def update_previous_timestep_fields(self, field_names):
for fn in field_names: for fn in field_names:
val = getattr(self, fn + '_1') val = getattr(self, f'{fn}_1')
val_c = np.copy(val) val_c = np.copy(val)
setattr(self, fn + '_0', val_c) setattr(self, f'{fn}_0', val_c)
def interpolate_to_sub_grid(self, field, coords): def interpolate_to_sub_grid(self, field, coords):
x, z, x_sg, z_sg = coords x, z, x_sg, z_sg = coords

4
gprMax/subgrids/updates.py
查看文件

@@ -33,10 +33,10 @@ def create_updates(G):
sg_type = type(sg) sg_type = type(sg)
if sg_type == SubGridHSG and sg.filter: if sg_type == SubGridHSG and sg.filter:
precursors = PrecursorNodesFiltered(G, sg) precursors = PrecursorNodesFiltered(G, sg)
elif sg_type == SubGridHSG and not sg.filter: elif sg_type == SubGridHSG:
precursors = PrecursorNodes(G, sg) precursors = PrecursorNodes(G, sg)
else: else:
logger.exception(str(sg) + ' is not a subgrid type') logger.exception(f'{str(sg)} is not a subgrid type')
raise ValueError raise ValueError
sgu = SubgridUpdater(sg, precursors, G) sgu = SubgridUpdater(sg, precursors, G)

6
gprMax/subgrids/user_objects.py
查看文件

@@ -45,7 +45,7 @@ class SubGridBase(UserObjectMulti):
elif isinstance(node, UserObjectGeometry): elif isinstance(node, UserObjectGeometry):
self.children_geometry.append(node) self.children_geometry.append(node)
else: else:
logger.exception(str(node) + ' this Object can not be added to a sub grid') logger.exception(f'{str(node)} this Object can not be added to a sub grid')
raise ValueError raise ValueError
def set_discretisation(self, sg, grid): def set_discretisation(self, sg, grid):
@@ -123,12 +123,12 @@ class SubGridBase(UserObjectMulti):
self.subgrid = sg self.subgrid = sg
# Copy over built in materials # Copy over built in materials
sg.materials = [copy(m) for m in grid.materials if m.numID in range(0, grid.n_built_in_materials + 1)] sg.materials = [copy(m) for m in grid.materials if m.type == 'builtin']
# Don't mix and match different subgrid types # Don't mix and match different subgrid types
for sg_made in grid.subgrids: for sg_made in grid.subgrids:
if type(sg) != type(sg_made): if type(sg) != type(sg_made):
logger.exception(self.__str__() + ' please only use one type of subgrid') logger.exception(f'{self.__str__()} please only use one type of subgrid')
raise ValueError raise ValueError
# Reference the subgrid under the main grid to which it belongs # Reference the subgrid under the main grid to which it belongs

10
gprMax/updates.py
查看文件

@@ -378,7 +378,7 @@ class CUDAUpdates:
for pml in self.grid.pmls['slabs']: for pml in self.grid.pmls['slabs']:
pml.htod_field_arrays() pml.htod_field_arrays()
pml.set_blocks_per_grid() pml.set_blocks_per_grid()
knl_name = 'order' + str(len(pml.CFS)) + '_' + pml.direction knl_name = f'order{len(pml.CFS)}_{pml.direction}'
self.subs_name_args['FUNC'] = knl_name self.subs_name_args['FUNC'] = knl_name
knl_electric = getattr(knl_pml_updates_electric, knl_name) knl_electric = getattr(knl_pml_updates_electric, knl_name)
@@ -880,21 +880,21 @@ class OpenCLUpdates:
for pml in self.grid.pmls['slabs']: for pml in self.grid.pmls['slabs']:
pml.set_queue(self.queue) pml.set_queue(self.queue)
pml.htod_field_arrays() pml.htod_field_arrays()
knl_name = 'order' + str(len(pml.CFS)) + '_' + pml.direction knl_name = f'order{len(pml.CFS)}_{pml.direction}'
knl_electric_name = getattr(knl_pml_updates_electric, knl_name) knl_electric_name = getattr(knl_pml_updates_electric, knl_name)
knl_magnetic_name = getattr(knl_pml_updates_magnetic, knl_name) knl_magnetic_name = getattr(knl_pml_updates_magnetic, knl_name)
pml.update_electric_dev = self.elwise(self.ctx, pml.update_electric_dev = self.elwise(self.ctx,
knl_electric_name['args_opencl'].substitute({'REAL': config.sim_config.dtypes['C_float_or_double']}), knl_electric_name['args_opencl'].substitute({'REAL': config.sim_config.dtypes['C_float_or_double']}),
knl_electric_name['func'].substitute(subs), knl_electric_name['func'].substitute(subs),
'pml_updates_electric_' + knl_name, f'pml_updates_electric_{knl_name}',
preamble=self.knl_common, preamble=self.knl_common,
options=config.sim_config.devices['compiler_opts']) options=config.sim_config.devices['compiler_opts'],)
pml.update_magnetic_dev = self.elwise(self.ctx, pml.update_magnetic_dev = self.elwise(self.ctx,
knl_magnetic_name['args_opencl'].substitute({'REAL': config.sim_config.dtypes['C_float_or_double']}), knl_magnetic_name['args_opencl'].substitute({'REAL': config.sim_config.dtypes['C_float_or_double']}),
knl_magnetic_name['func'].substitute(subs), knl_magnetic_name['func'].substitute(subs),
'pml_updates_magnetic_' + knl_name, f'pml_updates_magnetic_{knl_name}',
preamble=self.knl_common, preamble=self.knl_common,
options=config.sim_config.devices['compiler_opts']) options=config.sim_config.devices['compiler_opts'])

2
gprMax/user_inputs.py
查看文件

@@ -65,7 +65,7 @@ class UserInput:
except ValueError as err: except ValueError as err:
v = ['x', 'y', 'z'] v = ['x', 'y', 'z']
# Discretisation # Discretisation
dl = getattr(self.grid, 'd' + err.args[0]) dl = getattr(self.grid, f'd{err.args[0]}')
# Incorrect index # Incorrect index
i = p[v.index(err.args[0])] i = p[v.index(err.args[0])]
if name: if name:

2
gprMax/utilities/host_info.py
查看文件

@@ -515,7 +515,7 @@ def print_opencl_info(devs):
if i == 0: if i == 0:
platform = dev.platform.name platform = dev.platform.name
logger.basic(f' |--->Platform: {platform}') logger.basic(f' |--->Platform: {platform}')
if not platform == dev.platform.name: if platform != dev.platform.name:
logger.basic(f' |--->Platform: {dev.platform.name}') logger.basic(f' |--->Platform: {dev.platform.name}')
types = cl.device_type.to_string(dev.type) types = cl.device_type.to_string(dev.type)
if 'CPU' in types: if 'CPU' in types:

2
gprMax/utilities/utilities.py
查看文件

@@ -91,7 +91,7 @@ def logo(version):
v""" + version + '\n\n' v""" + version + '\n\n'
str = f"{description} {'=' * (get_terminal_width() - len(description) - 1)}\n\n" str = f"{description} {'=' * (get_terminal_width() - len(description) - 1)}\n\n"
str += Fore.CYAN + f'{logo}' str += f'{Fore.CYAN}{logo}'
str += Style.RESET_ALL + textwrap.fill(copyright, str += Style.RESET_ALL + textwrap.fill(copyright,
width=get_terminal_width() - 1, width=get_terminal_width() - 1,
initial_indent=' ') + '\n' initial_indent=' ') + '\n'

21
gprMax/waveforms.py
查看文件

@@ -54,18 +54,19 @@ class Waveform:
waveforms. waveforms.
""" """
if (self.type == 'gaussian' or self.type == 'gaussiandot' or if self.type in ['gaussian',
self.type == 'gaussiandotnorm' or self.type == 'gaussianprime' or 'gaussiandot',
self.type == 'gaussiandoubleprime'): 'gaussiandotnorm',
'gaussianprime',
'gaussiandoubleprime']:
self.chi = 1 / self.freq self.chi = 1 / self.freq
self.zeta = 2 * np.pi**2 * self.freq**2 self.zeta = 2 * np.pi**2 * self.freq**2
elif (self.type == 'gaussiandotdot' or elif self.type in ['gaussiandotdot', 'gaussiandotdotnorm', 'ricker']:
self.type == 'gaussiandotdotnorm' or self.type == 'ricker'):
self.chi = np.sqrt(2) / self.freq self.chi = np.sqrt(2) / self.freq
self.zeta = np.pi**2 * self.freq**2 self.zeta = np.pi**2 * self.freq**2
def calculate_value(self, time, dt): def calculate_value(self, time, dt):
"""Calculates value of the waveform at a specific time. """Calculates the value of the waveform at a specific time.
Args: Args:
time: float for absolute time. time: float for absolute time.
@@ -82,7 +83,7 @@ class Waveform:
delay = time - self.chi delay = time - self.chi
ampvalue = np.exp(-self.zeta * delay**2) ampvalue = np.exp(-self.zeta * delay**2)
elif self.type == 'gaussiandot' or self.type == 'gaussianprime': elif self.type in ['gaussiandot', 'gaussianprime']:
delay = time - self.chi delay = time - self.chi
ampvalue = -2 * self.zeta * delay * np.exp(-self.zeta * delay**2) ampvalue = -2 * self.zeta * delay * np.exp(-self.zeta * delay**2)
@@ -91,7 +92,7 @@ class Waveform:
normalise = np.sqrt(np.exp(1) / (2 * self.zeta)) normalise = np.sqrt(np.exp(1) / (2 * self.zeta))
ampvalue = -2 * self.zeta * delay * np.exp(-self.zeta * delay**2) * normalise ampvalue = -2 * self.zeta * delay * np.exp(-self.zeta * delay**2) * normalise
elif self.type == 'gaussiandotdot' or self.type == 'gaussiandoubleprime': elif self.type in ['gaussiandotdot', 'gaussiandoubleprime']:
delay = time - self.chi delay = time - self.chi
ampvalue = (2 * self.zeta * (2 * self.zeta * delay**2 - 1) * ampvalue = (2 * self.zeta * (2 * self.zeta * delay**2 - 1) *
np.exp(-self.zeta * delay**2)) np.exp(-self.zeta * delay**2))
@@ -116,8 +117,8 @@ class Waveform:
elif self.type == 'contsine': elif self.type == 'contsine':
rampamp = 0.25 rampamp = 0.25
ramp = rampamp * time * self.freq ramp = rampamp * time * self.freq
if ramp > 1: ramp = min(ramp, 1)
ramp = 1
ampvalue = ramp * np.sin(2 * np.pi * self.freq * time) ampvalue = ramp * np.sin(2 * np.pi * self.freq * time)
elif self.type == 'impulse': elif self.type == 'impulse':

38
setup.py
查看文件

@@ -37,8 +37,9 @@ if sys.version_info[:2] < MIN_PYTHON_VERSION:
# Importing gprMax _version__.py before building can cause issues. # Importing gprMax _version__.py before building can cause issues.
with open('gprMax/_version.py', 'r') as fd: with open('gprMax/_version.py', 'r') as fd:
version = re.search(r'^__version__\s*=\s*[\'"]([^\'"]*)[\'"]', version = re.search(
fd.read(), re.MULTILINE).group(1) r'^__version__\s*=\s*[\'"]([^\'"]*)[\'"]', fd.read(), re.MULTILINE
)[1]
def build_dispersive_material_templates(): def build_dispersive_material_templates():
"""Function to generate Cython .pyx file for dispersive media update. """Function to generate Cython .pyx file for dispersive media update.
@@ -46,7 +47,7 @@ def build_dispersive_material_templates():
functions. functions.
""" """
iswin = True if sys.platform == 'win32' else False iswin = (sys.platform == 'win32')
env = Environment(loader = FileSystemLoader(os.path.join('gprMax', 'cython')), ) env = Environment(loader = FileSystemLoader(os.path.join('gprMax', 'cython')), )
@@ -133,12 +134,12 @@ if 'cleanall' in sys.argv:
for file in cythonfiles: for file in cythonfiles:
filebase = os.path.splitext(file)[0] filebase = os.path.splitext(file)[0]
# Remove Cython C files # Remove Cython C files
if os.path.isfile(filebase + '.c'): if os.path.isfile(f'{filebase}.c'):
try: try:
os.remove(filebase + '.c') os.remove(f'{filebase}.c')
print(f'Removed: {filebase + ".c"}') print(f'Removed: {filebase}.c')
except OSError: except OSError:
print(f'Could not remove: {filebase + ".c"}') print(f'Could not remove: {filebase}.c')
# Remove compiled Cython modules # Remove compiled Cython modules
libfile = (glob.glob(os.path.join(os.getcwd(), libfile = (glob.glob(os.path.join(os.getcwd(),
os.path.splitext(file)[0]) + '*.pyd') + os.path.splitext(file)[0]) + '*.pyd') +
@@ -176,10 +177,6 @@ else:
linker_args = [] linker_args = []
libraries = [] libraries = []
# Compiler options - macOS - needs gcc (usually via HomeBrew) because the
# default compiler LLVM (clang) does not
# support OpenMP. With gcc -fopenmp option
# implies -pthread
elif sys.platform == 'darwin': elif sys.platform == 'darwin':
# Check for Intel or Apple M series CPU # Check for Intel or Apple M series CPU
cpuID = subprocess.check_output("sysctl -n machdep.cpu.brand_string", cpuID = subprocess.check_output("sysctl -n machdep.cpu.brand_string",
@@ -219,12 +216,14 @@ else:
pass pass
os.environ['MIN_SUPPORTED_MACOSX_DEPLOYMENT_TARGET'] = MIN_MACOS_VERSION os.environ['MIN_SUPPORTED_MACOSX_DEPLOYMENT_TARGET'] = MIN_MACOS_VERSION
# Sometimes worth testing with '-fstrict-aliasing', '-fno-common' # Sometimes worth testing with '-fstrict-aliasing', '-fno-common'
compile_args = ['-O3', '-w', '-fopenmp', '-march=native', compile_args = ['-O3',
'-mmacosx-version-min=' + MIN_MACOS_VERSION] '-w',
linker_args = ['-fopenmp', '-mmacosx-version-min=' + MIN_MACOS_VERSION] '-fopenmp',
'-march=native',
f'-mmacosx-version-min={MIN_MACOS_VERSION}']
linker_args = ['-fopenmp', f'-mmacosx-version-min={MIN_MACOS_VERSION}']
libraries=['gomp'] libraries=['gomp']
# Compiler options - Linux
elif sys.platform == 'linux': elif sys.platform == 'linux':
compile_args = ['-O3', '-w', '-fopenmp', '-march=native'] compile_args = ['-O3', '-w', '-fopenmp', '-march=native']
linker_args = ['-fopenmp'] linker_args = ['-fopenmp']
@@ -261,12 +260,12 @@ else:
version=version, version=version,
author='Craig Warren, Antonis Giannopoulos, and John Hartley', author='Craig Warren, Antonis Giannopoulos, and John Hartley',
url='http://www.gprmax.com', url='http://www.gprmax.com',
description='Electromagnetic Modelling Software based on the ' + description='Electromagnetic Modelling Software based on the '
'Finite-Difference Time-Domain (FDTD) method', + 'Finite-Difference Time-Domain (FDTD) method',
long_description=long_description, long_description=long_description,
long_description_content_type="text/x-rst", long_description_content_type="text/x-rst",
license='GPLv3+', license='GPLv3+',
python_requires='>' + str(MIN_PYTHON_VERSION[0]) + '.' + str(MIN_PYTHON_VERSION[1]), python_requires=f'>{str(MIN_PYTHON_VERSION[0])}.{str(MIN_PYTHON_VERSION[1])}',
install_requires=['colorama', install_requires=['colorama',
'cython', 'cython',
'h5py', 'h5py',
@@ -289,4 +288,5 @@ else:
'Operating System :: POSIX :: Linux', 'Operating System :: POSIX :: Linux',
'Programming Language :: Cython', 'Programming Language :: Cython',
'Programming Language :: Python :: 3', 'Programming Language :: Python :: 3',
'Topic :: Scientific/Engineering']) 'Topic :: Scientific/Engineering'],
)

4
testing/models_pmls/plot_pml_comparison.py
查看文件

@@ -93,8 +93,8 @@ for x, model in enumerate(testmodels):
# Diffs # Diffs
datadiffs = np.zeros(datatest.shape, dtype=np.float64) datadiffs = np.zeros(datatest.shape, dtype=np.float64)
for i in range(len(outputstest)): for i in range(len(outputstest)):
max = np.amax(np.abs(dataref[:, i])) maxi = np.amax(np.abs(dataref[:, i]))
datadiffs[:, i] = np.divide(np.abs(datatest[:, i] - dataref[:, i]), max, out=np.zeros_like(dataref[:, i]), where=max != 0) # Replace any division by zero with zero datadiffs[:, i] = np.divide(np.abs(datatest[:, i] - dataref[:, i]), maxi, out=np.zeros_like(dataref[:, i]), where=maxi != 0) # Replace any division by zero with zero
# Calculate power (ignore warning from taking a log of any zero values) # Calculate power (ignore warning from taking a log of any zero values)
with np.errstate(divide='ignore'): with np.errstate(divide='ignore'):

14
testing/models_pmls/pml_3D_pec_plate/pml_3D_pec_plate.py
查看文件

@@ -29,10 +29,16 @@ rx = gprMax.Rx(p1=(0.038, 0.114, 0.013))
plate = gprMax.Plate(p1=(0.013, 0.013, 0.013), plate = gprMax.Plate(p1=(0.013, 0.013, 0.013),
p2=(0.038, 0.113, 0.013), material_id='pec') p2=(0.038, 0.113, 0.013), material_id='pec')
gv1 = gprMax.GeometryView(p1=(0, 0, 0), p2=(x, y, z), dl=(dl, dl, dl), gv1 = gprMax.GeometryView(p1=(0, 0, 0),
filename=Path(*parts[:-1], parts[-1] + '_n'), output_type='n') p2=(x, y, z),
gv2 = gprMax.GeometryView(p1=(0, 0, 0), p2=(x, y, z), dl=(dl, dl, dl), dl=(dl, dl, dl),
filename=Path(*parts[:-1], parts[-1] + '_f'), output_type='f') filename=Path(*parts[:-1], f'{parts[-1]}_n'),
output_type='n',)
gv2 = gprMax.GeometryView(p1=(0, 0, 0),
p2=(x, y, z),
dl=(dl, dl, dl),
filename=Path(*parts[:-1], f'{parts[-1]}_f'),
output_type='f',)
pmls = {'CFS-PML': {'pml': gprMax.PMLProps(formulation='HORIPML', thickness=10), pmls = {'CFS-PML': {'pml': gprMax.PMLProps(formulation='HORIPML', thickness=10),
# Parameters from http://dx.doi.org/10.1109/TAP.2018.2823864 # Parameters from http://dx.doi.org/10.1109/TAP.2018.2823864

14
testing/models_pmls/pml_3D_pec_plate/pml_3D_pec_plate_ref.py
查看文件

@@ -29,10 +29,16 @@ rx = gprMax.Rx(p1=(0.113, 0.189, 0.088))
plate = gprMax.Plate(p1=(0.088, 0.088, 0.088), plate = gprMax.Plate(p1=(0.088, 0.088, 0.088),
p2=(0.113, 0.188, 0.088), material_id='pec') p2=(0.113, 0.188, 0.088), material_id='pec')
gv1 = gprMax.GeometryView(p1=(0, 0, 0), p2=(x, y, z), dl=(dl, dl, dl), gv1 = gprMax.GeometryView(p1=(0, 0, 0),
filename=Path(*parts[:-1], parts[-1] + '_n'), output_type='n') p2=(x, y, z),
gv2 = gprMax.GeometryView(p1=(0, 0, 0), p2=(x, y, z), dl=(dl, dl, dl), dl=(dl, dl, dl),
filename=Path(*parts[:-1], parts[-1] + '_f'), output_type='f') filename=Path(*parts[:-1], f'{parts[-1]}_n'),
output_type='n',)
gv2 = gprMax.GeometryView(p1=(0, 0, 0),
p2=(x, y, z),
dl=(dl, dl, dl),
filename=Path(*parts[:-1], f'{parts[-1]}_f'),
output_type='f',)
pml = gprMax.PMLProps(formulation='HORIPML', thickness=10) pml = gprMax.PMLProps(formulation='HORIPML', thickness=10)

8
testing/test_experimental.py
查看文件

@@ -81,8 +81,10 @@ realmax = np.where(np.abs(real[:, 1]) == 1)[0][0]
difftime = - (timemodel[modelmax] - real[realmax, 0]) difftime = - (timemodel[modelmax] - real[realmax, 0])
# Plot modelled and real data # Plot modelled and real data
fig, ax = plt.subplots(num=modelfile.stem + '_vs_' + realfile.stem, fig, ax = plt.subplots(num=f'{modelfile.stem}_vs_{realfile.stem}',
figsize=(20, 10), facecolor='w', edgecolor='w') figsize=(20, 10),
facecolor='w',
edgecolor='w',)
ax.plot(timemodel + difftime, model, 'r', lw=2, label='Model') ax.plot(timemodel + difftime, model, 'r', lw=2, label='Model')
ax.plot(real[:, 0], real[:, 1], 'r', ls='--', lw=2, label='Experiment') ax.plot(real[:, 0], real[:, 1], 'r', ls='--', lw=2, label='Experiment')
ax.set_xlabel('Time [s]') ax.set_xlabel('Time [s]')
@@ -93,7 +95,7 @@ ax.legend()
ax.grid() ax.grid()
# Save a PDF/PNG of the figure # Save a PDF/PNG of the figure
savename = modelfile.stem + '_vs_' + realfile.stem savename = f'{modelfile.stem}_vs_{realfile.stem}'
savename = modelfile.parent / savename savename = modelfile.parent / savename
# fig.savefig(savename.with_suffix('.pdf'), dpi=None, format='pdf', # fig.savefig(savename.with_suffix('.pdf'), dpi=None, format='pdf',
# bbox_inches='tight', pad_inches=0.1) # bbox_inches='tight', pad_inches=0.1)

24
testing/test_models.py
查看文件

@@ -112,11 +112,9 @@ for i, model in enumerate(testmodels):
filetest.attrs['dt'], filetest.attrs['dt'],
filetest.attrs['dx_dy_dz'], rxposrelative) filetest.attrs['dx_dy_dz'], rxposrelative)
filetest.close()
else: else:
# Get output for model and reference files # Get output for model and reference files
fileref = file.stem + '_ref' fileref = f'{file.stem}_ref'
fileref = file.parent / Path(fileref) fileref = file.parent / Path(fileref)
fileref = h5py.File(fileref.with_suffix('.h5'), 'r') fileref = h5py.File(fileref.with_suffix('.h5'), 'r')
filetest = h5py.File(file.with_suffix('.h5'), 'r') filetest = h5py.File(file.with_suffix('.h5'), 'r')
@@ -164,10 +162,10 @@ for i, model in enumerate(testmodels):
# Diffs # Diffs
datadiffs = np.zeros(datatest.shape, dtype=np.float64) datadiffs = np.zeros(datatest.shape, dtype=np.float64)
for i in range(len(outputstest)): for i in range(len(outputstest)):
max = np.amax(np.abs(dataref[:, i])) maxi = np.amax(np.abs(dataref[:, i]))
datadiffs[:, i] = np.divide(np.abs(dataref[:, i] - datatest[:, i]), max, datadiffs[:, i] = np.divide(np.abs(dataref[:, i] - datatest[:, i]), maxi,
out=np.zeros_like(dataref[:, i]), out=np.zeros_like(dataref[:, i]),
where=max != 0) # Replace any division by zero with zero where=maxi != 0) # Replace any division by zero with zero
# Calculate power (ignore warning from taking a log of any zero values) # Calculate power (ignore warning from taking a log of any zero values)
with np.errstate(divide='ignore'): with np.errstate(divide='ignore'):
@@ -188,17 +186,17 @@ for i, model in enumerate(testmodels):
facecolor='w', facecolor='w',
edgecolor='w') edgecolor='w')
ex1.plot(timetest, datatest[:, 0], 'r', lw=2, label=model) ex1.plot(timetest, datatest[:, 0], 'r', lw=2, label=model)
ex1.plot(timeref, dataref[:, 0], 'g', lw=2, ls='--', label=model + '(Ref)') ex1.plot(timeref, dataref[:, 0], 'g', lw=2, ls='--', label=f'{model}(Ref)')
ey1.plot(timetest, datatest[:, 1], 'r', lw=2, label=model) ey1.plot(timetest, datatest[:, 1], 'r', lw=2, label=model)
ey1.plot(timeref, dataref[:, 1], 'g', lw=2, ls='--', label=model + '(Ref)') ey1.plot(timeref, dataref[:, 1], 'g', lw=2, ls='--', label=f'{model}(Ref)')
ez1.plot(timetest, datatest[:, 2], 'r', lw=2, label=model) ez1.plot(timetest, datatest[:, 2], 'r', lw=2, label=model)
ez1.plot(timeref, dataref[:, 2], 'g', lw=2, ls='--', label=model + '(Ref)') ez1.plot(timeref, dataref[:, 2], 'g', lw=2, ls='--', label=f'{model}(Ref)')
hx1.plot(timetest, datatest[:, 3], 'r', lw=2, label=model) hx1.plot(timetest, datatest[:, 3], 'r', lw=2, label=model)
hx1.plot(timeref, dataref[:, 3], 'g', lw=2, ls='--', label=model + '(Ref)') hx1.plot(timeref, dataref[:, 3], 'g', lw=2, ls='--', label=f'{model}(Ref)')
hy1.plot(timetest, datatest[:, 4], 'r', lw=2, label=model) hy1.plot(timetest, datatest[:, 4], 'r', lw=2, label=model)
hy1.plot(timeref, dataref[:, 4], 'g', lw=2, ls='--', label=model + '(Ref)') hy1.plot(timeref, dataref[:, 4], 'g', lw=2, ls='--', label=f'{model}(Ref)')
hz1.plot(timetest, datatest[:, 5], 'r', lw=2, label=model) hz1.plot(timetest, datatest[:, 5], 'r', lw=2, label=model)
hz1.plot(timeref, dataref[:, 5], 'g', lw=2, ls='--', label=model + '(Ref)') hz1.plot(timeref, dataref[:, 5], 'g', lw=2, ls='--', label=f'{model}(Ref)')
ylabels = ['$E_x$, field strength [V/m]', '$H_x$, field strength [A/m]', ylabels = ['$E_x$, field strength [V/m]', '$H_x$, field strength [A/m]',
'$E_y$, field strength [V/m]', '$H_y$, field strength [A/m]', '$E_y$, field strength [V/m]', '$H_y$, field strength [A/m]',
'$E_z$, field strength [V/m]', '$H_z$, field strength [A/m]'] '$E_z$, field strength [V/m]', '$H_z$, field strength [A/m]']
@@ -232,7 +230,7 @@ for i, model in enumerate(testmodels):
ax.grid() ax.grid()
# Save a PDF/PNG of the figure # Save a PDF/PNG of the figure
filediffs = file.stem + '_diffs' filediffs = f'{file.stem}_diffs'
filediffs = file.parent / Path(filediffs) filediffs = file.parent / Path(filediffs)
# fig1.savefig(file.with_suffix('.pdf'), dpi=None, format='pdf', # fig1.savefig(file.with_suffix('.pdf'), dpi=None, format='pdf',
# bbox_inches='tight', pad_inches=0.1) # bbox_inches='tight', pad_inches=0.1)