# Nullspace Function Vs QR decomposition

**URL:** https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656
**Category:** Numerics
**Created:** [February 23, 2024, 5:01pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656 "2024-02-23T17:01:36Z")
**Posts on this page:** 9
**Page:** 1

<div class="post-metadata">

### Author: ![fimiller](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fimiller/32/51226_2.png) [@fimiller](https://discourse.julialang.org/u/fimiller)
#### Post date: [February 23, 2024, 5:01pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/1 "2024-02-23T17:01:37Z")

</div>

I am attempting to compute a nullspace for a matrix A. I have two potential techniques of doing this.

One is to use a permuted QR decomposition, and the other is to use the nullspace function.

These two approaches give different nullspace basis that drastically impact some experiments further down the line.

I was wondering if someone could provide the documentation needed for me to understand what the LinearAlgebra.nullspace function algorithm is, so I can compare against my algorithm with a permuted QR to see why they differ so greatly.

Thank you!

---

<div class="post-metadata">

### Author: ![Oscar\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oscar_smith/32/25343_2.png) [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)
#### Post date: [February 23, 2024, 5:24pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/2 "2024-02-23T17:24:55Z")

</div>

The basic problem is that there is no fully sensible way of computing the rank of numerical matrices. `nullspace` uses the SVD which is slightly more robust, but the real problem here is that fundamentally there is no numerically stable way to do this reliably.

---

<div class="post-metadata">

### Author: ![fimiller](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fimiller/32/51226_2.png) [@fimiller](https://discourse.julialang.org/u/fimiller)
#### Post date: [February 23, 2024, 5:52pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/3 "2024-02-23T17:52:24Z")

</div>

When you say “no sensible way of computing rank”, I was under the impression that numerical rank is looking for singular values above a certain tolerance. Also, there are numerically stable ways (to my understanding) of doing these matrix factorizations (modified gram schmidt for QR, for example).

One extra thought is that, as my matrix gets large, storing it in memory becomes problematic. Are there any matrix free / iterative procedures currently out there for computing a nullspace? Some looking online did not have a lot of results.

---

<div class="post-metadata">

### Author: ![stabbles](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stabbles/32/946_2.png) [@stabbles](https://discourse.julialang.org/u/stabbles)
#### Post date: [February 23, 2024, 9:09pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/4 "2024-02-23T21:09:33Z")

</div>

> One extra thought is that, as my matrix gets large, storing it in memory becomes problematic. Are there any matrix free / iterative procedures currently out there for computing a nullspace? Some looking online did not have a lot of results.

You can solve Ax = 0 with conjugate gradients actually, provided A is symmetric and does have a zero eigenvalue. I don’t know how you would find a full basis though without knowing its dimension ahead of time.

---

<div class="post-metadata">

### Author: ![fimiller](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fimiller/32/51226_2.png) [@fimiller](https://discourse.julialang.org/u/fimiller)
#### Post date: [February 23, 2024, 9:16pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/5 "2024-02-23T21:16:41Z")

</div>

This will happen. Do you know if ts possible if I have a bound on the rank? Could I compute at least n nullspace vectors with conjugate gradient?

---

<div class="post-metadata">

### Author: ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [February 23, 2024, 9:40pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/6 "2024-02-23T21:40:22Z")

</div>

> [@fimiller](#):
>
> Could I compute at least n nullspace vectors with conjugate gradient?

In theory yes, by starting with different random vectors. In practice, my recollection is that CG can sometimes experience breakdown for semidefinite matrices, due to floatingh-point issues, so care is required. You might be better off with iterative eigensolver methods like LOBPCG or Lanczos.

> [@fimiller](#):
>
> These two approaches give different nullspace basis that drastically impact some experiments further down the line.

When you say “different nullspace bases”, what do you mean? The basis is not unique, of course, since you can multiply it by any unitary matrix.

As @Oscar_Smith alludes to, because numerical computations of nullspace bases and ranks involve an arbitrary threshold on the smallest singular value (or the smallest R diagonal for QR), so even running the same algorithm on different machines (much less entirely different algorithms like QR vs SVD) will sometimes give different dimensions/ranks. You’ll need to think carefully about your problem if you don’t want your algorithm to be sensitive to this issue.

There isn’t any builtin function to compute the nullspace from a QR factorization (though perhaps there should be?), so without seeing your code it’s not clear whether you might have a bug there too.

See also the discussions in [How to find the linearly independent columns (rows) of a matrix - #14 by stevengj](https://discourse.julialang.org/t/how-to-find-the-linearly-independent-columns-rows-of-a-matrix/109772/14) and [How to find the linearly independent columns (rows) of a matrix - #25 by mstewart](https://discourse.julialang.org/t/how-to-find-the-linearly-independent-columns-rows-of-a-matrix/109772/25) on rank from QR.

---

<div class="post-metadata">

### Author: ![stabbles](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stabbles/32/946_2.png) [@stabbles](https://discourse.julialang.org/u/stabbles)
#### Post date: [February 23, 2024, 10:04pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/7 "2024-02-23T22:04:22Z")

</div>

Block methods for the eigenproblem are probably fastest, but only if the matrix isn’t (too) indefinite. And you still have to (over)estimate the dimension of the null space when picking the block size.

You could also try `partialschur` from ArnoldiMethod.jl with DoubleFloats.jl numbers if you need more accuracy at reasonable speed 🙂 with the caveat that it is not a block method, so may not find a full basis for the null space.

---

<div class="post-metadata">

### Author: ![fimiller](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fimiller/32/51226_2.png) [@fimiller](https://discourse.julialang.org/u/fimiller)
#### Post date: [February 24, 2024, 7:15pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/8 "2024-02-24T19:15:36Z")

</div>

I see. Thanks for the long response! I will definitely use some of these ideas for solving my problem.

---

<div class="post-metadata">

### Author: ![fimiller](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fimiller/32/51226_2.png) [@fimiller](https://discourse.julialang.org/u/fimiller)
#### Post date: [February 24, 2024, 7:16pm UTC](https://discourse.julialang.org/t/nullspace-function-vs-qr-decomposition/110656/9 "2024-02-24T19:16:54Z")

</div>

I will also look into this. Thanks!
