I want to quickly solve the following linear least-squares problem


$$\min\_{x \in \mathbb{R}^n} \left\| A x - b \right\|\_2^2$$


with a special sparse structure where each row in $A$ has only up to $4$ consecutive non-zero entries. This makes its normal matrix


$$C = A^T A$$


a positive definite $7$-*band* matrix with a condition number between $8^2$ and $400^2$. So, these condition numbers aren't too bad that solving the Gauß normal equation system


$$ C x = A^T b $$


instead wouldn't get me into trouble numerically, I think. What are my options? I could try conjugate gradient methods but I would prefer direct solvers that can deal with these kinds of special cases in $O(n)$ time independent of the condition. I'm aware of algorithms for the tridiagonal case and I guess I could try to adapt them for $5$ bands (?). But before reinventing the wheel and/or testing many different algorithms I wanted to ask you about what approach might be the most efficient in terms of time because I have lots of these problems (millions) with values for $n$ of around $5000$ where $A$ has about $4n$ rows.



Solving sparse linear least squares or a positive definite 5-band matrix system fast
The LU factorization of $C$ along with forward and backward substitution works well in this case. The factorization can still be done completely in-place. So, there is no need to touch or create other off-band elements.


