diff --git a/gprMax/config.py b/gprMax/config.py
index 2e8c3a4a..c15d1f1f 100644
--- a/gprMax/config.py
+++ b/gprMax/config.py
@@ -30,6 +30,7 @@ from scipy.constants import c
 from scipy.constants import epsilon_0 as e0
 from scipy.constants import mu_0 as m0
 
+from .exceptions import GeneralError
 from .utilities import detect_check_gpus
 from .utilities import get_host_info
 from .utilities import get_terminal_width
@@ -235,6 +236,9 @@ class SimulationConfig:
         if self.general['subgrid']:
             self.general['precision'] = 'double'
 
+        if self.general['subgrid'] and self.general['cuda']:
+            raise GeneralError('The CUDA-based solver cannot currently be used with models that contain sub-grids.')
+
         # Scenes parameter may not exist if user enters via CLI
         try:
             self.scenes = args.scenes
diff --git a/gprMax/contexts.py b/gprMax/contexts.py
index 4cc66833..b1c3d227 100644
--- a/gprMax/contexts.py
+++ b/gprMax/contexts.py
@@ -40,8 +40,8 @@ class Context:
 
     def __init__(self):
         self.model_range = range(config.sim_config.model_start, config.sim_config.model_end)
-        self.tsimend = 0
-        self.tsimstart = 1
+        self.tsimend = None
+        self.tsimstart = None
 
     def run(self):
         """Run the simulation in the correct context."""
@@ -91,13 +91,11 @@ class NoMPIContext(Context):
         for i in self.model_range:
             config.model_num = i
             write_model_config()
-            # Always create a solver for the first model.
-            # The next model to run only gets a new solver if the
-            # geometry is not re-used.
+
+            # Always create a grid for the first model. The next model to run
+            # only gets a new grid if the geometry is not re-used.
             if i != 0 and config.sim_config.args.geometry_fixed:
                 config.get_model_config().reuse_geometry = True
-                # Ensure re-used G is associated correctly with model
-                # G.model_num = i
             else:
                 G = create_G()
 
diff --git a/gprMax/grid.py b/gprMax/grid.py
index bd9d2652..edf2285b 100644
--- a/gprMax/grid.py
+++ b/gprMax/grid.py
@@ -340,21 +340,6 @@ class CUDAGrid(FDTDGrid):
         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()
-
-        # Copy arrays geometry and field component arrays
-        self.htod_geometry_arrays()
-        self.htod_field_arrays()
-        if config.get_model_config().materials['maxpoles'] > 0:
-            self.htod_dispersive_arrays()
-
-        # Copy arrays for fields in PML
-        for pml in self.pmls:
-            pml.initialise_field_arrays_gpu()
-
 
 def dispersion_analysis(G):
     """Analysis of numerical dispersion (Taflove et al, 2005, p112) -
diff --git a/gprMax/model_build_run.py b/gprMax/model_build_run.py
index a26b8d44..0de411a8 100644
--- a/gprMax/model_build_run.py
+++ b/gprMax/model_build_run.py
@@ -171,19 +171,20 @@ class ModelBuildRun:
             # Check to see if numerical dispersion might be a problem
             results = dispersion_analysis(gb.grid)
             if results['error']:
-                log.warning(Fore.RED + f"\nNumerical dispersion analysis ({gb.grid.name}) not carried out as {results['error']}" + Style.RESET_ALL)
+                log.warning(Fore.RED + f"\nNumerical dispersion analysis on {gb.grid.name} not carried out as {results['error']}" + Style.RESET_ALL)
             elif results['N'] < config.get_model_config().numdispersion['mingridsampling']:
                 raise GeneralError(f"\nNon-physical wave propagation in {gb.grid.name} detected. Material '{results['material'].ID}' has wavelength sampled by {results['N']} cells, less than required minimum for physical wave propagation. Maximum significant frequency estimated as {results['maxfreq']:g}Hz")
             elif (results['deltavp'] and np.abs(results['deltavp']) >
                   config.get_model_config().numdispersion['maxnumericaldisp']):
                 log.warning(Fore.RED + f"\n{gb.grid.name} has potentially significant numerical dispersion. Estimated largest physical phase-velocity error is {results['deltavp']:.2f}% in material '{results['material'].ID}' whose wavelength sampled by {results['N']} cells. Maximum significant frequency estimated as {results['maxfreq']:g}Hz" + Style.RESET_ALL)
             elif results['deltavp']:
-                log.info(f"\nNumerical dispersion analysis ({gb.grid.name}): estimated largest physical phase-velocity error is {results['deltavp']:.2f}% in material '{results['material'].ID}' whose wavelength sampled by {results['N']} cells. Maximum significant frequency estimated as {results['maxfreq']:g}Hz")
+                log.info(f"\nNumerical dispersion analysis on {gb.grid.name}: estimated largest physical phase-velocity error is {results['deltavp']:.2f}% in material '{results['material'].ID}' whose wavelength sampled by {results['N']} cells. Maximum significant frequency estimated as {results['maxfreq']:g}Hz")
 
     def reuse_geometry(self):
         # Reset iteration number
         self.G.iteration = 0
-        config.get_model_config().inputfilestr = f'\n--- Model {config.get_model_config().appendmodelnumber}/{config.sim_config.model_end}, input file (not re-processed, i.e. geometry fixed): {config.sim_config.input_file_path}'
+        s = f'\n--- Model {config.get_model_config().appendmodelnumber}/{config.sim_config.model_end}, input file (not re-processed, i.e. geometry fixed): {config.sim_config.input_file_path}'
+        config.get_model_config().inputfilestr = Fore.GREEN + f"{s} {'-' * (get_terminal_width() - 1 - len(s))}\n" + Style.RESET_ALL
         log.info(config.get_model_config().inputfilestr)
         for grid in [self.G] + self.G.subgrids:
             grid.reset_fields()
@@ -286,7 +287,7 @@ class GridBuilder:
 
     def build_pmls(self):
         pbar = tqdm(total=sum(1 for value in self.grid.pmlthickness.values() if value > 0),
-                    desc=f'Building PML boundaries ({self.grid.name})',
+                    desc=f'Building PML boundaries [{self.grid.name}]',
                     ncols=get_terminal_width() - 1, file=sys.stdout,
                     disable=not config.sim_config.general['progressbars'])
         for pml_id, thickness in self.grid.pmlthickness.items():
@@ -299,7 +300,7 @@ class GridBuilder:
         # Build the model, i.e. set the material properties (ID) for every edge
         # of every Yee cell
         log.info('')
-        pbar = tqdm(total=2, desc=f'Building Yee cells ({self.grid.name})',
+        pbar = tqdm(total=2, desc=f'Building Yee cells [{self.grid.name}]',
                     ncols=get_terminal_width() - 1, file=sys.stdout,
                     disable=not config.sim_config.general['progressbars'])
         build_electric_components(self.grid.solid, self.grid.rigidE, self.grid.ID, self.grid)
@@ -330,5 +331,5 @@ class GridBuilder:
         materialstable.outer_border = False
         materialstable.justify_columns[0] = 'right'
 
-        log.info(f'\nMaterials ({self.grid.name}):')
+        log.info(f'\nMaterials [{self.grid.name}]:')
         log.info(materialstable.table)