quatica.optiQ

quatica.optiQ

Optimization solvers for quaternion Hermitian matrices (pure quaternion algebra).

This module provides an equality-constrained Newton method for the log-det barrier formulation of quaternionic semidefinite programs (SDP):

minimize - mu * logdet(X) subject to A(X) = b, X ≻ 0, X ∈ H_n(H),

where = Re tr(U^H V) is the real inner product on H_n(H), and A is a linear map with Hermitian measurement matrices {H_i}.

All operations are quaternion-native using QuatIca utilities; no real/complex embeddings are used.

References/Notes: - Gradient/Hessian: ∇ logdet(X) = X^{-1}, ∇² logdet(X)[H] = - X^{-1} H X^{-1} - Barrier Newton KKT via Schur complement with H^{-1}[W] = (1/mu) X W X

BarrierParams dataclass

Parameters controlling the barrier Newton method.

Parameters

mu : float Initial barrier parameter. mu_decay : float Geometric decay factor for mu per outer iteration. mu_min : float Minimum barrier value to stop. newton_tol : float Tolerance for inner Newton convergence (on combined residual). newton_maxit : int Max Newton iterations per barrier value. backtrack_beta : float Backtracking multiplicative shrink. backtrack_sigma : float Armijo parameter for sufficient decrease. posdef_eps : float Minimal eigenvalue threshold to accept PD candidate.

QuaternionSDPOperator dataclass

Linear map A and its adjoint A* represented via Hermitian measurement matrices.

A(X)_i = _R = Re tr(H_i^H X) A*(y) = sum_i y_i H_i

Parameters

H_list : list[np.ndarray] List of Hermitian quaternion matrices (same shape n×n).

QuaternionSDPOperatorScaled dataclass

Scaled constraint operator ilde A with internal per-row scaling s_i.

( ilde A X)_i = / s_i, ilde A^*(y) = \sum_i (y_i / s_i) H_i.

A_unscaled(X)

Return original (unscaled) A(X) with raw H_i.

scale_b(b)

Return b̃ = b / s (elementwise), consistent with internal scaling.

SolverState dataclass

Lightweight state container for solver bookkeeping.

admm_barrier_fixed_mu(H_list, b, C, mu, rho=5.0, maxit=300, abstol=1e-06, reltol=1e-05, X0=None, Z0=None, U0=None, verbose=True, alpha=1.5, use_hat_projection=False)

Description

Scaled ADMM for the fixed-μ barrier subproblem: minimize - μ logdet(X) subject to A(X) = b, X ≻ 0, using a consensus split with Z carrying only the equality constraints.

Parameters

H_list : list[np.ndarray] Hermitian measurement matrices defining A and A*. b : np.ndarray Constraint right-hand side (R^m). C : np.ndarray Hermitian cost matrix. mu : float Barrier parameter (fixed for this solver call). rho : float ADMM penalty parameter. maxit : int Maximum ADMM iterations. abstol, reltol : float Absolute and relative stopping tolerances (Boyd criteria). X0, Z0, U0 : np.ndarray | None Optional initializations for primal/dual variables. verbose : bool If True, prints per-iteration diagnostics.

Returns

tuple (X, Z, U, history) where history contains norms and eigenvalue traces.

Notes/References

• Pure quaternion algebra via QuatIca utilities.
• Z-projection is exact in orthonormalized coordinates  (one line).
• X-update is closed-form SPD prox of -μ logdet with Frobenius metric.

build_canonical_orthonormal_basis(n)

Construct a canonical orthonormal basis of quaternion-Hermitian matrices for H_n(H).

Basis elements (orthonormal under =Re tr(U^H V)): - Diagonals: D_k = E_{kk} - Off-diagonals (i<j): for u in {1, i, j, k} S^{(u)}{ij} = (1/sqrt(2)) ( u E{ij} + conj(u) E_{ji} )

Returns

list[np.ndarray] List of length 2 n^2 - n with dtype quaternion, orthonormal.

build_central_mu_instance(n=3, m=4, mu=1.0, seed=0, basis='canonical')

Build a small, planted central-point instance for the fixed-μ barrier problem.

The planted solution is X_star with y_star=0 and C = μ X_star^{-1}.

Parameters

n : int Matrix size. m : int Number of constraints to use (subset of canonical basis if basis='canonical'). mu : float Barrier parameter used to plant C. seed : int | None RNG seed. basis : str 'canonical' for a subset of the canonical orthonormal basis; otherwise uses random Hermitian constraints orthonormalized via Gram–Cholesky.

Returns

(H_list, b, C, X_star, mu)

build_certificate_optimum(n=4, r=2, seed=0)

Regime B: Choose rank-r X_ and slack S_ on its complement. Set b = A(X_), C = S_, so that (X_, S_) is a certificate at optimum.

Returns

(H_list, b, C, X_star, S_star)

build_singleton_feasible(n=3, seed=0)

Regime A: Complete constraints that pin X uniquely.

Returns

(H_list, b, C, X_star)

build_singleton_feasible_canonical(n=3, m=None, seed=0)

Build a singleton-feasible instance using a canonical orthonormal basis for constraints.

Parameters

n : int Matrix size. m : int | None Number of constraints to use. If None, use full dimension (2 n^2 - n). seed : int | None RNG seed for X_star and for selecting a subset when m < d.

Returns

(H_list, b, C, X_star)

dim_hermitian_quat(n)

Real dimension of H_n(H) = {X ∈ H^{n×n} : X = X^H} → 2 n^2 − n.

eighH(X)

