Low-rank ADI for a continuous Lyapunov equation

Let $A\in F^{n\times n}$ be Hurwitz, let $B\in F^{n\times r}$, and let $p_1,\dots ,p_k\in \mathbb{C}$ be a finite shift list. The target equation is $AX+X{A}^{\ast }+B{B}^{\ast }=0$. Inputs also include the shift list and exact linear-solve operations. For real $A,B$, nonreal shifts must occur in consecutive complex-conjugate pairs; the procedure may use complex arithmetic for these paired shifts.

Before starting matrix updates, validate the real-data shift pairing by scanning the shift list once: at a real shift advance the scan index by one; at a nonreal shift require the next list entry to equal its complex conjugate, then advance the scan index by two. A missing or unequal partner returns unpaired-complex-shift. For complex data no pairing restriction is imposed. This validation scan consumes pairs; the update loop below still processes every individual shift in order.

Initialize $R_0=B$, let $Z_0$ have zero columns, and set $X_0=0$. For $j=1,\dots ,k$:

1. If $\mathrm{Re}{p}_j\ge 0$, return invalid-shift. If $A+p_{j}I$ is singular, return singular-shift.
2. Set $c_j=-2\mathrm{Re}{p}_j>0$ and solve $$V_j=(A+p_{j}I{)}^{-1}{R}_{j-1}\mathrm{.}$$
3. Set $$Z_j=[Z_{j-1},\sqrt{c_j}{V}_j],X_j=Z_j{Z}_j^{\ast }=X_{j-1}+c_j{V}_j{V}_j^{\ast },R_j=R_{j-1}+c_j{V}_j\mathrm{.}$$

Return completed with the low-rank factor $Z_k$, approximation $X_k$, and residual factor $R_k$. A finite shift list may stop before convergence; this procedure does not select shifts or assert unconditional convergence.