diff --git a/CHANGELOG.md b/CHANGELOG.md index d982e56..47d30a9 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,7 +4,7 @@ All notable changes to this project will be documented in this file. The format is based on [Keep a Changelog](http://keepachangelog.com/) and this project adheres to [Semantic Versioning](http://semver.org/). -## [Unreleased] - 2025-09-01 +## [Unreleased] ### Added @@ -12,6 +12,17 @@ and this project adheres to [Semantic Versioning](http://semver.org/). ### Fixed +## [0.0.7] - 2025-09-09 + +### Added +- linalg_tools.extract_lower_rank_solution: new algorithm to extract + lower-rank solution from a solution to an SDP. + +### Changed +- sdp_solvers.solve_low_rank_sdp: now returns the actuall success status of IPOPT. + +### Fixed + ## [0.0.6] - 2025-09-01 ### Added diff --git a/cert_tools/__init__.py b/cert_tools/__init__.py index 7e0abbd..fe2a767 100644 --- a/cert_tools/__init__.py +++ b/cert_tools/__init__.py @@ -4,4 +4,4 @@ from .linalg_tools import * from .sdp_solvers import * -__version__ = "0.0.6" +__version__ = "0.0.7" diff --git a/cert_tools/linalg_tools.py b/cert_tools/linalg_tools.py index 6d93bb1..2bee1dd 100644 --- a/cert_tools/linalg_tools.py +++ b/cert_tools/linalg_tools.py @@ -1,3 +1,4 @@ +import warnings from copy import deepcopy import numpy as np @@ -84,7 +85,7 @@ def project_so3(X): return U @ Vh -def rank_project(X, p=1, tolerance=1e-10): +def rank_project(X, p=1, tolerance=1e-10) -> tuple[np.ndarray, dict]: """Project matrix X to matrix of rank p.""" try: assert la.issymmetric(X, atol=tolerance) @@ -113,7 +114,52 @@ def rank_project(X, p=1, tolerance=1e-10): "error eigs": np.sum(np.abs(E[p:])), "EVR": abs(E[p - 1] / E[p]), # largest over second-largest } - return x, info + return np.asarray(x), info + + +def extract_lower_rank_solution(X, A_list, tol=1e-8, max_trials=10): + """Given a psd matrix X of rank r, returns a matrix X' of rank r-1. + + This is an implementation of Algorithm 2.1 by [1]. + + [1] Lemon et al. -- Low-Rank Semidefinite Programming Theory and Applications + Foundations and Trends in Optimization, 2015. + """ + r_start = 0 + for k in range(max_trials): + E, U = np.linalg.eigh(X) + r = int(np.sum(E > tol)) + if r == 1: + return X, {"success": True, "rank": 1, "iter": k} + elif r == r_start - 1: + return X, {"success": True, "rank": r, "iter": k} + if r_start == 0: + r_start = r + V = U[:, -r:] @ np.diag(np.sqrt(E[-r:])) + np.testing.assert_allclose(X, V @ V.T, atol=tol) + + AV = np.vstack([svec(V.T @ Ai @ V) for Ai in A_list]) # type: ignore + + # below, the rows of bi have the nullspace vectors. + # AV @ bi = 0 + b, info_null = get_nullspace(AV, tolerance=tol) + if b.shape[0] == 0: + return X, {"success": True, "rank": r, "iter": k} + break + + # chose one of the vectors + idx = np.random.choice(b.shape[0]) + bi = b[idx] + + assert np.all(AV @ bi < tol) + Delta = smat(bi) + + lambdas = np.abs(np.linalg.eigvalsh(Delta)) + # maximum magnitude eigenvalue + alpha = -1 / np.max(lambdas) + X = V @ (np.eye(r) + alpha * Delta) @ V.T + warnings.warn(f"Did not find a lower-rank solution in {max_trials} trials.") + return X, {"success": True, "rank": r, "iter": k} def find_dependent_columns(A_sparse, tolerance=1e-10, verbose=False, debug=False): diff --git a/cert_tools/sdp_solvers.py b/cert_tools/sdp_solvers.py index c04a352..8fecfef 100644 --- a/cert_tools/sdp_solvers.py +++ b/cert_tools/sdp_solvers.py @@ -170,13 +170,21 @@ def solve_low_rank_sdp( if x_cand is not None: sol_input["x0"] = x_cand.reshape((-1, 1)) r = S(**sol_input) - Y_opt = r["x"] + + # Get solver info + stats = S.stats() + success = stats.get("success", False) + info = {"solver_stats": stats, "success": success} + # Reshape and generate SDP solution + Y_opt = r["x"] Y_opt = np.array(Y_opt).reshape((n, rank), order="F") X_opt = Y_opt @ Y_opt.T + # Get cost scale, offset = adjust cost = np.trace(Q @ X_opt) * scale + offset + # Construct certificate mults = np.array(r["lam_g"]) H = Q @@ -186,7 +194,7 @@ def solve_low_rank_sdp( H = H + A * mults[i, 0] # Return - info = {"X": X_opt, "H": H, "cost": cost} + info = {"X": X_opt, "H": H, "cost": cost, "success": success} return Y_opt, info diff --git a/setup.cfg b/setup.cfg index 81583ce..e83fb68 100644 --- a/setup.cfg +++ b/setup.cfg @@ -1,6 +1,6 @@ [metadata] name = cert_tools -version = 0.0.6 +version = 0.0.7 authors = [ {name = "Frederike Dümbgen", email = "frederike.dumbgen@utoronto.ca" }, {name = "Connor Holmes", email = "connor.holmes@mail.utoronto.ca" }]