minor changes, QR presolve no, adjusted rtol;

plan; take current code prototype and refactor it in one
big Solver class; all memory allocation is done in init;
all methods call simple functions, many already written, that
do not allocate any memory; easy to replace them with C equivalents
for testing, so we can start plugging those in; possible improvements
in my opinion from preconditioning (exact diag precond can ~half the iters)
most probably diag plus rank-1, and from better line search, thinking to
to 1-dim Newton; Minamide is too complicated; QR makes code too heavy;
refinement is not usable unless active set already guessed, can only be
used as polishing in the end; infeas/unbound logic will all have to be designed;
other cones: only SOC for the foreseeable future, main addition with refactoring;
then take Cvxportfolio problems as benchmark for optimizing solver parameters

project_euromir/direction_calculator.py

@@ -251,8 +251,8 @@ def get_direction(
         diag_H = np.array(
             np.diag(_densify(self._hessian_function(current_point))))
 
-        diag_H[diag_H < 1e-4] = 1.
-        self._preconditioner = sp.diags(1./diag_H)
+        diag_H += 1e-6 #[diag_H < 1e-12] = 1.
+        self._preconditioner = sp.sparse.diags(1./diag_H)
         return super().get_direction(
             current_point=current_point, current_gradient=current_gradient)
 

project_euromir/loss_no_hsde.py

@@ -71,7 +71,9 @@ def dual_project(vec, result, zero, nonneg, soc):
     result[zero+nonneg:zero+nonneg+soc_cursor] *= -1.
 
 def make_derivative_project(vec, result, nonneg, soc):
-    """Make DPi_K. We skip zero cone.
+    """Make DPi_K.
+
+    We skip zero cone.
     """
     y_mask = np.ones(nonneg, dtype=float)
     y_mask[:] = result < 0.
@@ -80,8 +82,7 @@ def operator(dy):
     return operator
 
 def make_derivative_dual_project(vec, result, zero, nonneg, soc):
-    """Make DPi_K*.
-    """
+    """Make DPi_K*."""
     s_mask = np.ones(zero+nonneg, dtype=float)
     s_mask[zero:] = result[zero:] < 0.
     def operator(ds):

project_euromir/solver_new.py

@@ -42,7 +42,7 @@
 
 logger = logging.getLogger(__name__)
 
-QR_PRESOLVE = True
+QR_PRESOLVE = False
 
 def solve(matrix, b, c, zero, nonneg, soc=(),
         # xy = None, # need to import logic for equilibration
@@ -120,7 +120,8 @@ def _local_derivative_residual(xy):
 
     # line_searcher = LogSpaceLineSearcher(
     #     loss_function=_local_loss,
-    #     min_step=1e-12,
+    #     min_step=1e-16,
+    #     grid_len=1000,
     #     #gradient_function=_local_grad,
     #     #c_1=1e-4,
     #     #c_2=0.9,
@@ -133,14 +134,14 @@ def _local_derivative_residual(xy):
         loss_function=_local_loss,
         max_iters=1000)
 
-    direction_calculator = WarmStartedCGNewton(
-        # warm start causes issues if null space changes b/w iterations
-        hessian_function=_local_hessian,
-        rtol_termination=lambda x, g: min(0.5, np.linalg.norm(g)**0.5),
-        max_cg_iters=None,
-        minres=False,
-        regularizer=1e-8, # it seems 1e-10 is best, but it's too sensitive to it :(
-        )
+    # direction_calculator = WarmStartedCGNewton(
+    #     # warm start causes issues if null space changes b/w iterations
+    #     hessian_function=_local_hessian,
+    #     rtol_termination=lambda x, g: min(0.5, np.linalg.norm(g)**0.5),
+    #     max_cg_iters=None,
+    #     minres=False,
+    #     regularizer=1e-8, # it seems 1e-10 is best, but it's too sensitive to it :(
+    #     )
 
     # doesn't improve, yet; we can make many tests
     # direction_calculator = DiagPreconditionedCGNewton(
@@ -160,7 +161,7 @@ def _local_derivative_residual(xy):
     direction_calculator  = CGNewton(
         # warm start causes issues if null space changes b/w iterations
         hessian_function=_local_hessian,
-        rtol_termination=lambda x, g: min(0.5, np.linalg.norm(g)), #,**2),
+        rtol_termination=lambda x, g: min(0.5, np.linalg.norm(g)**0.5), #,**2),
         max_cg_iters=None,
         # minres=True, # less stable, needs more stringent termination,
         # and/or better logic to transition to refinement, but it is a bit faster