2
$\begingroup$

I'm looking for an algorithm to solve for the classic:

$$A\mathbf{x} = \mathbf{b}$$

I cannot compute $A$ directly, but rather can compute matrix-vector products $A\mathbf{v}$ for any $\mathbf{v}$.

At first, this appears to be a straightforward application of matrix-free solvers. Since I don't have any guarantees on symmetry or positive-definitness, this seems like a good application of the BICGSTAB algorithm.

The only problem is that I am interested in the special case where $\mathbf{b} = \mathbf{0}$, i.e. solving for a vector in the nullspace of $A$. Obviously, the trivial vector $\mathbf{x}$ = $\mathbf{0}$ is a solution, but not an interesting or useful one. Unfortunately, it's also one that matrix-free methods are prone to find. (It's also the min-norm solution, which I think is problematic as some solutions search for the min-norm solution.)

Is there any trick I can use to extract a non-trivial solution $\mathbf{x}$ with a matrix-free solver for general $A$?

(Cross-post from math.stackexchange.com; I think this problem is a little more research-y than I originally thought)

Improve this question
$\endgroup$
5
  • 4
    $\begingroup$ I wonder whether this would work: pick any nonzero vector $v$, calculate $Av$, let's say $Av=c$, then ask your solver to solve $Ax=c$. If it gives you some answer $w\ne v$, then $w-v$ is in the nullspace. $\endgroup$
    – Gerry Myerson
    Commented Dec 12, 2021 at 1:16
  • 2
    $\begingroup$ The m.se version is math.stackexchange.com/questions/4328461/… $\endgroup$
    – Gerry Myerson
    Commented Dec 12, 2021 at 11:38
  • $\begingroup$ I like this idea, it is simple, although it is entirely possible that we will receive an answer $w=v$. Is there some way we can modify this problem to guarantee they will not be equal to each other? Since we know $\mathbf{v} \cdot \mathbf{v} \neq 0$, if we concatenate 0 to the end of $\mathbf{c}$ (call this $\mathbf{d}$, and $\mathbf{v}$ as the final row of $A$ (call this $B$). we could force a vector $\mathbf{w} \perp \mathbf{v}$. Can we say for certain that we can solve $B\mathbf{x} = \mathbf{d}$ s.t. the vector from the first $n-1$ elements of solution $w$ are orthogonal to $v$? $\endgroup$
    – user650261
    Commented Dec 12, 2021 at 18:34
  • 1
    $\begingroup$ Have you considered the (admittedly, heavily unimaginative) approach to just compute the matrix $A$ by applying it to all canonical unit vectors (and then to let your favourite implementation of Gaussian elimination compute the rref of $[A | b]$)? $\endgroup$
    – Jochen Glueck
    Commented Dec 13, 2021 at 8:48
  • 2
    $\begingroup$ @JochenGlueck I have, it's a last-ditch fallback approach. It's not the most efficient solution, and will be hard to make scale without throwing a lot of hardware at the problem. $\endgroup$
    – user650261
    Commented Dec 13, 2021 at 18:06

3 Answers 3

Reset to default
2
$\begingroup$

Depending on the matrix free method: If it is iterative, you may just initialize it with a nonzero vector $x_0$. For example the Richardson iteration $x_{k+1} = x_k - A^T Ax_k$ does converge to the projection of $x_0$ onto the nullspace of $A$. Inpecting Krylov spaces one should see that something similaria true for Krylov methods like CG.

Improve this answer
$\endgroup$
1
  • $\begingroup$ This is a good idea we may try; I don't know enough about Richardson iteration's convergence rates to know if it's the best choice, but it seems simple enough to try. Other methods like BICGSTAB typically rely on the relative residual at some point, which is undefined when the RHS is all 0s and then causes problems. Some other CG-like methods may be suitable. $\endgroup$
    – user650261
    Commented Dec 13, 2021 at 16:51
1
$\begingroup$

This trick is often used for Markov chains: partition $$ A = \begin{bmatrix}A_{11} & a_{12} \\ a_{21} & \alpha_{22}\end{bmatrix}, $$ where the partition has sizes $(n-1,1)$, i.e., $\alpha_{22}$ is a scalar and $a_{12}$ and $a_{21}^T$ are single columns. Then, a vector $v = \begin{bmatrix}v_1 \\ 1 \end{bmatrix}$ (partitioned conformably) is in the kernel iff $A_{11}v_1 + a_{12} = 0$. This latter equation is a $(n-1) \times (n-1)$ linear system that you can solve with any matrix-free solver. It is easy to compute the action of $A_{11}$ from that of $A$ given in matrix-free form.

This method works if there exists a vector $v$ in the kernel with nonzero last entry, and if $A_{11}\in\mathbb{R}^{(n-1)\times(n-1)}$ is non-singular. In Markov chains, these properties often follow from connectedness arguments via the Perron-Frobenius theorem. Otherwise, one can permute rows and/or columns to use a different entry than the last.

In any case, the method will fail if $A$ has rank lower than $n-1$, since then $A_{11}$ is necessarily singular. If you know the rank $r$, you can do something similar by setting the last $r$ entries of $v$ to an arbitrary value and solving a system with $A_{11} \in \mathbb{R}^{(n-r)\times(n-r)}$.

Improve this answer
$\endgroup$
6
  • $\begingroup$ Is this method in some sense the same as @GerryMyerson's comment? $\endgroup$
    – LSpice
    Commented Dec 12, 2021 at 16:58
  • 1
    $\begingroup$ @LSpice I don't rule out that there could be a relationship, but at first sight they look quite different to me. One solves systems with $A$ and one with $A_{11}$, for instance. $\endgroup$
    – Federico Poloni
    Commented Dec 12, 2021 at 17:05
  • $\begingroup$ Maybe I am a little bit confused. You claim: "Then, a vector $v = \begin{bmatrix}v_1 \\ 1 \end{bmatrix}$ (partitioned conformably) is in the kernel iff $A_{11}v_1 + a_{12} = 0$" I with the forward part of this iff. But for the backward part of the claim, don't we also need that $a_{21}^T v_1 + a_{22} = 0$? What exactly is the proposed algorithm here to solve this problem? $\endgroup$
    – user650261
    Commented Dec 12, 2021 at 21:08
  • 1
    $\begingroup$ @user650261 I am assuming that $A$ is singular and $A_{11}$ is nonsingular. Then, the last row of $A$ must be a combination of the previous ones. (I can give a formal proof of this fact if you wish.) I have edited the post to better state an algorithm. $\endgroup$
    – Federico Poloni
    Commented Dec 13, 2021 at 7:14
  • 1
    $\begingroup$ I assume you meant "if $A_{11}$ is nonsingular", because that's the condition you really need, and the property on the last row is just a consequence of this. In general no, unfortunately, but it might be known theoretically like in the case of Markov chains. Or you could think of preconditioning $Ax=b$ with a random matrix to ensure it with high probability, see for instance "randomized Hadamard transforms" which is a similar technique from randomized linear algebra. (Just throwing around ideas here, nothing too definite unfortunately.) $\endgroup$
    – Federico Poloni
    Commented Dec 13, 2021 at 14:11
0
$\begingroup$

A colleague had an idea today that I think works. We'll test it soon but I'm curious if people here have feedback.

Construct the matrix $B$ where $B$ is $A$ with a row of $1$'s concatenated to the bottom, and vector $\mathbf{c}$ where $\mathbf{c}$ is $\mathbf{x}$ with a $1$ concatenated to the bottom. This "concatenation" can be done in a matrix-free way, since only the vector-vector product is needed.

Then $B\mathbf{x} = \mathbf{c}$ iff $A\mathbf{x} = \mathbf{b}$ and $\|\mathbf{x}\|_1 = 1$. The latter condition guarantees that the trivial vector will never be found, and it should always be possible to solve for an $\mathbf{x}$ that satisfies this because any nullspace vector can be arbitrarily rescaled and still remains in the nullspace. Thus, any solution to $B \mathbf{x} = \mathbf{c}$ is a nontrivial solution of $A\mathbf{x} = \mathbf{0}$.

Improve this answer
$\endgroup$
5
  • $\begingroup$ But the OP doesn't know $A$, this is the crux of the problem! $\endgroup$
    – Alex M.
    Commented Dec 13, 2021 at 17:07
  • $\begingroup$ As it’s currently written, your answer is unclear. Please edit to add additional details that will help others understand how this addresses the question asked. You can find more information on how to write good answers in the help center. $\endgroup$
    – Community Bot
    Commented Dec 13, 2021 at 17:08
  • 2
    $\begingroup$ Hi, I am OP :). We don't know $A$, but a vector-vector product can easily be concatenated to any matrix-vector product in a matrix-free way. $\endgroup$
    – user650261
    Commented Dec 13, 2021 at 17:08
  • $\begingroup$ (1) $B$ is not square, though. How do you plan to solve linear systems with it? (2) Actually your condition guarantees that $\sum x_i = 1$, which is a different thing from $\|x\|_1 = 1$. It could be possible that all vectors in the nullspace satisfy $\sum x_i = 0$, in which case the approach doesn't work. Imposing $\sum x_i = 1$ isn't too different from imposing $x_n = 1$ like I suggested; in both cases it's a linear condition (and a suitable change of bases can convert any linear condition into any other). $\endgroup$
    – Federico Poloni
    Commented Dec 13, 2021 at 19:57
  • $\begingroup$ (1) You can solve systems of non-square matrices, $B$ will be rank deficient but that doesn't mean it doesn't have a solution, and some matrix-free solvers will solve rectangular systems. In fact, in our setting, $A$ is rectangular. (2) That is a good point, that is an error on my part, but \|x\|_1 = 1 is at least one solution. If a solution satisfying a completely different constraint is found, that's fine - but it guarantees the trivial solution won't work. I think what's important here is that \|x\|_1 = 1 is a subset of the constraint, I'm unsure if a random vector works. $\endgroup$
    – user650261
    Commented Dec 13, 2021 at 21:12

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .