Aliased output matrix in LinearAlgebra.mul! does not throw error for matrix wrapper types

This came up here: https://www.reddit.com/r/Julia/comments/ugm4l7/working_with_linearalgebramul/

It looks like a bug to me, but I'm unable to test on 1.8. Perhaps someone can verify?

According to the docstring of LinearAlgebra.mul!:

Note that C must not be aliased with either A or B.

However:

julia> A = rand(5,5);

julia> B_ = rand(5,5);

julia> B = B_ + B_';

julia> C = similar(A); 

julia> mul!(C, A, Symmetric(B))  # works as expected
5×5 Matrix{Float64}:
 1.10149  1.93057  1.67467  1.39908  1.59191
 1.37769  1.07008  2.50417  1.21931  2.02994
 1.82134  2.73737  3.66532  1.5658   3.55056
 1.45685  1.44589  2.73083  1.17807  2.39641
 1.79866  2.58796  2.48849  1.33673  2.46281

julia> mul!(A, A, B)  # fails as expected
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 C:\Users\dnf\.julia\juliaup\julia-1.7.2+0~x64\share\julia\stdlib\v1.7\LinearAlgebra\src\matmul.jl:647
 [2] mul!
   @ C:\Users\dnf\.julia\juliaup\julia-1.7.2+0~x64\share\julia\stdlib\v1.7\LinearAlgebra\src\matmul.jl:169 [inlined]
 [3] mul!(C::Matrix{Float64}, A::Matrix{Float64}, B::Matrix{Float64})
   @ LinearAlgebra C:\Users\dnf\.julia\juliaup\julia-1.7.2+0~x64\share\julia\stdlib\v1.7\LinearAlgebra\src\matmul.jl:275
 [4] top-level scope
   @ REPL[91]:1

julia> mul!(A, A, Symmetric(B))  # does not fail as it should, and produces wrong result
5×5 Matrix{Float64}:
 0.201062  0.334617  0.478842  0.102553  1.21272
 0.32496   0.540814  0.773914  0.165748  1.96002
 0.575235  0.957335  1.36996   0.293402  3.46957
 0.362894  0.603946  0.864256  0.185096  2.18882
 0.439629  0.731652  1.04701   0.224235  2.65165

The same thing happens for other matrix wrapper types, like Diagonal etc.

julia> versioninfo()
Julia Version 1.7.2
Commit bf53498635 (2022-02-06 15:21 UTC)
Platform Info:
  OS: Windows (x86_64-w64-mingw32)
  CPU: Intel(R) Core(TM) i7-9750H CPU @ 2.60GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-12.0.1 (ORCJIT, skylake)
Environment:
  JULIA_EDITOR = code

Start by reproducing the examples with Symmetric and Diagonal wrappers, then read LinearAlgebra/src/matmul.jl around gemm_wrapper! at line 647 and the mul! methods. Done means aliased output is rejected with the expected ArgumentError for wrapper types instead of producing a wrong result.