diff --git a/src/prog_models/prognostics_model.py b/src/prog_models/prognostics_model.py
index f594dfbef..ee2a3410b 100644
--- a/src/prog_models/prognostics_model.py
+++ b/src/prog_models/prognostics_model.py
@@ -771,7 +771,7 @@ def simulate_to_threshold(self, future_loading_eqn : Callable, first_output : di
         t = config['t0']
         u = future_loading_eqn(t)
         if 'x' in config:
-            x = config['x']
+            x = deepcopy(config['x'])
         else:
             x = self.initialize(u, first_output)
         
@@ -1160,10 +1160,10 @@ def generate_surrogate(self, load_functions, method = 'dmd', **kwargs):
         config.update(kwargs)
 
         # List of user-define values to include in surrogate model: 
-        states_dmd = deepcopy(self.states)
-        inputs_dmd = deepcopy(self.inputs)
-        outputs_dmd = deepcopy(self.outputs)
-        events_dmd = deepcopy(self.events)
+        states_dmd = self.states.copy()
+        inputs_dmd = self.inputs.copy()
+        outputs_dmd = self.outputs.copy()
+        events_dmd = self.events.copy()
 
         # Validate user inputs 
         try:
@@ -1409,17 +1409,17 @@ def next_state(self, x, u, _):
 
             def simulate_to_threshold(self, future_loading_eqn, first_output = None, threshold_keys = None, **kwargs):
                 # Save keyword arguments same as DMD training for approximation 
-                kwargs_temp = {
-                    'save_freq': dmd_dt,
-                    'dt': dmd_dt
-                }
                 kwargs_sim = kwargs.copy()
-                kwargs_sim.update(kwargs_temp)
+                kwargs_sim['save_freq'] = dmd_dt
+                kwargs_sim['dt'] = dmd_dt
 
                 # Simulate to threshold at DMD time step
                 results = super().simulate_to_threshold(future_loading_eqn,first_output, threshold_keys, **kwargs_sim)
                 
                 # Interpolate results to be at user-desired time step
+                if 'dt' in kwargs:
+                    warn("dt is not used in DMD approximation")
+
                 # Default parameters 
                 config = {
                     'dt': None,
@@ -1427,20 +1427,29 @@ def simulate_to_threshold(self, future_loading_eqn, first_output = None, thresho
                     'save_pts': []
                 }
                 config.update(kwargs)
-                if config['dt'] != None:
-                    warn("dt is not used in DMD approximation")
 
-                if config['save_freq'] == dmd_dt and config['save_pts'] == []: 
+                if (config['save_freq'] == dmd_dt or
+                    (isinstance(config['save_freq'], tuple) and
+                        config['save_freq'][0]%dmd_dt < 1e-9 and
+                        config['save_freq'][1] == dmd_dt)
+                    ) and config['save_pts'] == []:
                     # In this case, the user wants what the DMD approximation returns 
                     return results 
-                
+
                 # In this case, the user wants something different than what the DMD approximation retuns, so we must interpolate 
                 # Define time vector based on user specifications
                 time_basic = [results.times[0], results.times[-1]]
                 time_basic.extend(config['save_pts'])                       
                 if config['save_freq'] != None:
-                    # Add Save Frequency
-                    time_array = np.arange(results.times[0]+config['save_freq'],results.times[-1],config['save_freq'])
+                    if isinstance(config['save_freq'], tuple):
+                        # Tuple used to specify start and frequency
+                        t_step = config['save_freq'][1]
+                        # Use starting time or the next multiple
+                        t_start = config['save_freq'][0]
+                        start = max(t_start, results.times[0] - (results.times[0]-t_start)%t_step)
+                        time_array = np.arange(start+t_step,results.times[-1],t_step)
+                    else: 
+                        time_array = np.arange(results.times[0]+config['save_freq'],results.times[-1],config['save_freq'])
                     time_basic.extend(time_array.tolist())
                 time_interp = sorted(time_basic)
 
diff --git a/src/prog_models/sim_result.py b/src/prog_models/sim_result.py
index 77ff7db3e..26b133910 100644
--- a/src/prog_models/sim_result.py
+++ b/src/prog_models/sim_result.py
@@ -20,7 +20,7 @@ class SimResult(UserList):
     __slots__ = ['times', 'data']  # Optimization 
     
     def __init__(self, times : list = [], data : list = []):
-        self.times = deepcopy(times)
+        self.times = times.copy()
         self.data = deepcopy(data)
 
     def __eq__(self, other : "SimResult") -> bool:
@@ -174,7 +174,7 @@ def __init__(self, fcn : Callable, times : list = [], states : list = []) -> Non
             data (array(dict)): Data points where data[n] corresponds to times[n]
         """
         self.fcn = fcn
-        self.times = deepcopy(times)
+        self.times = times.copy()
         self.states = deepcopy(states)
         self.__data = None
 
@@ -207,7 +207,7 @@ def extend(self, other : "LazySimResult") -> None:
 
         """
         if (isinstance(other, self.__class__)):
-            self.times.extend(deepcopy(other.times))  # lgtm [py/modification-of-default-value]
+            self.times.extend(other.times)  # lgtm [py/modification-of-default-value]
             self.states.extend(deepcopy(other.states))  # lgtm [py/modification-of-default-value]
             if self.__data is None or not other.is_cached():
                 self.__data = None
diff --git a/src/prog_models/utils/containers.py b/src/prog_models/utils/containers.py
index a732231b5..324dbc73d 100644
--- a/src/prog_models/utils/containers.py
+++ b/src/prog_models/utils/containers.py
@@ -52,7 +52,9 @@ def __len__(self) -> int:
         return len(self._keys)
 
     def __eq__(self, other : "DictLikeMatrixWrapper") -> bool:
-        return self._keys == other._keys and (self.matrix == other.matrix).all()
+        if isinstance(other, dict):
+            return list(self.keys()) == list(other.keys()) and (self.matrix == np.array([[other[key]] for key in self._keys])).all()
+        return self.keys() == other.keys() and (self.matrix == other.matrix).all()
 
     def __hash__(self):
         return hash(self.keys) + hash(self.matrix)
@@ -60,6 +62,9 @@ def __hash__(self):
     def __str__(self) -> str:
         return self.__repr__()
 
+    def copy(self) -> "DictLikeMatrixWrapper":
+        return DictLikeMatrixWrapper(self._keys, self.matrix.copy())
+
     def keys(self) -> list:
         return self._keys