Handle NaNs in objective values explicitly

Author: lindonroberts

dfols/model.py

@@ -306,7 +306,7 @@ def solve_geom_system(self, rhs):
             return col_scale(LA.lstsq(W, col_scale(rhs * left_scaling))[0], right_scaling)
 
     def interpolate_mini_models_svd(self, verbose=False, make_full_rank=False, min_sing_val=1e-6, sing_val_frac=1.0, max_jac_cond=1e8,
-                                    get_chg_J=False):
+                                    get_chg_J=False, throw_error_on_nans=False):
         W, left_scaling, right_scaling = self.interpolation_matrix()
         self.factorise_geom_system()
         ls_interp_cond_num = np.linalg.cond(W) if verbose else 0.0  # scipy.linalg does not have condition number!
@@ -327,12 +327,18 @@ def interpolate_mini_models_svd(self, verbose=False, make_full_rank=False, min_s
             self.model_jac = np.dot(self.model_jac, np.dot(Qhat, Qhat.T))
 
         rhs = self.fval_v[fval_row_idx, :]  # size npt * m
+        if np.any(np.isnan(rhs)) and throw_error_on_nans:
+            if self.do_logging:
+                logging.warning("model.interpolate_mini_models_svd: NaNs encountered in objective evaluations, raising error")
+            raise np.linalg.LinAlgError("NaN encountered in objective evaluations")
         try:
             dg = self.solve_geom_system(rhs)  # size (n+1)*m
         except LA.LinAlgError:
             return False, None, None, None, None  # flag error
         except ValueError:
             return False, None, None, None, None  # flag error (e.g. inf or NaN encountered)
+        if not np.all(np.isfinite(dg)):  # another check for inf or NaN
+            return False, None, None, None, None
         J_old = self.model_jac.copy()
         self.model_jac = dg[1:,:].T
         self.model_const = dg[0,:] - np.dot(self.model_jac, xopt)  # shift base to xbase

dfols/params.py

@@ -44,6 +44,7 @@ def __init__(self, n, npt, maxfun, objfun_has_noise=False):
         self.params["init.random_directions_make_orthogonal"] = True  # although random > orthogonal, avoid for init
         # Interpolation
         self.params["interpolation.precondition"] = True
+        self.params["interpolation.throw_error_on_nans"] = False  # throw numpy.linalg.LinAlgError if interpolating to nan data?
         # Logging
         self.params["logging.n_to_print_whole_x_vector"] = 6
         self.params["logging.save_diagnostic_info"] = False
@@ -142,6 +143,8 @@ def param_type(self, key, npt):
             type_str, nonetype_ok, lower, upper = 'bool', False, None, None
         elif key == "interpolation.precondition":
             type_str, nonetype_ok, lower, upper = 'bool', False, None, None
+        elif key == "interpolation.throw_error_on_nans":
+            type_str, nonetype_ok, lower, upper = 'bool', False, None, None
         elif key == "logging.n_to_print_whole_x_vector":
             type_str, nonetype_ok, lower, upper = 'int', False, 0, None
         elif key == "logging.save_diagnostic_info":

dfols/solver.py

@@ -247,7 +247,8 @@ def solve_main(objfun, x0, args, xl, xu, npt, rhobeg, rhoend, maxfun, nruns_so_f
                 min_sing_val=params("growing.full_rank.min_sing_val"),
                 sing_val_frac=params("growing.full_rank.svd_scale_factor"),
                 max_jac_cond=params("growing.full_rank.svd_max_jac_cond"),
-                get_chg_J=params("restarts.use_restarts") and params("restarts.auto_detect"))
+                get_chg_J=params("restarts.use_restarts") and params("restarts.auto_detect"),
+                throw_error_on_nans=params("interpolation.throw_error_on_nans"))
         if not interp_ok:
             if params("restarts.use_restarts") and params("restarts.use_soft_restarts"):
                 number_of_samples = max(nsamples(control.delta, control.rho, current_iter, nruns_so_far), 1)

dfols/tests/test_solver.py

@@ -41,6 +41,20 @@ def rosenbrock_jacobian(x):
     return np.array([[-20.0*x[0], 10.0], [-1.0, 0.0]])
 
 
+class TestNans(unittest.TestCase):
+    # Generic objective that only returns NaNs (like optclim code)
+    # Verify get a sensible termination
+    def runTest(self):
+        x0 = np.array([-1.2, 1.0])
+        r_error = lambda x: np.array([np.nan, np.nan, np.nan])
+        # First attempt: exit gracefully
+        soln = dfols.solve(r_error, x0)
+        self.assertEqual(soln.flag, soln.EXIT_LINALG_ERROR, "Wrong error message")
+        # Second attempt: throw error when trying to interpolate
+        with self.assertRaises(np.linalg.LinAlgError):
+            soln = dfols.solve(r_error, x0, user_params={"interpolation.throw_error_on_nans": True})
+
+
 class TestRosenbrockGeneric(unittest.TestCase):
     # Minimise the (2d) Rosenbrock function
     def runTest(self):