LU Decomposition - return permutation matrix instead of vector

New to Julia
,
1

Dear Community,

Is there any way to return in Julia full permutation matrix instead of vector?

using LinearAlgebra
A = [0 4 2; 10 2 1; 1 1 1]
L,U,P = lu(A)

P is returned as [2,1,3].
However, I would like to have in the following form: [0,1,0;1,0,0;0,0,1]

Is there any way to enforce this in Julia?

Thank you in advance for any help

2 

julia> p = [2, 1, 3];

julia> n = length(p);

julia> II = Matrix{Int}(I, n, n)
3×3 Array{Int64,2}:
 1  0  0
 0  1  0
 0  0  1

julia> II[p, :]
3×3 Array{Int64,2}:
 0  1  0
 1  0  0
 0  0  1
5 Likes
3

Alternatively you can also use I(n) to get the identity matrix.

julia> II = Matrix{Int}(I, 3, 3)
3×3 Array{Int64,2}:
 1  0  0
 0  1  0
 0  0  1

julia> I(3)
3×3 Diagonal{Bool,Array{Bool,1}}:
 1  ⋅  ⋅
 ⋅  1  ⋅
 ⋅  ⋅  1
2 Likes
4

Except then indexing into it with the vector gives you a sparse matrix which may or may not be what you want / need.

julia> I(3)[p, :]
3×3 SparseArrays.SparseMatrixCSC{Bool,Int64} with 3 stored entries:
  [2, 1]  =  1
  [1, 2]  =  1
  [3, 3]  =  1
5

Yeah, I see the issue here…

6

Or simply

julia> using LinearAlgebra

julia> A = rand(5,5);

julia> lu(A).P
5×5 Array{Float64,2}:
 0.0  1.0  0.0  0.0  0.0
 0.0  0.0  1.0  0.0  0.0
 0.0  0.0  0.0  1.0  0.0
 0.0  0.0  0.0  0.0  1.0
 1.0  0.0  0.0  0.0  0.0
3 Likes
7

Why is that a Float64 matrix instead of an Int matrix?! The permutation vector is an Int vector.

8

I wonder if you read the docs of lu, since you can get it directly as a property:

julia> lu(A).P
3×3 Array{Float64,2}:
 0.0  1.0  0.0
 1.0  0.0  0.0
 0.0  0.0  1.0
1 Like