Numerical stability of the symplectic $LL^T$ factorization

Abstract

In this paper we give the detailed error analysis of two algorithms $W_1$ and $W_2$ for computing the symplectic factorization of a symmetric positive definite and symplectic matrix $A \in \mathbb R^{2n \times 2n}$ in the form $A=LL^T$, where $L \in \mathbb R^{2n \times 2n}$ is a symplectic block lower triangular matrix. We prove that Algorithm $W_2$ is numerically stable for a broader class of symmetric positive definite matrices $A \in \mathbb R^{2n \times 2n}$. It means that Algorithm $W_2$ is producing the computed factors $\tilde L$ in floating-point arithmetic with machine precision $\mathcal{u}$ such that $||A-\tilde L {\tilde L}^T||_{2} = {\cal O}(\mathcal{u} ||{A}||_{2})$. On the other hand, Algorithm $W_1$ is unstable, in general, for symmetric positive definite and symplectic matrix $A$. In this paper we also give corresponding bounds for Algorithm $W_1$ that are weaker. We show that the factorization error depends on the condition number $\kappa_2(A_{11})$ of the principal submatrix $A_{11}$. Bounds for the loss of symplecticity of the lower block triangular matrices $L$ for both Algorithms $W_1$ and $W_2$ that hold in exact arithmetic for a broader class of symmetric positive definite matrices $A$ (but not necessarily symplectic) are also given. The tests performed in \textsl{MATLAB} illustrate that our error bounds for considered algorithms are reasonably sharp.