diff --git a/gprMax/cmds_multiple.py b/gprMax/cmds_multiple.py index 9fc47bdc..90fff1e3 100644 --- a/gprMax/cmds_multiple.py +++ b/gprMax/cmds_multiple.py @@ -664,7 +664,7 @@ class Snapshot(UserObjectMulti): #else: s = SnapshotUser(xs, ys, zs, xf, yf, zf, dx, dy, dz, iterations, filename) - log.info(f'Snapshot from {xs * grid.dx:g}m, {ys * grid.dy:g}m, {zs * grid.dz:g}m, to {xf * grid.dx:g}m, {yf * grid.dy:g}m, {zf * grid.dz:g}m, discretisation {dx * grid.dx:g}m, {dy * grid.dy:g}m, {dz * grid.dz:g}m, at {s.time * grid.dt:g} secs with filename {s.basefilename} created.') + log.info(f'Snapshot from {xs * grid.dx:g}m, {ys * grid.dy:g}m, {zs * grid.dz:g}m, to {xf * grid.dx:g}m, {yf * grid.dy:g}m, {zf * grid.dz:g}m, discretisation {dx * grid.dx:g}m, {dy * grid.dy:g}m, {dz * grid.dz:g}m, at {s.time * grid.dt:g} secs with filename {s.filename} created.') grid.snapshots.append(s) diff --git a/gprMax/cmds_single_use.py b/gprMax/cmds_single_use.py index 85d169e9..b273c1d7 100644 --- a/gprMax/cmds_single_use.py +++ b/gprMax/cmds_single_use.py @@ -504,7 +504,7 @@ class OutputDir(UserObjectSingle): self.order = 11 def create(self, grid, uip): - grid.outputdirectory = self.kwargs['dir'] + config.model_configs[grid.model_num].set_output_file_path(self.kwargs['dir']) class NumberOfModelRuns(UserObjectSingle): diff --git a/gprMax/config.py b/gprMax/config.py index b89eb268..eab33301 100644 --- a/gprMax/config.py +++ b/gprMax/config.py @@ -73,22 +73,20 @@ class ModelConfig: self.reuse_geometry = False + # String to print at start of each model run + inputfilestr = f'\n--- Model {self.i + 1}/{sim_config.model_end}, input file: {sim_config.input_file_path}' + self.set_inputfilestr(inputfilestr) + if not sim_config.single_model: self.appendmodelnumber = str(self.i + 1) # Indexed from 1 else: self.appendmodelnumber = '' # Output file path for specific model - parts = sim_config.output_file_path.with_suffix('').parts - self.output_file_path = Path(*parts[:-1], parts[-1] + self.appendmodelnumber) - self.output_file_path_ext = self.output_file_path.with_suffix('.out') + self.set_output_file_path() - # Make a snapshot directory - self.snapshot_dir = '_snaps' - - # String to print at start of each model run - inputfilestr = f'\n--- Model {self.i + 1}/{sim_config.model_end}, input file: {sim_config.input_file_path}' - self.set_inputfilestr(inputfilestr) + # Specify a snapshot directory + self.set_snapshots_file_path() # Numerical dispersion analysis parameters # highestfreqthres: threshold (dB) down from maximum power (0dB) of main frequency used @@ -129,6 +127,34 @@ class ModelConfig: """ self.inputfilestr = Fore.GREEN + f"{inputfilestr} {'-' * (get_terminal_width() - 1 - len(inputfilestr))}\n" + Style.RESET_ALL + def set_output_file_path(self, outputdir=None): + """Output file path can be provided by the user via the API or an input file + command. If they haven't provided one use the input file path instead. + + Args: + outputdir (str): Output file directory given from input file command. + """ + + if not outputdir: + try: + self.output_file_path = Path(self.args.outputfile) + except AttributeError: + self.output_file_path = sim_config.input_file_path.with_suffix('') + else: + try: + Path(outputdir).mkdir(exist_ok=True) + self.output_file_path = Path(outputdir, sim_config.input_file_path.stem) + except AttributeError: + self.output_file_path = sim_config.input_file_path.with_suffix('') + + parts = self.output_file_path.parts + self.output_file_path = Path(*parts[:-1], parts[-1] + self.appendmodelnumber) + self.output_file_path_ext = self.output_file_path.with_suffix('.out') + + def set_snapshots_file_path(self): + """Set directory to store any snapshots.""" + parts = self.output_file_path.with_suffix('').parts + self.snapshot_file_path = Path(*parts[:-1], parts[-1] + '_snaps') class SimulationConfig: """Configuration parameters for a standard simulation. @@ -200,8 +226,8 @@ class SimulationConfig: # Set more complex parameters self.set_precision() + self.get_byteorder() self.set_input_file_path() - self.set_output_file_path() self.set_model_start_end() self.set_single_model() @@ -237,14 +263,22 @@ class SimulationConfig: 'cython_float_or_double': cython.float, 'cython_complex': cython.floatcomplex, 'C_float_or_double': 'float', - 'C_complex': 'pycuda::complex'} + 'C_complex': 'pycuda::complex', + 'vtk_float': 'Float32'} elif self.general['precision'] == 'double': self.dtypes = {'float_or_double': np.float64, 'complex': np.complex128, 'cython_float_or_double': cython.double, 'cython_complex': cython.doublecomplex, 'C_float_or_double': 'double', - 'C_complex': 'pycuda::complex'} + 'C_complex': 'pycuda::complex', + 'vtk_float': 'Float64'} + + def get_byteorder(self): + """Check 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_single_model(self): if self.model_start == 0 and self.model_end == 1: @@ -269,20 +303,12 @@ class SimulationConfig: self.model_end = modelend def set_input_file_path(self): - """If the API is in use an id for the simulation must be provided.""" + """Set input file path for CLI or API.""" if self.args.inputfile is None: self.input_file_path = Path(self.args.outputfile) else: self.input_file_path = Path(self.args.inputfile) - def set_output_file_path(self): - """Output file path can be provided by the user. If they havent provided one - use the inputfile file path instead.""" - try: - self.output_file_path = Path(self.args.outputfile) - except AttributeError: - self.output_file_path = Path(self.input_file_path) - class SimulationConfigMPI(SimulationConfig): """Configuration parameters for a MPI simulation. diff --git a/gprMax/fields_outputs.py b/gprMax/fields_outputs.py index 87d5cb4d..5aa5e823 100644 --- a/gprMax/fields_outputs.py +++ b/gprMax/fields_outputs.py @@ -31,7 +31,7 @@ def store_outputs(G): iteration (int): Current iteration number. Ex, Ey, Ez, Hx, Hy, Hz (memory view): Current electric and magnetic field values. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ iteration = G.iteration @@ -105,7 +105,7 @@ def write_hdf5_main_grid_outputfile(outputfile, G): Args: outputfile (str): Name of the output file. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ write_data(outputfile, G) @@ -116,7 +116,7 @@ def write_hdf5_sub_grid_outputfile(outputfile, G): Args: outputfile (str): Name of the output file. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ stem = outputfile.stem @@ -145,7 +145,7 @@ def write_data(outputfile, G): Args: outputfile (str): Name of the output file. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ f = h5py.File(outputfile, 'w') diff --git a/gprMax/geometry_outputs.py b/gprMax/geometry_outputs.py index eb35b7e4..dc5c7561 100644 --- a/gprMax/geometry_outputs.py +++ b/gprMax/geometry_outputs.py @@ -16,6 +16,7 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . +import logging import os from pathlib import Path import sys @@ -31,16 +32,13 @@ from .cython.geometry_outputs import define_normal_geometry from .cython.geometry_outputs import define_fine_geometry from .utilities import round_value +log = logging.getLogger(__name__) + class GeometryView: """Views of the geometry of the model.""" - if sys.byteorder == 'little': - byteorder = 'LittleEndian' - else: - byteorder = 'BigEndian' - - def __init__(self, xs=None, ys=None, zs=None, xf=None, yf=None, zf=None, dx=None, dy=None, dz=None, filename=None, fileext=None, grid=None): + def __init__(self, xs=None, ys=None, zs=None, xf=None, yf=None, zf=None, dx=None, dy=None, dz=None, filename=None, fileext=None, G=None): """ Args: xs, xf, ys, yf, zs, zf (int): Extent of the volume in cells. @@ -48,6 +46,7 @@ class GeometryView: filename (str): Filename to save to. fileext (str): File extension of VTK file - either '.vti' for a per cell geometry view, or '.vtp' for a per cell edge geometry view. + G (FDTDGrid): Parameters describing a grid in a model. """ self.xs = xs @@ -64,6 +63,7 @@ class GeometryView: self.dz = dz self.basefilename = filename self.fileext = fileext + self.G = G if self.fileext == '.vti': # Calculate number of cells according to requested sampling for geometry view @@ -77,9 +77,9 @@ class GeometryView: self.vtk_nycells = round_value(self.ny / self.dy) self.vtk_nzcells = round_value(self.nz / self.dz) self.vtk_ncells = self.vtk_nxcells * self.vtk_nycells * self.vtk_nzcells - self.datawritesize = (np.dtype(np.uint32).itemsize * self.vtk_ncells - + 2 * np.dtype(np.int8).itemsize * self.vtk_ncells - + 3 * np.dtype(np.uint32).itemsize) + self.datawritesize = (np.dtype(np.uint32).itemsize * self.vtk_ncells + + 2 * np.dtype(np.int8).itemsize * self.vtk_ncells + + 3 * np.dtype(np.uint32).itemsize) elif self.fileext == '.vtp': self.vtk_numpoints = (self.nx + 1) * (self.ny + 1) * (self.nz + 1) @@ -89,45 +89,44 @@ class GeometryView: self.vtk_nylines = self.ny * (self.nx + 1) * (self.nz + 1) self.vtk_nzlines = self.nz * (self.nx + 1) * (self.ny + 1) self.vtk_numlines = self.vtk_nxlines + self.vtk_nylines + self.vtk_nzlines - self.vtk_connectivity_offset = round_value((self.vtk_numpoints - * self.vtk_numpoint_components - * np.dtype(np.float32).itemsize) - + np.dtype(np.uint32).itemsize) - self.vtk_offsets_offset = round_value(self.vtk_connectivity_offset - + (self.vtk_numlines * self.vtk_numline_components * np.dtype(np.uint32).itemsize) - + np.dtype(np.uint32).itemsize) - self.vtk_materials_offset = round_value(self.vtk_offsets_offset - + (self.vtk_numlines * np.dtype(np.uint32).itemsize) - + np.dtype(np.uint32).itemsize) - vtk_cell_offsets = ((self.vtk_numline_components * self.vtk_numlines) - + self.vtk_numline_components - self.vtk_numline_components - 1) // self.vtk_numline_components + 1 - self.datawritesize = (np.dtype(np.float32).itemsize * self.vtk_numpoints * self.vtk_numpoint_components - + np.dtype(np.uint32).itemsize * self.vtk_numlines * self.vtk_numline_components - + np.dtype(np.uint32).itemsize * self.vtk_numlines - + np.dtype(np.uint32).itemsize * vtk_cell_offsets - + np.dtype(np.uint32).itemsize * 4) + self.vtk_connectivity_offset = round_value((self.vtk_numpoints * + self.vtk_numpoint_components * + np.dtype(np.float32).itemsize) + + np.dtype(np.uint32).itemsize) + self.vtk_offsets_offset = round_value(self.vtk_connectivity_offset + + (self.vtk_numlines * + self.vtk_numline_components * + np.dtype(np.uint32).itemsize) + + np.dtype(np.uint32).itemsize) + self.vtk_materials_offset = round_value(self.vtk_offsets_offset + + (self.vtk_numlines * + np.dtype(np.uint32).itemsize) + + np.dtype(np.uint32).itemsize) + vtk_cell_offsets = (((self.vtk_numline_components * self.vtk_numlines) + + self.vtk_numline_components - self.vtk_numline_components - 1) // + self.vtk_numline_components + 1) + self.datawritesize = (np.dtype(np.float32).itemsize * self.vtk_numpoints * + self.vtk_numpoint_components + np.dtype(np.uint32).itemsize * + self.vtk_numlines * self.vtk_numline_components + + np.dtype(np.uint32).itemsize * self.vtk_numlines + + np.dtype(np.uint32).itemsize * vtk_cell_offsets + + np.dtype(np.uint32).itemsize * 4) - def set_filename(self, appendmodelnumber): - """ - Construct filename from user-supplied name and model run number. - - Args: - appendmodelnumber (str): Text to append to filename. - """ - - parts = config.sim_config.input_file_path.parts - self.filename = Path(*parts[:-1], parts[-1] + appendmodelnumber) + def set_filename(self): + """Construct filename from user-supplied name and model run number.""" + parts = config.model_configs[self.G.model_num].output_file_path.with_suffix('').parts + self.filename = Path(*parts[:-1], parts[-1]) self.filename = self.filename.with_suffix(self.fileext) def write_vtk(self, G, pbar): - """ - Writes the geometry information to a VTK file. Either ImageData (.vti) for a - per-cell geometry view, or PolygonalData (.vtp) for a per-cell-edge geometry view. + """Writes the geometry information to a VTK file. + Either ImageData (.vti) for a per-cell geometry view, or + PolygonalData (.vtp) for a per-cell-edge geometry view. N.B. No Python 3 support for VTK at time of writing (03/2015) Args: - G (class): Grid class instance - holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. pbar (class): Progress bar class instance. """ @@ -143,18 +142,29 @@ class GeometryView: for index, rx in enumerate(G.rxs): self.rxs[rx.xcoord, rx.ycoord, rx.zcoord] = index + 1 - vtk_srcs_pml_offset = round_value((np.dtype(np.uint32).itemsize * self.vtk_nxcells * self.vtk_nycells * self.vtk_nzcells) + np.dtype(np.uint32).itemsize) - vtk_rxs_offset = round_value((np.dtype(np.uint32).itemsize * self.vtk_nxcells * self.vtk_nycells * self.vtk_nzcells) + np.dtype(np.uint32).itemsize + (np.dtype(np.int8).itemsize * self.vtk_nxcells * self.vtk_nycells * self.vtk_nzcells) + np.dtype(np.uint32).itemsize) + vtk_srcs_pml_offset = round_value((np.dtype(np.uint32).itemsize * + self.vtk_nxcells * + self.vtk_nycells * + self.vtk_nzcells) + + np.dtype(np.uint32).itemsize) + vtk_rxs_offset = round_value((np.dtype(np.uint32).itemsize * + self.vtk_nxcells * self.vtk_nycells * + self.vtk_nzcells) + + np.dtype(np.uint32).itemsize + + (np.dtype(np.int8).itemsize * + self.vtk_nxcells * self.vtk_nycells * + self.vtk_nzcells) + + np.dtype(np.uint32).itemsize) with open(self.filename, 'wb') as f: f.write('\n'.encode('utf-8')) - f.write('\n'.format(GeometryView.byteorder).encode('utf-8')) - f.write('\n'.format(self.vtk_xscells, self.vtk_xfcells, self.vtk_yscells, self.vtk_yfcells, self.vtk_zscells, self.vtk_zfcells, self.dx * G.dx, self.dy * G.dy, self.dz * G.dz).encode('utf-8')) - f.write('\n'.format(self.vtk_xscells, self.vtk_xfcells, self.vtk_yscells, self.vtk_yfcells, self.vtk_zscells, self.vtk_zfcells).encode('utf-8')) + f.write(f'\n'.encode('utf-8')) + f.write(f'\n'.encode('utf-8')) + f.write(f'\n'.encode('utf-8')) f.write('\n'.encode('utf-8')) f.write('\n'.encode('utf-8')) - f.write('\n'.format(vtk_srcs_pml_offset).encode('utf-8')) - f.write('\n'.format(vtk_rxs_offset).encode('utf-8')) + f.write(f'\n'.encode('utf-8')) + f.write(f'\n'.encode('utf-8')) f.write('\n'.encode('utf-8')) f.write('\n\n\n_'.encode('utf-8')) @@ -207,33 +217,30 @@ class GeometryView: elif self.fileext == '.vtp': with open(self.filename, 'wb') as f: f.write('\n'.encode('utf-8')) - f.write('\n'.format(GeometryView.byteorder).encode('utf-8')) - f.write('\n\n'.format(self.vtk_numpoints, self.vtk_numlines).encode('utf-8')) - + f.write(f'\n'.encode('utf-8')) + f.write(f'\n\n'.encode('utf-8')) f.write('\n\n\n'.encode('utf-8')) - f.write('\n\n'.format(self.vtk_connectivity_offset).encode('utf-8')) - f.write('\n\n'.format(self.vtk_offsets_offset).encode('utf-8')) - + f.write(f'\n\n'.encode('utf-8')) + f.write(f'\n\n'.encode('utf-8')) f.write('\n'.encode('utf-8')) f.write('\n'.format(self.vtk_materials_offset).encode('utf-8')) f.write('\n'.encode('utf-8')) - f.write('\n\n\n_'.encode('utf-8')) # Coordinates of each point points = np.zeros((self.vtk_numpoints, 3), dtype=np.float32) - # Number of x components - # Node connectivity. Each index contains a pair of connected x nodes x_lines = np.zeros((self.vtk_nxlines, 2), dtype=np.uint32) # Material at Ex location in Yee cell. x_materials = np.zeros((self.vtk_nxlines), dtype=np.uint32) - + # Node connectivity. Each index contains a pair of connected y nodes y_lines = np.zeros((self.vtk_nylines, 2), dtype=np.uint32) + # Material at Ey location in Yee cell. y_materials = np.zeros((self.vtk_nylines), dtype=np.uint32) - + # Node connectivity. Each index contains a pair of connected z nodes z_lines = np.zeros((self.vtk_nzlines, 2), dtype=np.uint32) + # Material at Ez location in Yee cell. z_materials = np.zeros((self.vtk_nzlines), dtype=np.uint32) define_fine_geometry(self.nx, @@ -263,7 +270,8 @@ class GeometryView: pbar.update(n=points.nbytes) # Write connectivity data - f.write(pack('I', np.dtype(np.uint32).itemsize * self.vtk_numlines * self.vtk_numline_components)) + f.write(pack('I', np.dtype(np.uint32).itemsize * + self.vtk_numlines * self.vtk_numline_components)) pbar.update(n=4) f.write(x_lines) pbar.update(n=x_lines.nbytes) @@ -275,7 +283,9 @@ class GeometryView: # Write cell type (line) offsets f.write(pack('I', np.dtype(np.uint32).itemsize * self.vtk_numlines)) pbar.update(n=4) - for vtk_offsets in range(self.vtk_numline_components, (self.vtk_numline_components * self.vtk_numlines) + self.vtk_numline_components, self.vtk_numline_components): + for vtk_offsets in range(self.vtk_numline_components, + (self.vtk_numline_components * self.vtk_numlines) + + self.vtk_numline_components, self.vtk_numline_components): f.write(pack('I', vtk_offsets)) pbar.update(n=4) @@ -293,25 +303,24 @@ class GeometryView: self.write_gprmax_info(f, G, materialsonly=True) def write_gprmax_info(self, f, G, materialsonly=False): - """ - Writes gprMax specific information relating material, source, - and receiver names to numeric identifiers. + """Writes gprMax specific information relating material, source, + and receiver names to numeric identifiers. Args: f (filehandle): VTK file. - G (class): Grid class instance - holds essential parameters describing the model. - materialsonly (boolean): Only write information on materials + G (FDTDGrid): Parameters describing a grid in a model. + materialsonly (bool): Only write information on materials """ f.write('\n\n\n'.encode('utf-8')) for material in G.materials: - f.write('{}\n'.format(material.ID, material.numID).encode('utf-8')) + f.write(f'{material.numID}\n'.encode('utf-8')) if not materialsonly: f.write('1\n'.encode('utf-8')) for index, src in enumerate(G.hertziandipoles + G.magneticdipoles + G.voltagesources + G.transmissionlines): - f.write('{}\n'.format(src.ID, index + 2).encode('utf-8')) + f.write(f'{index + 2}\n'.encode('utf-8')) for index, rx in enumerate(G.rxs): - f.write('{}\n'.format(rx.ID, index + 1).encode('utf-8')) + f.write(f'{index + 1}\n'.encode('utf-8')) f.write('\n'.encode('utf-8')) @@ -334,8 +343,13 @@ class GeometryObjects: self.nx = self.xf - self.xs self.ny = self.yf - self.ys self.nz = self.zf - self.zs - self.filename = basefilename + '.h5' - self.materialsfilename = basefilename + '_materials.txt' + + # Set filenames + parts = config.sim_config.input_file_path.with_suffix('').parts + self.filename_hdf5 = Path(*parts[:-1], basefilename) + self.filename_hdf5 = self.filename_hdf5.with_suffix('.h5') + self.filename_materials = Path(*parts[:-1], basefilename + '_materials') + self.filename_materials = self.filename_materials.with_suffix('.txt') # Sizes of arrays to write necessary to update progress bar self.solidsize = (self.nx + 1) * (self.ny + 1) * (self.nz + 1) * np.dtype(np.int16).itemsize @@ -347,12 +361,12 @@ class GeometryObjects: """Write a geometry objects file in HDF5 format. Args: - G (class): Grid class instance - holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. pbar (class): Progress bar class instance. """ # Write the geometry objects to a HDF5 file - fdata = h5py.File(os.path.abspath(os.path.join(G.inputdirectory, self.filename)), 'w') + fdata = h5py.File(self.filename_hdf5, 'w') fdata.attrs['gprMax'] = __version__ fdata.attrs['Title'] = G.title fdata.attrs['dx_dy_dz'] = (G.dx, G.dy, G.dz) @@ -370,27 +384,28 @@ class GeometryObjects: # Write materials list to a text file # This includes all materials in range whether used in volume or not - fmaterials = open(os.path.abspath(os.path.join(G.inputdirectory, self.materialsfilename)), 'w') + fmaterials = open(self.filename_materials, 'w') for numID in range(minmat, maxmat + 1): for material in G.materials: if material.numID == numID: - fmaterials.write('#material: {:g} {:g} {:g} {:g} {}\n'.format(material.er, material.se, material.mr, material.sm, material.ID)) + fmaterials.write(f'#material: {material.er:g} {material.se:g} {material.mr:g} {material.sm:g} {material.ID}\n') if material.poles > 0: if 'debye' in material.type: - dispersionstr = '#add_dispersion_debye: {:g} '.format(material.poles) + dispersionstr = f'#add_dispersion_debye: {material.poles:g} ' for pole in range(material.poles): - dispersionstr += '{:g} {:g} '.format(material.deltaer[pole], material.tau[pole]) + dispersionstr += f'{material.deltaer[pole]:g} {material.tau[pole]:g} ' elif 'lorenz' in material.type: - dispersionstr = '#add_dispersion_lorenz: {:g} '.format(material.poles) + dispersionstr = f'#add_dispersion_lorenz: {material.poles:g} ' for pole in range(material.poles): - dispersionstr += '{:g} {:g} {:g} '.format(material.deltaer[pole], material.tau[pole], material.alpha[pole]) + dispersionstr += f'{material.deltaer[pole]:g} {material.tau[pole]:g} {material.alpha[pole]:g} ' elif 'drude' in material.type: - dispersionstr = '#add_dispersion_drude: {:g} '.format(material.poles) + dispersionstr = f'#add_dispersion_drude: {material.poles:g} ' for pole in range(material.poles): - dispersionstr += '{:g} {:g} '.format(material.tau[pole], material.alpha[pole]) + dispersionstr += f'{material.tau[pole]:g} {material.alpha[pole]:g} ' dispersionstr += material.ID fmaterials.write(dispersionstr + '\n') + class GeometryViewFineMultiGrid: """Geometry view for all grids in the simulation. @@ -400,11 +415,6 @@ class GeometryViewFineMultiGrid: view them at once in Paraview. """ - if sys.byteorder == 'little': - byteorder = 'LittleEndian' - else: - byteorder = 'BigEndian' - def __init__(self, xs, ys, zs, xf, yf, zf, dx, dy, dz, filename, fileext, G): """ Args: @@ -432,7 +442,6 @@ class GeometryViewFineMultiGrid: # total additional lines required for subgrid self.additional_lines += sg_gv.n_total_lines # total additional points required for subgrid - self.additional_points += sg_gv.n_total_points self.vtk_numpoints = self.additional_points + (self.nx + 1) * (self.ny + 1) * (self.nz + 1) @@ -444,24 +453,21 @@ class GeometryViewFineMultiGrid: self.vtk_materials_offset = round_value(int(self.vtk_offsets_offset + (self.vtk_numlines * np.dtype(np.uint32).itemsize) + np.dtype(np.uint32).itemsize)) self.datawritesize = np.dtype(np.float32).itemsize * self.vtk_numpoints * self.vtk_numpoint_components + np.dtype(np.uint32).itemsize * self.vtk_numlines * self.vtk_numline_components + np.dtype(np.uint32).itemsize * self.vtk_numlines + np.dtype(np.uint32).itemsize * self.vtk_numlines - def set_filename(self, appendmodelnumber): - """ - Construct filename from user-supplied name and model run number. - - Args: - appendmodelnumber (str): Text to append to filename. - """ - - self.filename = os.path.abspath(os.path.join(os.path.dirname(os.path.abspath(config.general['inputfilepath'])), self.basefilename + appendmodelnumber)) - self.filename += self.fileext + def set_filename(self): + """Construct filename from user-supplied name and model run number.""" + parts = config.model_configs[self.G.model_num].input_file_path.parts + self.filename = Path(*parts[:-1], parts[-1] + config.model_configs[self.G.model_num].appendmodelnumber) + self.filename = self.filename.with_suffix(self.fileext) def write_vtk(self, *args): - """ - Writes the geometry information to a VTK file. Either ImageData (.vti) for a - per-cell geometry view, or PolygonalData (.vtp) for a per-cell-edge geometry view. + """Writes the geometry information to a VTK file. + Either ImageData (.vti) for a per-cell geometry view, or + PolygonalData (.vtp) for a per-cell-edge geometry view. + N.B. No Python 3 support for VTK at time of writing (03/2015) + Args: - G (class): Grid class instance - holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. """ G = self.G @@ -471,7 +477,7 @@ class GeometryViewFineMultiGrid: # refine parameters for subgrid f.write('\n'.encode('utf-8')) - f.write('\n'.format(GeometryViewFineMultiGrid.byteorder).encode('utf-8')) + f.write('\n'.format(config.sim_config.vtk_byteorder).encode('utf-8')) f.write('\n\n'.format(self.vtk_numpoints, self.vtk_numlines).encode('utf-8')) f.write('\n\n\n'.encode('utf-8')) @@ -484,7 +490,7 @@ class GeometryViewFineMultiGrid: f.write('\n\n\n_'.encode('utf-8')) # Write points - print('writing points main grid') + log.info('writing points main grid') datasize = np.dtype(np.float32).itemsize * self.vtk_numpoints * self.vtk_numpoint_components f.write(pack('I', datasize)) for i in range(0, G.nx + 1): @@ -493,7 +499,7 @@ class GeometryViewFineMultiGrid: f.write(pack('fff', i * G.dx, j * G.dy, k * G.dz)) for sg_v in self.sg_views: - print('writing points subgrid') + log.info('writing points subgrid') sg_v.write_points(f, G) n_x_lines = self.nx * (self.ny + 1) * (self.nz + 1) @@ -508,7 +514,7 @@ class GeometryViewFineMultiGrid: z_lines = np.zeros((n_z_lines, 2), dtype=np.uint32) z_materials = np.zeros((n_z_lines), dtype=np.uint32) - print('calculate connectivity main grid') + log.info('calculate connectivity main grid') label = 0 counter_x = 0 counter_y = 0 @@ -540,7 +546,7 @@ class GeometryViewFineMultiGrid: label = label + 1 - print('calculate connectivity subgrids') + log.info('calculate connectivity subgrids') for sg_v in self.sg_views: sg_v.populate_connectivity_and_materials(label) # use the last subgrids label for the next view @@ -584,13 +590,13 @@ class GeometryViewFineMultiGrid: #self.write_gprmax_info(f, G, materialsonly=True) def write_gprmax_info(self, f, G, materialsonly=False): - """ - Writes gprMax specific information relating material, source, - and receiver names to numeric identifiers. + """Writes gprMax specific information relating material, source, + and receiver names to numeric identifiers. + Args: f (filehandle): VTK file. - G (class): Grid class instance - holds essential parameters describing the model. - materialsonly (boolean): Only write information on materials + G (FDTDGrid): Parameters describing a grid in a model. + materialsonly (bool): Only write information on materials. """ f.write('\n\n\n'.encode('utf-8')) @@ -607,6 +613,7 @@ class GeometryViewFineMultiGrid: class SubgridGeometryView: + def __init__(self, sg): self.sg = sg @@ -646,8 +653,10 @@ class SubgridGeometryView: f.write(pack('fff', p_x, p_y, p_z)) def populate_connectivity_and_materials(self, label): - """Label is the starting label. 0 if no other grids are present but +1 the last label used - for a multigrid view""" + """Label is the starting label. 0 if no other grids are present but + +1 the last label used for a multigrid view. + """ + sg = self.sg self.label = label diff --git a/gprMax/grid.py b/gprMax/grid.py index e880f9ab..53fe8e19 100644 --- a/gprMax/grid.py +++ b/gprMax/grid.py @@ -174,21 +174,19 @@ class FDTDGrid: self.updatecoeffsH = np.zeros((len(self.materials), 5), dtype=config.sim_config.dtypes['float_or_double']) - def initialise_dispersive_arrays(self, dtype): - """Initialise arrays for storing coefficients when there are dispersive materials present. - - Args: - dtype (dtype): Dtype to use for dispersive arrays. + def initialise_dispersive_arrays(self): + """Initialise arrays for storing coefficients when there are dispersive + materials present. """ self.Tx = np.zeros((config.model_configs[self.model_num].materials['maxpoles'], - self.nx + 1, self.ny + 1, self.nz + 1), dtype=dtype) + self.nx + 1, self.ny + 1, self.nz + 1), dtype=config.model_configs[self.model_num].materials['dispersivedtype']) self.Ty = np.zeros((config.model_configs[self.model_num].materials['maxpoles'], - self.nx + 1, self.ny + 1, self.nz + 1), dtype=dtype) + self.nx + 1, self.ny + 1, self.nz + 1), dtype=config.model_configs[self.model_num].materials['dispersivedtype']) self.Tz = np.zeros((config.model_configs[self.model_num].materials['maxpoles'], - self.nx + 1, self.ny + 1, self.nz + 1), dtype=dtype) + self.nx + 1, self.ny + 1, self.nz + 1), dtype=config.model_configs[self.model_num].materials['dispersivedtype']) self.updatecoeffsdispersive = np.zeros((len(self.materials), 3 * config.model_configs[self.model_num].materials['maxpoles']), - dtype=dtype) + dtype=config.model_configs[self.model_num].materials['dispersivedtype']) def mem_est_basic(self): """Estimate the amount of memory (RAM) required for grid arrays. @@ -275,6 +273,8 @@ class FDTDGrid: """Clear arrays for field components and PMLs.""" # Clear arrays for field components self.initialise_field_arrays() + if config.model_configs[self.model_num].materials['maxpoles'] != 0: + self.initialise_dispersive_arrays() # Clear arrays for fields in PML for pml in self.pmls: @@ -313,12 +313,19 @@ class CUDAGrid(FDTDGrid): self.bpg = (int(np.ceil(((self.nx + 1) * (self.ny + 1) * (self.nz + 1)) / self.tpb[0])), 1, 1) - def initialise_arrays(self): + def initialise_geometry_arrays(self): + """Initialise an array for cell edge IDs (ID) on GPU.""" + import pycuda.gpuarray as gpuarray + + super().initialise_geometry_arrays() + self.ID_gpu = gpuarray.to_gpu(self.ID) + + def initialise_field_arrays(self): """Initialise geometry and field arrays on GPU.""" import pycuda.gpuarray as gpuarray - self.ID_gpu = gpuarray.to_gpu(self.ID) + super().initialise_field_arrays() self.Ex_gpu = gpuarray.to_gpu(self.Ex) self.Ey_gpu = gpuarray.to_gpu(self.Ey) self.Ez_gpu = gpuarray.to_gpu(self.Ez) @@ -326,22 +333,31 @@ class CUDAGrid(FDTDGrid): self.Hy_gpu = gpuarray.to_gpu(self.Hy) self.Hz_gpu = gpuarray.to_gpu(self.Hz) - def initialise_dispersive_arrays(self, dtype): - """Initialise dispersive material coefficient arrays on GPU. - - Args: - dtype (dtype): Dtype to use for dispersive arrays. - """ + def initialise_dispersive_arrays(self): + """Initialise dispersive material coefficient arrays on GPU.""" import pycuda.gpuarray as gpuarray - super().initialise_dispersive_arrays(dtype) - + super().initialise_dispersive_arrays() self.Tx_gpu = gpuarray.to_gpu(self.Tx) self.Ty_gpu = gpuarray.to_gpu(self.Ty) self.Tz_gpu = gpuarray.to_gpu(self.Tz) self.updatecoeffsdispersive_gpu = gpuarray.to_gpu(self.updatecoeffsdispersive) + def reset_fields(self): + """Clear arrays for field components and PMLs.""" + + super().reset_fields() + + # Clear arrays for field components + self.initialise_field_arrays() + if config.model_configs[self.model_num].materials['maxpoles'] != 0: + self.initialise_dispersive_arrays() + + # Clear arrays for fields in PML + for pml in self.pmls: + pml.initialise_field_arrays() + def memory_check(self, snapsmemsize=0): """Check if model can be run on specified GPU.""" @@ -487,7 +503,7 @@ def Ix(x, y, z, Hx, Hy, Hz, G): Args: x, y, z (float): Coordinates of position in grid. Hx, Hy, Hz (memory view): numpy array of magnetic field values. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ if y == 0 or z == 0: @@ -504,7 +520,7 @@ def Iy(x, y, z, Hx, Hy, Hz, G): Args: x, y, z (float): Coordinates of position in grid. Hx, Hy, Hz (memory view): numpy array of magnetic field values. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ if x == 0 or z == 0: @@ -521,7 +537,7 @@ def Iz(x, y, z, Hx, Hy, Hz, G): Args: x, y, z (float): Coordinates of position in grid. Hx, Hy, Hz (memory view): numpy array of magnetic field values. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ if x == 0 or y == 0: diff --git a/gprMax/hash_cmds_file.py b/gprMax/hash_cmds_file.py index 8e317002..27487a5b 100644 --- a/gprMax/hash_cmds_file.py +++ b/gprMax/hash_cmds_file.py @@ -19,6 +19,7 @@ from io import StringIO import logging import os +from pathlib import Path import sys import gprMax.config as config @@ -157,23 +158,22 @@ def process_include_files(hashcmds, inputfile): return processedincludecmds -def write_processed_file(processedlines, appendmodelnumber): +def write_processed_file(processedlines, G): """Writes an input file after any Python code and include commands in the original input file have been processed. Args: processedlines (list): Input commands after after processing any Python code and include commands. - appendmodelnumber (str): Text to append to filename. + G (FDTDGrid): Parameters describing a grid in a model. """ - processedfile = (os.path.join(config.general['outputfilepath'], - os.path.splitext(config.general['inputfilepath'])[0] + - appendmodelnumber + '_processed.in')) + parts = config.model_configs[G.model_num].output_file_path.parts + processedfile = (Path(*parts[:-1], parts[-1] + '_processed.in')) with open(processedfile, 'w') as f: for item in processedlines: - f.write('{}'.format(item)) + f.write(f'{item}') log.info(f'Written input commands, after processing any Python code and include commands, to file: {processedfile}\n') @@ -282,13 +282,12 @@ def get_user_objects(processedlines, check=True): return user_objs -def parse_hash_commands(model_config, G, scene): +def parse_hash_commands(scene, G): """Parse user hash commands and add them to the scene. Args: - model_config (ModelConfig): Model level configuration object. - G (FDTDGrid): Holds essential parameters describing a model. scene (Scene): Scene object. + G (FDTDGrid): Parameters describing a grid in a model. Returns: scene (Scene): Scene object. @@ -296,7 +295,7 @@ def parse_hash_commands(model_config, G, scene): with open(config.sim_config.input_file_path) as inputfile: - usernamespace = model_config.get_usernamespace() + usernamespace = config.model_configs[G.model_num].get_usernamespace() # Read input file and process any Python and include file commands processedlines = process_python_include_code(inputfile, usernamespace) @@ -311,7 +310,7 @@ def parse_hash_commands(model_config, G, scene): # Write a file containing the input commands after Python or include # file commands have been processed if config.sim_config.args.write_processed: - write_processed_file(processedlines, model_config.appendmodelnumber, G) + write_processed_file(processedlines, G) user_objs = get_user_objects(processedlines, check=True) for user_obj in user_objs: diff --git a/gprMax/materials.py b/gprMax/materials.py index b12d35e0..4fafc21f 100644 --- a/gprMax/materials.py +++ b/gprMax/materials.py @@ -49,7 +49,7 @@ class Material: """Calculates the magnetic update coefficients of the material. Args: - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ HA = (config.m0 * self.mr / G.dt) + 0.5 * self.sm @@ -64,7 +64,7 @@ class Material: """Calculates the electric update coefficients of the material. Args: - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ EA = (config.sim_config.em_consts['e0'] * self.er / G.dt) + 0.5 * self.se @@ -126,7 +126,7 @@ class DispersiveMaterial(Material): """Calculates the electric update coefficients of the material. Args: - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ # The implementation of the dispersive material modelling comes from the @@ -216,7 +216,7 @@ def process_materials(G): store in arrays, and build text list of materials/properties Args: - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. Returns: materialsdata (list): List of material IDs, names, and properties to @@ -315,7 +315,7 @@ class PeplinskiSoil: Args: nbins (int): Number of bins to use to create the different materials. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ # Debye model properties of water @@ -382,7 +382,7 @@ def create_built_in_materials(G): """Create pre-defined (built-in) materials. Args: - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ G.n_built_in_materials = len(G.materials) diff --git a/gprMax/model_build_run.py b/gprMax/model_build_run.py index b9a37089..2bde60a4 100644 --- a/gprMax/model_build_run.py +++ b/gprMax/model_build_run.py @@ -103,6 +103,8 @@ class ModelBuildRun: # Normal model reading/building process; bypassed if geometry information to be reused self.reuse_geometry() if config.model_configs[G.model_num].reuse_geometry else self.build_geometry() + log.info(f'\nOutput path: {config.model_configs[G.model_num].output_file_path.parent}') + # Adjust position of simple sources and receivers if required if G.srcsteps[0] != 0 or G.srcsteps[1] != 0 or G.srcsteps[2] != 0: for source in itertools.chain(G.hertziandipoles, G.magneticdipoles): @@ -135,7 +137,8 @@ class ModelBuildRun: if not (G.geometryviews or G.geometryobjectswrite) and config.sim_config.args.geometry_only: log.warning(Fore.RED + f'\nNo geometry views or geometry objects found.' + Style.RESET_ALL) for i, geometryview in enumerate(G.geometryviews): - geometryview.set_filename(config.model_configs[G.model_num].appendmodelnumber) + log.info('') + geometryview.set_filename() pbar = tqdm(total=geometryview.datawritesize, unit='byte', unit_scale=True, desc=f'Writing geometry view file {i + 1}/{len(G.geometryviews)}, {geometryview.filename.name}', ncols=get_terminal_width() - 1, file=sys.stdout, @@ -143,6 +146,7 @@ class ModelBuildRun: geometryview.write_vtk(G, pbar) pbar.close() for i, geometryobject in enumerate(G.geometryobjectswrite): + log.info('') pbar = tqdm(total=geometryobject.datawritesize, unit='byte', unit_scale=True, desc=f'Writing geometry object file {i + 1}/{len(G.geometryobjectswrite)}, {geometryobject.filename.name}', ncols=get_terminal_width() - 1, file=sys.stdout, @@ -157,22 +161,11 @@ class ModelBuildRun: scene = self.build_scene() - # Combine available grids and check memory requirements + # Combine available grids and check basic memory requirements grids = [G] + G.subgrids for grid in grids: config.model_configs[G.model_num].mem_use += grid.mem_est_basic() mem_check(config.model_configs[G.model_num].mem_use) - log.info(f'\nMemory (RAM) required: ~{human_size(config.model_configs[G.model_num].mem_use)}') - - gridbuilders = [GridBuilder(grid) for grid in grids] - - for gb in gridbuilders: - pml_information(gb.grid) - gb.build_pmls() - gb.build_components() - gb.tm_grid_update() - gb.update_voltage_source_materials() - gb.grid.initialise_std_update_coeff_arrays() # Set datatype for dispersive arrays if there are any dispersive materials. if config.model_configs[G.model_num].materials['maxpoles'] != 0: @@ -187,10 +180,6 @@ class ModelBuildRun: # Update estimated memory (RAM) usage config.model_configs[G.model_num].mem_use += G.mem_est_dispersive() mem_check(config.model_configs[G.model_num].mem_use) - log.info(f'Memory (RAM) required - updated (dispersive): ~{human_size(config.model_configs[G.model_num].mem_use)}') - - for gb in gridbuilders: - gb.grid.initialise_dispersive_arrays(config.model_configs[G.model_num].materials['dispersivedtype']) # Check there is sufficient memory to store any snapshots if G.snapshots: @@ -203,10 +192,20 @@ class ModelBuildRun: mem_check(config.model_configs[G.model_num].mem_use) if config.sim_config.general['cuda']: mem_check_gpu_snaps(G.model_num, snaps_mem) - log.info(f'Memory (RAM) required - updated (snapshots): ~{human_size(config.model_configs[G.model_num].mem_use)}') - # Build materials + log.info(f'\nMemory (RAM) required: ~{human_size(config.model_configs[G.model_num].mem_use)}') + + # Build grids + gridbuilders = [GridBuilder(grid) for grid in grids] for gb in gridbuilders: + pml_information(gb.grid) + gb.build_pmls() + gb.build_components() + gb.tm_grid_update() + gb.update_voltage_source_materials() + gb.grid.initialise_std_update_coeff_arrays() + if config.model_configs[G.model_num].materials['maxpoles'] != 0: + gb.grid.initialise_dispersive_arrays() gb.build_materials() # Check to see if numerical dispersion might be a problem @@ -237,25 +236,13 @@ class ModelBuildRun: if not scene: scene = Scene() # Parse the input file into user objects and add them to the scene - scene = parse_hash_commands(config.model_configs[self.G.model_num], self.G, scene) + scene = parse_hash_commands(scene, self.G) # Creates the internal simulation objects scene.create_internal_objects(self.G) return scene - def create_output_directory(self): - log.debug('Fix output directory path setting') - # if self.G.outputdirectory: - # # Check and set output directory and filename - # try: - # os.mkdir(self.G.outputdirectory) - # log.info(f'\nCreated output directory: {self.G.outputdirectory}') - # except FileExistsError: - # pass - # # Modify the output path (hack) - # config.model_configs[G.model_num].output_file_path_ext = Path(self.G.outputdirectory, config.model_configs[G.model_num].output_file_path_ext) - def write_output_data(self): """Write output data, i.e. field data for receivers and snapshots to file(s). @@ -266,19 +253,17 @@ class ModelBuildRun: # Write any snapshots to file if self.G.snapshots: - # Create directory and construct filename from user-supplied name - # and model run number - snapshotdir = config.model_configs[self.G.model_num].snapshot_dir - if not os.path.exists(snapshotdir): - os.mkdir(snapshotdir) + # Create directory for snapshots + config.model_configs[self.G.model_num].set_snapshots_file_path() + snapshotdir = config.model_configs[self.G.model_num].snapshot_file_path + snapshotdir.mkdir(exist_ok=True) log.info('') for i, snap in enumerate(self.G.snapshots): - fn = snapshotdir / Path(config.model_configs[self.G.model_num].output_file_path.stem + '_' + snap.basefilename) + fn = snapshotdir / Path(snap.filename) snap.filename = fn.with_suffix('.vti') pbar = tqdm(total=snap.vtkdatawritesize, leave=True, unit='byte', - unit_scale=True, desc=f'Writing snapshot file {i + 1} of {len(self.G.snapshots)}, {os.path.split(snap.filename)[1]}', ncols=get_terminal_width() - 1, file=sys.stdout, - disable=not config.general['progressbars']) + unit_scale=True, desc=f'Writing snapshot file {i + 1} of {len(self.G.snapshots)}, {snap.filename.name}', ncols=get_terminal_width() - 1, file=sys.stdout, disable=not config.sim_config.general['progressbars']) snap.write_vtk_imagedata(pbar, self.G) pbar.close() log.info('') @@ -308,7 +293,6 @@ class ModelBuildRun: tsolve (float): time taken to execute solving (seconds). """ - self.create_output_directory() log.info(f'\nOutput file: {config.model_configs[self.G.model_num].output_file_path_ext}') # Check number of OpenMP threads diff --git a/gprMax/pml.py b/gprMax/pml.py index 08cb8408..d948718b 100644 --- a/gprMax/pml.py +++ b/gprMax/pml.py @@ -339,11 +339,13 @@ class CUDAPML(PML): solving on GPU using CUDA. """ - def initialise_field_arrays_gpu(self): + def initialise_field_arrays(self): """Initialise PML field and coefficient arrays on GPU.""" import pycuda.gpuarray as gpuarray + super().initialise_field_arrays() + self.ERA_gpu = gpuarray.to_gpu(self.ERA) self.ERB_gpu = gpuarray.to_gpu(self.ERB) self.ERE_gpu = gpuarray.to_gpu(self.ERE) @@ -352,52 +354,16 @@ class CUDAPML(PML): self.HRB_gpu = gpuarray.to_gpu(self.HRB) self.HRE_gpu = gpuarray.to_gpu(self.HRE) self.HRF_gpu = gpuarray.to_gpu(self.HRF) - - if self.direction[0] == 'x': - self.EPhi1_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx + 1, self.ny, self.nz + 1), - dtype=config.sim_config.dtypes['float_or_double'])) - self.EPhi2_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx + 1, self.ny + 1, self.nz), - dtype=config.sim_config.dtypes['float_or_double'])) - self.HPhi1_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx, self.ny + 1, self.nz), - dtype=config.sim_config.dtypes['float_or_double'])) - self.HPhi2_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx, self.ny, self.nz + 1), - dtype=config.sim_config.dtypes['float_or_double'])) - elif self.direction[0] == 'y': - self.EPhi1_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx, self.ny + 1, self.nz + 1), - dtype=config.sim_config.dtypes['float_or_double'])) - self.EPhi2_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx + 1, self.ny + 1, self.nz), - dtype=config.sim_config.dtypes['float_or_double'])) - self.HPhi1_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx + 1, self.ny, self.nz), - dtype=config.sim_config.dtypes['float_or_double'])) - self.HPhi2_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx, self.ny, self.nz + 1), - dtype=config.sim_config.dtypes['float_or_double'])) - elif self.direction[0] == 'z': - self.EPhi1_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx, self.ny + 1, self.nz + 1), - dtype=config.sim_config.dtypes['float_or_double'])) - self.EPhi2_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx + 1, self.ny, self.nz + 1), - dtype=config.sim_config.dtypes['float_or_double'])) - self.HPhi1_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx + 1, self.ny, self.nz), - dtype=config.sim_config.dtypes['float_or_double'])) - self.HPhi2_gpu = gpuarray.to_gpu(np.zeros((len(self.CFS), - self.nx, self.ny + 1, self.nz), - dtype=config.sim_config.dtypes['float_or_double'])) + self.EPhi1_gpu = gpuarray.to_gpu(self.EPhi1) + self.EPhi2_gpu = gpuarray.to_gpu(self.EPhi2) + self.HPhi1_gpu = gpuarray.to_gpu(self.HPhi1) + self.HPhi2_gpu = gpuarray.to_gpu(self.HPhi2) def set_blocks_per_grid(self, G): """Set the blocks per grid size used for updating the PML field arrays on a GPU. Args: - G (FDTDGrid): Holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. """ self.bpg = (int(np.ceil(((self.EPhi1_gpu.shape[1] + 1) * @@ -422,7 +388,7 @@ class CUDAPML(PML): correction on the GPU. Args: - G (FDTDGrid): Holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. """ self.update_electric_gpu(np.int32(self.xs), np.int32(self.xf), @@ -449,7 +415,7 @@ class CUDAPML(PML): correction on the GPU. Args: - G (FDTDGrid): Holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. """ self.update_magnetic_gpu(np.int32(self.xs), np.int32(self.xf), np.int32(self.ys), np.int32(self.yf), @@ -474,7 +440,7 @@ def pml_information(G): """Information about PMLs. Args: - G (FDTDGrid): Holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. """ # No PML if all(value == 0 for value in G.pmlthickness.values()): @@ -497,7 +463,7 @@ def build_pml(G, key, value): (based on underlying material er and mr from solid array). Args: - G (FDTDGrid): Holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. key (str): Identifier of PML slab. value (int): Thickness of PML slab in cells. """ diff --git a/gprMax/snapshots.py b/gprMax/snapshots.py index b2a42d9a..0a968a3e 100644 --- a/gprMax/snapshots.py +++ b/gprMax/snapshots.py @@ -39,11 +39,6 @@ class Snapshot: # GPU - blocks per grid - set according to largest requested snapshot bpg = None - # Set string for byte order - byteorder = 'LittleEndian' if sys.byteorder == 'little' else 'BigEndian' - - - def __init__(self, xs=None, ys=None, zs=None, xf=None, yf=None, zf=None, dx=None, dy=None, dz=None, time=None, filename=None): """ @@ -54,14 +49,6 @@ class Snapshot: filename (str): Filename to save to. """ - # Set format text and string depending on float type - if config.sim_config.dtypes['float_or_double'] == np.float32: - self.floatname = 'Float32' - self.floatstring = 'f' - elif config.sim_config.dtypes['float_or_double'] == np.float64: - self.floatname = 'Float64' - self.floatstring = 'd' - self.fieldoutputs = {'electric': True, 'magnetic': True} self.xs = xs self.ys = ys @@ -79,21 +66,20 @@ class Snapshot: self.sy = slice(self.ys, self.yf + self.dy, self.dy) self.sz = slice(self.zs, self.zf + self.dz, self.dz) self.ncells = self.nx * self.ny * self.nz - self.datasizefield = (3 * np.dtype(config.dtypes['float_or_double']).itemsize + self.datasizefield = (3 * np.dtype(config.sim_config.dtypes['float_or_double']).itemsize * self.ncells) - self.vtkdatawritesize = ((self.fieldoutputs['electric'] - + self.fieldoutputs['magnetic']) * self.datasizefield - + (self.fieldoutputs['electric'] - + self.fieldoutputs['magnetic']) - * np.dtype(np.uint32).itemsize) + self.vtkdatawritesize = ((self.fieldoutputs['electric'] + + self.fieldoutputs['magnetic']) * + self.datasizefield + (self.fieldoutputs['electric'] + + self.fieldoutputs['magnetic']) * np.dtype(np.uint32).itemsize) self.time = time - self.basefilename = filename + self.filename = filename def store(self, G): """Store (in memory) electric and magnetic field values for snapshot. Args: - G (class): Grid class instance - holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. """ # Memory views of field arrays to dimensions required for the snapshot @@ -105,12 +91,12 @@ class Snapshot: Hzslice = np.ascontiguousarray(G.Hz[self.sx, self.sy, self.sz]) # Create arrays to hold the field data for snapshot - Exsnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.dtypes['float_or_double']) - Eysnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.dtypes['float_or_double']) - Ezsnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.dtypes['float_or_double']) - Hxsnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.dtypes['float_or_double']) - Hysnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.dtypes['float_or_double']) - Hzsnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.dtypes['float_or_double']) + Exsnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.sim_config.dtypes['float_or_double']) + Eysnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.sim_config.dtypes['float_or_double']) + Ezsnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.sim_config.dtypes['float_or_double']) + Hxsnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.sim_config.dtypes['float_or_double']) + Hysnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.sim_config.dtypes['float_or_double']) + Hzsnap = np.zeros((self.nx, self.ny, self.nz), dtype=config.sim_config.dtypes['float_or_double']) # Calculate field values at points (comes from averaging field components in cells) calculate_snapshot_fields( @@ -141,28 +127,28 @@ class Snapshot: Args: pbar (class): Progress bar class instance. - G (class): Grid class instance - holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. """ - hfield_offset = (3 * np.dtype(config.dtypes['float_or_double']).itemsize + hfield_offset = (3 * np.dtype(config.sim_config.dtypes['float_or_double']).itemsize * self.ncells + np.dtype(np.uint32).itemsize) self.filehandle = open(self.filename, 'wb') self.filehandle.write('\n'.encode('utf-8')) - self.filehandle.write(f'\n'.encode('utf-8')) + self.filehandle.write(f'\n'.encode('utf-8')) self.filehandle.write(f'\n'.encode('utf-8')) self.filehandle.write(f'\n'.encode('utf-8')) if self.fieldoutputs['electric'] and self.fieldoutputs['magnetic']: self.filehandle.write('\n'.encode('utf-8')) - self.filehandle.write(f'\n'.encode('utf-8')) - self.filehandle.write(f'\n'.encode('utf-8')) + self.filehandle.write(f"""\n""".encode('utf-8')) + self.filehandle.write(f"""\n""".encode('utf-8')) elif self.fieldoutputs['electric']: self.filehandle.write('\n'.encode('utf-8')) - self.filehandle.write(f'\n'.encode('utf-8')) + self.filehandle.write(f"""\n""".encode('utf-8')) elif self.fieldoutputs['magnetic']: self.filehandle.write('\n'.encode('utf-8')) - self.filehandle.write(f'\n'.encode('utf-8')) + self.filehandle.write(f"""\n""".encode('utf-8')) self.filehandle.write('\n\n\n\n_'.encode('utf-8')) @@ -188,7 +174,7 @@ def initialise_snapshot_array_gpu(G): """Initialise array on GPU for to store field data for snapshots. Args: - G (class): Grid class instance - holds essential parameters describing the model. + G (FDTDGrid): Parameters describing a grid in a model. Returns: snapE_gpu, snapH_gpu (float): numpy arrays of snapshot data on GPU. @@ -206,24 +192,26 @@ def initialise_snapshot_array_gpu(G): Snapshot.nz_max = snap.nz # GPU - blocks per grid - according to largest requested snapshot - Snapshot.bpg = (int(np.ceil(((Snapshot.nx_max) * (Snapshot.ny_max) * (Snapshot.nz_max)) / Snapshot.tpb[0])), 1, 1) + Snapshot.bpg = (int(np.ceil(((Snapshot.nx_max) * + (Snapshot.ny_max) * + (Snapshot.nz_max)) / Snapshot.tpb[0])), 1, 1) # 4D arrays to store snapshots on GPU, e.g. snapEx(time, x, y, z); # if snapshots are not being stored on the GPU during the simulation then # they are copied back to the host after each iteration, hence numsnaps = 1 - numsnaps = 1 if config.cuda['snapsgpu2cpu'] else len(G.snapshots) + numsnaps = 1 if config.model_configs[G.model_num].cuda['snapsgpu2cpu'] else len(G.snapshots) snapEx = np.zeros((numsnaps, Snapshot.nx_max, Snapshot.ny_max, Snapshot.nz_max), - dtype=config.dtypes['float_or_double']) + dtype=config.sim_config.dtypes['float_or_double']) snapEy = np.zeros((numsnaps, Snapshot.nx_max, Snapshot.ny_max, Snapshot.nz_max), - dtype=config.dtypes['float_or_double']) + dtype=config.sim_config.dtypes['float_or_double']) snapEz = np.zeros((numsnaps, Snapshot.nx_max, Snapshot.ny_max, Snapshot.nz_max), - dtype=config.dtypes['float_or_double']) + dtype=config.sim_config.dtypes['float_or_double']) snapHx = np.zeros((numsnaps, Snapshot.nx_max, Snapshot.ny_max, Snapshot.nz_max), - dtype=config.dtypes['float_or_double']) + dtype=config.sim_config.dtypes['float_or_double']) snapHy = np.zeros((numsnaps, Snapshot.nx_max, Snapshot.ny_max, Snapshot.nz_max), - dtype=config.dtypes['float_or_double']) + dtype=config.sim_config.dtypes['float_or_double']) snapHz = np.zeros((numsnaps, Snapshot.nx_max, Snapshot.ny_max, Snapshot.nz_max), - dtype=config.dtypes['float_or_double']) + dtype=config.sim_config.dtypes['float_or_double']) # Copy arrays to GPU snapEx_gpu = gpuarray.to_gpu(snapEx) diff --git a/gprMax/subgrids/updates.py b/gprMax/subgrids/updates.py index ae77c0d4..fb299771 100644 --- a/gprMax/subgrids/updates.py +++ b/gprMax/subgrids/updates.py @@ -75,7 +75,7 @@ class SubgridUpdater(CPUUpdates): subgrid (SubGrid3d): Subgrid to be updated. precursors (PrecursorNodes): Precursor nodes associated with the subgrid - contain interpolated fields. - G (FDTDGrid): Holds essential parameters describing a model. + G (FDTDGrid): Parameters describing a grid in a model. """ super().__init__(subgrid) self.precursors = precursors