Matrix decomposition for a linear combination of matrices

Question

I want to solve the linear system, $$(\mathbf{A}+\alpha \mathbf{B}) \mathbf{x} = \mathbf{b},$$ where $\mathbf{A,B}\in\mathbb{R}^{n\times n}$, $\mathbf{x,b}\in\mathbb{R}^{n\times 1}$ and $\alpha\in\mathbb{R}$.

Let say that I want to solve that linear system for different $\alpha$ values. Is there any $\mathbf{A}+\alpha \mathbf{B}$ decomposition that could be easily updated for different $\alpha$ values so I could solve fast each system? If not, is there any way to solve this set of systems fast?

PD: I do not need a full answer but just some hints in case it can be solved fast.

Answer 1

Suppose that

$$(A+\alpha B)x = b$$ $A,B \in \mathbb{R}^{n \times n} , x ,b \in \mathbb{R}^{n} , \alpha \in \mathbb{R} $

I'm not sure how to preserve the original matrices if you just want the solution

Let $C= A+\alpha B$

now $C =U \Sigma V^{T} $

$$ Cx=b \implies U\Sigma V^{T}x =b$$ $$ C^{+} = V \Sigma^{-1} U^{T} $$ $$ x = V \Sigma^{-1}U^{T} b $$

You should understand that even slight changes in $\alpha$ could change the condition number $(A+\alpha B)$ for some matrices which would destroy your precision. The computational time for solving this is very little to update it.

Answer 2

This is a very partial answer.

First, I am working over an arbitrary algebraically closed field $\mathbb{K}$. Let $\mathbf{A}$ and $\mathbf{B}$ be $n$-by-$n$ matrices over $\mathbb{K}$. The only situation in which you can invert $$\mathbf{X}(\alpha):=\mathbf{A}+\alpha\,\mathbf{B}$$ for any $\alpha \in\mathbb{K}$ is when $\mathbf{A}$ is invertible and $\mathbf{A}^{-1}\,\mathbf{B}$ is nilpotent.

Let $\mathbf{b}\in\mathbb{K}^n$ be given. Thus, the only $\mathbf{x}(\alpha)\in\mathbb{K}^n$ satisfying $$\mathbf{X}(\alpha)\,\mathbf{x}(\alpha)=\mathbf{b}$$ must satisfy $$\left(\mathbf{I}+\alpha \,\mathbf{A}^{-1}\,\mathbf{B}\right)\,\mathbf{x}(\alpha)=\mathbf{A}^{-1}\,\mathbf{b}\,,$$ where $\mathbf{I}$ is the $n$-by-$n$ identity matrix. Because $\mathbf{A}^{-1}\,\mathbf{B}$ is a nilpotent $n$-by-$n$ matrix, we conclude that $\left(\mathbf{A}^{-1}\,\mathbf{B}\right)^n$ is the zero matrix. Ergo, $$\mathbf{x}(\alpha)=\big(\mathbf{I}+\alpha\,\mathbf{A}^{-1}\,\mathbf{B}\big)\,\mathbf{A}^{-1}\,\mathbf{b}=\sum_{k=0}^{n-1}\,(-1)^k\,\alpha^k\,\left(\mathbf{A}^{-1}\,\mathbf{B}\right)^{k}\,\mathbf{A}^{-1}\,\mathbf{b}\,,$$ where we use the convention $t^0:=1$ and $\mathbf{M}^0:=\mathbf{I}$ for any $t\in\mathbb{K}$ and any $n$-by-$n$ matrix $\mathbf{M}$ over $\mathbb{K}$. Thus, the only lengthy computation is the computation of $\mathbf{A}^{-1}$ (and $\mathbf{A}^{-1}\,\mathbf{B}$).

---

If you work over $\mathbb{R}$ or $\mathbb{C}$, then you can fix an operator norm $\|\_\|$. The easiest choice is perhaps the Frobenius norm. When $\mathbf{A}$ is invertible and $$|\alpha|\,\left\|\mathbf{A}^{-1}\,\mathbf{B}\right\|<1\,,$$ then $\mathbf{X}(\alpha)$ is invertible and you can approximate $\mathbf{x}(\alpha)$ using the infinite sum below: $$\mathbf{x}(\alpha)=\sum_{k=0}^\infty\,(-1)^k\,\alpha^k\,\left(\mathbf{A}^{-1}\,\mathbf{B}\right)^k\,\mathbf{A}^{-1}\,\mathbf{b}\,.\tag{*}$$ The convergence will be painfully slow when $\delta:=|\alpha |\,\left\|\mathbf{A}^{-1}\,\mathbf{B}\right\|$ is very close to $1$. Nevertheless, if you want a certain amount of accuracy (say, within $\epsilon$), you can pre-compute $\left(\mathbf{A}^{-1}\,\mathbf{B}\right)^k$ up to $$k\approx \dfrac{\ln(1-\delta)+\ln(\epsilon)}{\ln(\delta)}$$ and be done with the job for all $\alpha$ with $|\alpha|\leq \dfrac{\delta}{\left\|\mathbf{A}^{-1}\,\mathbf{B}\right\|}$.

---

The field $\mathbb{K}$ is either $\mathbb{R}$ or $\mathbb{C}$. For some special $\alpha$, we can use a similar power-series approach as in (*).

If $\mathbf{A}$ is singular but $\mathbf{B}$ is invertible, then you can just look at $$\mathbf{Y}(\beta):=\mathbf{B}+\beta\,\mathbf{A}\,.$$ For $\alpha\neq 0$, we see that $$\alpha\,\mathbf{Y}\left(\frac{1}\alpha\right)=\mathbf{X}(\alpha)$$ and we can do a similar trick to (*) with $\mathbf{Y}(\beta)$ instead.

In the case where both $\mathbf{A}$ and $\mathbf{B}$ are singular, you simply have to find an element $\lambda\in\mathbb{K}$ such that $\mathbf{C}:=\mathbf{A}+\lambda\,\mathbf{B}$ is nonsingular. If $\lambda$ does not exist, then you will be solving degenerate linear system for all $\lambda$, which is obviously no fun. Nonetheless, if you are lucky, this $\lambda$ exists, and we can use $$\mathbf{Z}(\gamma):=\mathbf{C}+\gamma\,\mathbf{B}$$ instead, noting that $\mathbf{X}(\alpha)=\mathbf{Z}(\alpha-\lambda)$. Then, you can work with $\mathbf{Z}(\gamma)$ just like before.