A Family of Iteration Functions for General Linear Systems

We introduce innovative algorithms for computing exact or approximate (minimum-norm) solutions to $Ax=b$ or the {\it normal equation} $A^TAx=ATb$, where $A$ is an $m \times n$ real matrix of arbitrary rank. We present more efficient algorithms when $A$ is symmetric PSD. First, we introduce the {\it Triangle Algorithm} (TA), a {\it convex-hull membership} algorithm that given $bk=Axk$ in the ellipsoid $E{A,\rho}={Ax: \Vert x \Vert \leq \rho}$, it either computes an improved approximation $b{k+1}=Ax{k+1}$ or proves $b \not \in E{A,\rho}$. We then give a dynamic variant of TA, the {\it Centering Triangle Algorithm} (CTA), generating residual, $rk=b -Axk$ via the iteration of $F1(r)=r-(r^{THr/rTH2r)Hr$,} where $H=AA^T$. If $A$ is symmetric PSD, $H$ can be taken as $A$. Next, for each $t=1, \dots, m$, we derive $Ft(r)=r- \sum{i=1}^t \alpha{t,i}(r) Hⁱ r$ whose iterations correspond to a Krylov subspace method with restart. If $\kappa^+(H)$ is the ratio of the largest to smallest positive eigenvalues of $H$, when $Ax=b$ is consistent, in $k=O({\kappa^+(H)}{t{-1}} \ln \varepsilon^{-1})$ iterations of $Ft$, $\Vert rk \Vert \leq \varepsilon$. Each iteration takes $O(tN+t^3)$ operations, $N$ the number of nonzero entries in $A$. By directly applying $Ft$ to the normal equation, we get $\Vert A^TAxk - A^Tb \Vert \leq \varepsilon$ in $O({\kappa^{+(AAT)}{t}{-1}} \ln \varepsilon^{-1})$ iterations. On the other hand, given any residual $r$, we compute $s$, the degree of its minimal polynomial with respect to $H$ in $O(sN+s^3)$ operations. Then $F_s(r)$ gives the minimum-norm solution of $Ax=b$ or an exact solution of $A^TAx=ATb$. The proposed algorithms are simple to implementation and theoretically robust. We present sample computational results, comparing the performance of CTA with CG and GMRES methods. The results support CTA as a highly competitive option.