Quaternion Hermitian eigendecomposition: X = V diag(lam) V^H.

Returns

lam : np.ndarray of shape (n,) Real eigenvalues in unspecified order. V : np.ndarray (dtype=quaternion) Unitary quaternion eigenvectors.

eigvalsH(X)

Eigenvalues (real) of a Hermitian quaternion matrix X (with symmetrization).

hess_inv_apply(X, mu, G)

Apply H^{-1}[G] for H[W] = mu * X^{-1} W X^{-1}.

Parameters

X : np.ndarray Current PD iterate. mu : float Barrier parameter. G : np.ndarray Right-hand quaternion matrix.

Returns

np.ndarray (1/mu) X G X

inner_real(U, V)

Real inner product _R = Re tr(U^H V) using quaternion adjoint.

Notes

Must use the adjoint U^H, not the Hermitian part, to preserve adjointness and accurate Gram/Schur assembly.

invH(X)

Inverse of a Hermitian positive definite quaternion matrix.

Parameters

X : np.ndarray Hermitian PD matrix.

Returns

np.ndarray X^{-1} computed via eigen-decomposition in quaternion algebra.

Notes

This routine assumes X is strictly PD; it raises if any eigenvalue ≤ 0.

logdetH(X)

Return sum(log(lam_i)) for lam_i eigenvalues of PD quaternion X.

qadj(A)

Quaternion adjoint (conjugate transpose).

qeye(n)

Return the n×n quaternion identity matrix.

qherm(A)

Hermitian projection: (A + A^H)/2.

qmm(A, B)

Quaternion matrix multiply using QuatIca optimized matmul.

qzeros(n, m)

Return an n×m quaternion zero matrix.

random_H_basis(n, seed=None)

Build a random (numerically independent) set of Hermitian directions H_i. For small n, sampling random Hermitians is sufficient for toy problems.

random_hermitian(n, seed=None)

Random dense quaternion Hermitian matrix.

random_spd(n, lam_min=0.5, lam_max=2.0, seed=None)

Generate a random quaternion SPD matrix via unitary diagonalization.

Parameters

n : int Matrix size. lam_min, lam_max : float Bounds for positive eigenvalues. seed : int | None RNG seed.

Returns

np.ndarray Quaternion SPD matrix.

save_admm_history_plot(hist, out_path)

Description

Save a semilog plot of ADMM residuals across iterations.

Parameters

hist : dict History dictionary from admm_barrier_fixed_mu with keys 'pr','dr','eq'. out_path : str Path to save the figure (PNG). Parent directories will be created.

Notes/References

Figures are saved; not displayed interactively. Use validation_output directory.

save_barrier_history_plot(history, out_path)

Save iteration metrics for primal barrier Newton history.

Plots semilog curves of ||r_p||, ||r_d||, and step t over iterations.

save_pdipm_history_plot(history, out_path)

Save iteration metrics for the barrier-path solver history.

Plots semilog curves of ||r_p||, ||r_d||, gap/n, and step t over iterations.

solve_barrier(A_list, b, C, X0=None, params=BarrierParams(), verbose=True)

Equality-constrained log-det barrier Newton method (quaternion-native), in hat space.

Solves (for each mu stage): minimize _R - mu logdet(X) subject to A_hat(X) = b_hat, X ≻ 0,

with KKT residuals

r_p = b_hat - A_hat(X), r_d = C + A_hat^*(y) - mu X^{-1}.

solve_logdet_barrier_newton(H_list, b, C, *, X0=None, y0=None, mu=1.0, max_iter=50, eps=1e-08, verbose=True, schur_solver='dense', schur_precond='none', schur_precond_rank=None, schur_precond_ridge_scale=1e-06, schur_precond_seed=0, cg_tol=1e-10, cg_maxit=500, ops=None, assume_hat=False, return_ops=False)

Solve the fixed-μ equality-constrained logdet barrier KKT system via Newton steps in hat-space.

solve_logdet_barrier_path(H_list, b, C, *, X0=None, y0=None, mu_init=1.0, beta_mu=0.3, mu_min=1e-08, max_iter=200, eps=1e-08, verbose=True, schur_solver='dense', schur_precond='none', schur_precond_rank=None, schur_precond_ridge_scale=1e-06, schur_precond_seed=0, cg_tol=1e-10, cg_maxit=500, ops=None, assume_hat=False, return_ops=False)

Barrier path-following (μ-continuation) for the logdet barrier SDP subproblem in hat-space.

solve_pd_mehrotra(H_list, b, C, X0=None, S0=None, y0=None, mu_init=1.0, beta_mu=0.3, eps_p=1e-08, eps_d=1e-08, eps_gap=1e-08, max_iter=50, cg_tol=1e-10, cg_maxit=500, verbose=True, fixed_mu=False, *, schur_solver='dense', schur_precond='none', schur_precond_rank=None, schur_precond_ridge_scale=1e-06, schur_precond_seed=0, ops=None, assume_hat=False, return_ops=False)

Barrier-path (μ-continuation) logdet solver in hat space. (Historical name: solve_pd_mehrotra)

Current design (robust core): Solve the fixed-μ barrier KKT system (in hat space): Ahat(X) = b_hat, C + Ahat^*(y) - μ X^{-1} = 0.

If fixed_mu=True: solve only for μ=mu_init.

If fixed_mu=False: continuation on μ: μ <- max(beta_mu*μ, eps_gap), warm-starting from previous (X,y).

sqrtH(X)

Hermitian square root via quaternion eigendecomposition.