quatica.utils
quatica.utils
A2A0123(A)
Extract component matrices from real matrix A = [A0 A2 A1 A3].
Parameters:
A : numpy.ndarray Real matrix with columns arranged as [A0 A2 A1 A3] where each component matrix has n columns
Returns:
tuple : (A0, A1, A2, A3) component matrices
Notes:
Based on MATLAB implementation by Zhigang Jia (Aug 14, 2014)
GRSGivens(g1, g2=None, g3=None, g4=None)
Generate quaternion Givens rotation matrix for the last column.
Parameters:
g1 : float or numpy.ndarray First quaternion component or 4-vector g2, g3, g4 : float, optional Additional quaternion components (only used if g1 is scalar)
Returns:
numpy.ndarray : 4x4 real Givens rotation matrix
Notes:
Based on MATLAB implementation by Zhigang Jia
Hess_QR_ggivens(Hess)
QR decomposition of quaternion Hessenberg matrix using Givens rotations.
Parameters:
Hess : numpy.ndarray Real block matrix representation of quaternion Hessenberg matrix Hess = [A0; A1; A2; A3] where A0, A1, A2, A3 are the quaternion components
Returns:
tuple : (W, Hess) where W*Hess gives the QR factorization
Notes:
Based on MATLAB implementation by Zhigang Jia (Apr 21, 2020)
Realp(A1, A2, A3, A4)
Convert quaternion matrix components to real block matrix representation.
Parameters:
A1, A2, A3, A4 : numpy.ndarray or scalar Quaternion matrix components (can be scalars or matrices)
Returns:
numpy.ndarray : Real block matrix AR = [A1 -A2 -A3 -A4; A2 A1 -A4 A3; A3 A4 A1 -A2; A4 -A3 A2 A1]
Notes:
Based on MATLAB implementation by Zhigang Jia
UtriangleQsparse(R0, R1, R2, R3, b0, b1, b2, b3, tol=1e-14)
Solve quaternion upper triangular system R * x = b using backward substitution.
Parameters:
R0, R1, R2, R3 : numpy.ndarray Upper triangular matrix components b0, b1, b2, b3 : numpy.ndarray Right-hand side vector components tol : float, optional Tolerance for checking zero elements (default: 1e-14)
Returns:
tuple : (b0, b1, b2, b3) - Solution vector components (overwrites input b)
Notes:
Based on MATLAB implementation by Zhigang Jia & Xuan Liu (Apr 22, 2020) Implements backward substitution: Algorithm 3.1.2, page 89, Matrix Computations, Golub and Van Loan, 3rd ed.
absQsparse(A0, A1, A2, A3)
Compute absolute value/norm of quaternion in component format.
Parameters:
A0, A1, A2, A3 : float or numpy.ndarray Quaternion components
Returns:
tuple : (r, s0, s1, s2, s3) where r is the norm and s0,s1,s2,s3 are normalized components
Notes:
Based on MATLAB implementation by Zhigang Jia
compute_real_svd_pinv(X_real)
Compute pseudoinverse using SVD in the real domain
det(X, d)
Compute determinant of a quaternion matrix.
Parameters:
X : numpy.ndarray with dtype=quaternion Square quaternion matrix d : str Determinant type: - 'Moore': Product of eigenvalues (requires Hermitian matrix) - 'Dieudonné' or 'Dieudonne': Product of singular values - 'Study': Determinant of the adjoint matrix
Returns:
complex or float : The computed determinant
Notes:
- Moore determinant can be negative or complex, but requires Hermitian matrix
- Dieudonné determinant is always real
- Study determinant is the square of Dieudonné determinant
dotinvQsparse(A0, A1, A2, A3)
Compute quaternion inverse in component format.
Parameters:
A0, A1, A2, A3 : float or numpy.ndarray Quaternion components
Returns:
tuple : (inv0, inv1, inv2, inv3) - Inverse quaternion components
Notes:
Based on MATLAB implementation by Zhigang Jia
ggivens(x1, x2)
Generate quaternion Givens rotation matrix.
Parameters:
x1, x2 : numpy.ndarray Quaternion vectors (4-component arrays)
Returns:
numpy.ndarray : 8x8 real Givens rotation matrix G such that G'*[x1;x2] = [||x||_2; 0]
Notes:
Based on MATLAB implementation by Zhigang Jia, Musheng Wei, Meixiang Zhao and Yong Chen (Mar 24, 2017)
induced_matrix_norm_1(A)
Matrix 1-norm induced by vector 1-norm: max column sum of |A_ij|.
Notes
- Supports dense quaternion ndarrays. For sparse, convert to dense first or use component-space helpers.
induced_matrix_norm_inf(A)
Matrix infinity-norm induced by vector infinity-norm: max row sum of |A_ij|.
Notes
- Supports dense quaternion ndarrays. For sparse, convert to dense first or use component-space helpers.
ishermitian(A, tol=None)
Check if a quaternion matrix is Hermitian to within the given tolerance.
Parameters:
A : numpy.ndarray with dtype=quaternion Quaternion matrix to test tol : float, optional Tolerance for comparison (default: machine epsilon)
Returns:
bool : True if matrix is Hermitian within tolerance
Notes:
A matrix A is Hermitian if A = A^H where A^H is the conjugate transpose.
matrix_norm(A, ord=None)
Compute common matrix norms for quaternion matrices.
ord supported
- None or 'fro' or 'F': Frobenius norm
- 1: Induced 1-norm (max column sum)
- np.inf or 'inf': Induced infinity-norm (max row sum)
- 2: Spectral norm (largest singular value)
normQ(A, opt=None)
Compute norm of quaternion matrix A.
Parameters:
A : numpy.ndarray with dtype=quaternion Quaternion matrix opt : str, optional Norm type: 'd' (dual) or None (Frobenius)
Returns:
float : The computed norm
Notes:
Based on MATLAB implementation by Zhigang Jia (Jan 24, 2018)
normQsparse(A0, A1, A2, A3, opt=None)
Compute norm of quaternion matrix in component format (A0, A1, A2, A3).
Parameters:
A0, A1, A2, A3 : numpy.ndarray or scipy.sparse matrix Quaternion matrix components opt : str, optional Norm type: 'd' (dual), '2' (2-norm), '1' (1-norm), or None (Frobenius)
Returns:
float : The computed norm
Notes:
Based on MATLAB implementation by Zhigang Jia (Apr 27, 2020)
power_iteration(A, max_iterations=100, tol=1e-10, return_eigenvalue=False, verbose=False)
Compute the dominant eigenvector of a quaternion matrix using power iteration.
Parameters:
A : numpy.ndarray with dtype=quaternion Square quaternion matrix of size n×n max_iterations : int, optional Maximum number of iterations (default: 100) tol : float, optional Convergence tolerance for vector norm difference (default: 1e-10) return_eigenvalue : bool, optional If True, also return an eigenvalue estimate (default: False) verbose : bool, optional If True, print convergence information and non-Hermitian warning (default: False)
Returns:
numpy.ndarray or tuple: - If return_eigenvalue=False: dominant eigenvector (n×1 quaternion vector) - If return_eigenvalue=True: (eigenvector, eigenvalue_estimate) tuple
Notes:
- Intended primarily for HERMITIAN quaternion matrices, where eigenvalues are real. In that case this routine returns a meaningful dominant eigenpair.
- For general (non-Hermitian) quaternion matrices, the returned scalar when
return_eigenvalue=True is a magnitude-based Rayleigh-quotient estimate (real),
which can be interpreted as a real-part/magnitude heuristic, not a true complex
eigenvalue. For non-Hermitian matrices, use
power_iteration_nonhermitianwhich returns complex eigenvalues (in a fixed complex subfield) and quaternion eigenvectors.
Algorithm: 1) Start with a random vector 2) Iterate v_{k+1} = A v_k / ||A v_k|| 3) Stop when ||v_{k+1} - v_k|| < tol or max_iterations reached
power_iteration_nonhermitian(A, max_iterations=5000, eig_tol=1e-12, res_tol=1e-10, seed=0, return_vector=True, eigenvalue_format='complex', subfield_axis='x', block_purify=True)
Power iteration for NON-Hermitian quaternion matrices (experimental).
Returns a complex eigenvalue estimate (in a chosen complex subfield) and a corresponding quaternion eigenvector approximation.
Parameters
A : np.ndarray Square quaternion matrix (non-Hermitian allowed) max_iterations : int Maximum iterations for the complex power iteration eig_tol : float Tolerance for eigenvalue stabilization res_tol : float | None Residual tolerance on ||Mv - lambda v|| to declare convergence seed : int RNG seed for initialization return_vector : bool If True, also return the quaternion eigenvector approximation eigenvalue_format : {"complex", "quaternion"} Output format for the eigenvalue. Default is "complex". subfield_axis : {"x"} Imaginary axis defining complex subfield for the adjoint (currently only "x").
Notes
- For Hermitian matrices, prefer
power_iteration, which returns a quaternion eigenvector and a real-valued eigenvalue magnitude. - Eigenvalues appear in conjugate pairs in the adjoint; this returns one complex root.
quat_abs_scalar(q)
Return the modulus |q| of a quaternion scalar q.
Computes the absolute value (modulus) of a quaternion q = w + xi + yj + zk as |q| = sqrt(w² + x² + y² + z²).
Parameters:
q : quaternion.quaternion Input quaternion scalar
Returns:
float Modulus (absolute value) of the quaternion
Notes:
The quaternion modulus satisfies |q₁ * q₂| = |q₁| * |q₂| and is used in defining norms and distances in quaternion space.
quat_eye(n)
Create an n×n identity quaternion matrix.
Generates the quaternion identity matrix I where I_ij = δ_ij * (1 + 0i + 0j + 0k), i.e., ones on the diagonal and zeros elsewhere.
Parameters:
n : int Size of the square identity matrix
Returns:
np.ndarray An n×n quaternion identity matrix
Notes:
The quaternion identity matrix satisfies A @ I = I @ A = A for any n×n quaternion matrix A.
quat_frobenius_norm(A)
Compute the Frobenius norm of a quaternion matrix (dense or sparse).
Calculates ||A||_F = sqrt(sum(|A_ij|^2)) where |A_ij| is the modulus of the quaternion at position (i,j).
Parameters:
A : np.ndarray or SparseQuaternionMatrix Input quaternion matrix
Returns:
float Frobenius norm of the matrix
Notes:
For sparse matrices, the norm is computed efficiently by summing the squared norms of each component separately.
quat_hermitian(A)
Return the conjugate transpose (Hermitian) of quaternion matrix A (dense or sparse).
Computes A^H = (A)^T where A is the complex conjugate and T is transpose. For quaternions q = w + xi + yj + zk, the conjugate is q* = w - xi - yj - zk.
Parameters:
A : np.ndarray or SparseQuaternionMatrix Input quaternion matrix
Returns:
np.ndarray or SparseQuaternionMatrix Conjugate transpose A^H of the input matrix
Notes:
The Hermitian (conjugate transpose) is fundamental in quaternion linear algebra and appears in definitions of unitary matrices, eigenvalue problems, and norms.
quat_kernel(A, side='right', rtol=1e-10)
Compute the kernel (null space) of a quaternion matrix.
Alias for quat_null_space() with identical functionality.
Parameters:
A : numpy.ndarray with dtype=quaternion Input quaternion matrix of shape (m, n) side : str, optional 'right' for ker(A), 'left' for ker(A^H) Default: 'right' rtol : float, optional Relative tolerance for determining rank Default: 1e-10
Returns:
N : numpy.ndarray with dtype=quaternion Kernel matrix such that A @ N ≈ 0 (right) or N^H @ A ≈ 0 (left)
quat_matmat(A, B)
Multiply two quaternion matrices (supports dense × dense, sparse × dense, dense × sparse, sparse × sparse).
Performs quaternion matrix multiplication using optimized algorithms based on the input types. Handles mixed sparse/dense operations efficiently.
Parameters:
A : np.ndarray or SparseQuaternionMatrix First quaternion matrix B : np.ndarray or SparseQuaternionMatrix Second quaternion matrix
Returns:
np.ndarray or SparseQuaternionMatrix Result of quaternion matrix multiplication A @ B
Notes:
The function automatically selects the appropriate multiplication algorithm: - Dense × Dense: Component-wise quaternion multiplication - Sparse × Dense/Sparse: Uses sparse matrix multiplication routines - Dense × Sparse: Uses left multiplication method
quat_null_left(A, rtol=1e-10)
Compute the left null space of a quaternion matrix: null(A^H).
Parameters:
A : numpy.ndarray with dtype=quaternion Input quaternion matrix of shape (m, n) rtol : float, optional Relative tolerance for determining rank Default: 1e-10
Returns:
N : numpy.ndarray with dtype=quaternion Left null space matrix of shape (m, m-rank) such that N^H @ A ≈ 0
Notes:
Convenience function equivalent to quat_null_space(A, side='left', rtol=rtol)
quat_null_right(A, rtol=1e-10)
Compute the right null space of a quaternion matrix: null(A).
Parameters:
A : numpy.ndarray with dtype=quaternion Input quaternion matrix of shape (m, n) rtol : float, optional Relative tolerance for determining rank Default: 1e-10
Returns:
N : numpy.ndarray with dtype=quaternion Right null space matrix of shape (n, n-rank) such that A @ N ≈ 0
Notes:
Convenience function equivalent to quat_null_space(A, side='right', rtol=rtol)
quat_null_space(A, side='right', rtol=1e-10)
Compute the null space (kernel) of a quaternion matrix using Q-SVD.
Parameters:
A : numpy.ndarray with dtype=quaternion Input quaternion matrix of shape (m, n) side : str, optional 'right' for right null space (null(A)), 'left' for left null space (null(A^H)) Default: 'right' rtol : float, optional Relative tolerance for determining rank (singular values <= rtol * max(s) are zero) Default: 1e-10
Returns:
N : numpy.ndarray with dtype=quaternion For side='right': null space matrix of shape (n, n-rank) such that A @ N ≈ 0 For side='left': null space matrix of shape (m, m-rank) such that N^H @ A ≈ 0
Notes:
Uses Q-SVD: A = U @ Σ @ V^H - Right null space: columns of V corresponding to zero singular values - Left null space: columns of U corresponding to zero singular values
Examples:
A = create_test_matrix(5, 3) # 5x3 matrix, rank ≤ 3 N_right = quat_null_space(A, side='right') # null(A) N_left = quat_null_space(A, side='left') # null(A^H) print(f"Right null space: {N_right.shape}") # (3, 3-rank) print(f"Left null space: {N_left.shape}") # (5, 5-rank)
quaternion_to_complex_adjoint(A, axis='x')
Map quaternion matrix A to 2n×2n complex adjoint matrix for a chosen complex subfield.
By default the complex subfield is tied to the x-axis (unit i): A = W + X i + Y j + Z k → C = W + i X, D = Y + i Z, Adj(A) = [[C, D], [-D, C]]
Parameters
A : np.ndarray (n×n, dtype=np.quaternion) axis : {"x"} Imaginary axis defining the complex subfield. Currently only "x" is supported.
Notes
Extending to other axes requires a consistent symplectic embedding relative to that axis.
rank(X, tol=None)
Compute the rank of a quaternion matrix by counting non-zero singular values.
Parameters:
X : numpy.ndarray with dtype=quaternion Quaternion matrix of any shape (m, n) tol : float, optional Tolerance for considering singular values as non-zero (default: machine epsilon * max(m,n) * max singular value)
Returns:
int : The rank of the matrix (number of non-zero singular values)
Notes:
The rank is computed by: 1. Computing the SVD of the matrix: X = U @ S @ V^H 2. Counting singular values above the tolerance threshold 3. The rank equals the number of non-zero singular values
The tolerance is automatically adjusted based on matrix size and magnitude to handle numerical precision issues.
real_contract(R, m, n)
Convert real block matrix R back to quaternion matrix. Invert real_expand back into an m×n quaternion array
real_expand(Q)
Convert quaternion matrix Q to real block matrix representation. Given an m×n quaternion array Q, return a (4m)×(4n) real block matrix [[Qw, -Qx, -Qy, -Qz], [Qx, Qw, -Qz, Qy], [Qy, Qz, Qw, -Qx], [Qz, -Qy, Qx, Qw]]
spectral_norm_2(A)
Matrix 2-norm (spectral norm): largest singular value of A.
Delegates to quaternion SVD implementation.
timesQsparse(B0, B1, B2, B3, C0, C1, C2, C3)
Quaternion matrix multiplication in component format.
Parameters:
B0, B1, B2, B3 : numpy.ndarray or scipy.sparse matrix First quaternion matrix components C0, C1, C2, C3 : numpy.ndarray or scipy.sparse matrix Second quaternion matrix components
Returns:
tuple : (A0, A1, A2, A3) - Result quaternion matrix components
Notes:
Implements quaternion multiplication: A = B * C Based on MATLAB implementation by Zhigang Jia