Reposting a reddit post here to solicit comments from the Julia community. I may have chanced upon a bug? According to reddit users `No-Distribution4263`

and `Foura5`

both `foo!`

and `bar!`

should result in an error. But the case below, `foo!`

works while `bar!`

does.

Also, n my use case `m`

and `n`

will be large thus `Z`

will be a large array. What is the best way of mutating the values in `Z`

without allocating another array?

Thanks very much.

```
using LinearAlgebra
using Random
using Statistics
Sigma = Symmetric(
[1.00 0.25 0.50;
0.25 1.00 0.15;
0.50 0.15 1.00]
)
mu = [1.50; 2.00; 2.50]
seed = 42
m = 10
n = 3
# foo!
C = cholesky(Sigma)
Z = randn(Xoshiro(seed), (m, n))
function foo!(F::Cholesky,
Z::Matrix{Float64},
mu::Vector{Float64})
Z .= mul!(Z, Z, F.U) .+ mu'
return Z
end
julia> foo!(C, Z, mu)
10×3 Matrix{Float64}:
1.13664 1.11003 1.09999
1.75174 2.87654 1.38161
1.18501 2.34136 3.47988
1.18875 1.5393 4.29595
2.31631 2.70479 2.74265
1.97674 3.52045 2.50674
0.640445 1.69499 2.03278
0.0307118 0.400502 1.82319
-0.614335 1.86276 0.617501
1.54378 1.91743 3.35531
# bar!
E = eigen(Sigma)
Z = randn(Xoshiro(seed), (m, n))
function bar!(F::Eigen,
Z::Matrix{Float64},
mu::Vector{Float64})
mul!(F.vectors, Diagonal(Base.sqrt.(F.values)), F.vectors')
Z .= mul!(Z, Z, F.vectors) .+ mu'
return Z
end
julia> bar!(E, Z, mu)
ERROR: ArgumentError: output matrix must not be aliased with input matrix
Stacktrace:
[1] gemm_wrapper!(C::Matrix{Float64}, tA::Char, tB::Char, A::Matrix{Float64}, B::Matrix{Float64}, _add::LinearAlgebra.MulAddMul{true, true, Bool, Bool})
@ LinearAlgebra ~/.julia/share/julia/stdlib/v1.7/LinearAlgebra/src/matmul.jl:647
[2] mul!
@ ~/.julia/share/julia/stdlib/v1.7/LinearAlgebra/src/matmul.jl:169 [inlined]
[3] mul!
@ ~/.julia/share/julia/stdlib/v1.7/LinearAlgebra/src/matmul.jl:275 [inlined]
[4] bar!(F::Eigen{Float64, Float64, Matrix{Float64}, Vector{Float64}}, Z::Matrix{Float64}, mu::Vector{Float64})
@ Main ./REPL[15]:5
[5] top-level scope
@ REPL[16]:1
```