LU decomposition in Julia

I think the confusion here is that, even though this is selecting the k-th row of U, the result of a 1d slice in Julia is always treated as a “column vector” (i.e. a 1d slice is a 1d array, no matter what dimension is sliced).

So, I think what you want is:

U[k+1:n,:] = U[k+1:n,:] - L[k+1:n,k]*transpose(U[k,:])

to make sure that U[k,:] is treated as a row vector for your rank-1 update here.

See also Unintuitive Julia result: selecting a row from a matrix and Problem: extracting a row from an array, returns a column and Trouble Understanding Slicing and RFC: Drop dimensions indexed by scalars by mbauman · Pull Request #13612 · JuliaLang/julia · GitHub (which first implemented this) and Taking vector transposes seriously · Issue #4774 · JuliaLang/julia · GitHub (where it was proposed as “APL-style indexing”).

1 Like