A = [1.0 0.0;0.0 1.0;0.0 0.0]

A1 = pinv(A):

[1.0 0.0 1.0;0.0 0.0 0.0] wrong

A1*A

[1.0 0.0;0.0 0.0] wrong

A2 = inv(A’*A)*A’:
[1.0 0.0 0.0;0.0 1.0 0.0] correct
A2*A

[1.0 0.0;0.0 1.0] correct

A = [1.0 0.0;0.0 1.0;0.0 0.0]

A1 = pinv(A):

[1.0 0.0 1.0;0.0 0.0 0.0] wrong

A1*A

[1.0 0.0;0.0 0.0] wrong

A2 = inv(A’*A)*A’:
[1.0 0.0 0.0;0.0 1.0 0.0] correct
A2*A

[1.0 0.0;0.0 1.0] correct

1 Like

Agreed, it looks like a bug:

```
julia> SVD = svd(A)
LinearAlgebra.SVD{Float64, Float64, Matrix{Float64}}
U factor:
3×2 Matrix{Float64}:
1.00000e+00 0.00000e+00
0.00000e+00 1.00000e+00
0.00000e+00 0.00000e+00
singular values:
2-element Vector{Float64}:
1.00000e+00
1.00000e+00
Vt factor:
2×2 Matrix{Float64}:
1.00000e+00 0.00000e+00
0.00000e+00 1.00000e+00
julia> (SVD.Vt * Diagonal(SVD.S)) * SVD.U' * A
2×2 Matrix{Float64}:
1.00000e+00 0.00000e+00
0.00000e+00 1.00000e+00
```

Can you spot it? https://github.com/JuliaLang/julia/blob/bf534986350a991e4a1b29126de0342ffd76205e/stdlib/LinearAlgebra/src/dense.jl#L1427

OK, a hint: is `A`

diagonal?

We should fix this bug, but it’s worth noting that you should pretty much never use `pinv`

. Doing linear solves directly is pretty much always better.

Right, but I am actually computing a projection matrix which will be used many times.

then call `factorize(A)`

and reuse that.

The internal part of `pinv`

is just the SVD, as the above posts show.

So you might want to look at using `svd(A)`

and using the SVD factors to construct your projection.

Thank you.