and the rook-pivoted Bunch–Kaufman factorization is returning the wrong permutation.

Note, however, than most of the time you use the factorization object simply by x = F \ b to solve Ax=b, and the \ solver for rook-pivoted Bunch–Kaufman appears to be correct.

I noticed that \ appears to return the correct value (presumably because it calls LAPACK directly instead of relying on the permutation). However, in my use case I specifically need access to L and D.

In case you haven’t already derived the fix from the LAPACK docs, I think the following may do the job:

# adapted from LinearAlgebra stdlib
function _ipiv2p_rook(v::AbstractVector{T}, maxi::Integer, uplo::AbstractChar) where T
p = T[1:maxi;]
uploL = uplo == 'L'
i = uploL ? 1 : maxi
# if uplo == 'U' we construct the permutation backwards
@inbounds while 1 <= i <= length(v)
vi = v[i]
if vi > 0 # the 1x1 blocks
p[i], p[vi] = p[vi], p[i]
i += uploL ? 1 : -1
else # the 2x2 blocks
if uploL
p[i], p[-vi] = p[-vi], p[i]
vp = -v[i+1]
p[i + 1], p[vp] = p[vp], p[i + 1]
i += 2
else # 'U'
p[i], p[-vi] = p[-vi], p[i]
vp = -v[i-1]
p[i - 1], p[vp] = p[vp], p[i - 1]
i -= 2
end
end
end
return p
end
function permvec(B::BunchKaufman{T}) where {T}
B.rook || return B.p
n = size(B,1)
return _ipiv2p_rook(getfield(B, :ipiv), n, getfield(B, :uplo))
end

Then

F = bunchkaufman(A,rook)
p = permvec(F)
F.L * F.D * tfunc(F.L) ≈ A[p,p]

where tfunc is adjoint or transpose as appropriate.