Repository: JuliaLang/julia
File: stdlib/LinearAlgebra/src/symmetric.jl (Julia)

Code before:
*(A::Transpose{<:Any,<:AbstractTriangular}, B::Adjoint{<:Any,<:RealHermSymComplexHerm}) = A * B.parent
*(A::Transpose{<:Any,<:AbstractTriangular}, B::Transpose{<:Any,<:RealHermSymComplexSym}) = A * B.parent
*(A::Adjoint{<:Any,<:RealHermSymComplexHerm}, B::Adjoint{<:Any,<:AbstractTriangular}) = A.parent * B
*(A::Adjoint{<:Any,<:RealHermSymComplexHerm}, B::Transpose{<:Any,<:AbstractTriangular}) = A.parent * B
*(A::Transpose{<:Any,<:RealHermSymComplexSym}, B::Adjoint{<:Any,<:AbstractTriangular}) = A.parent * B
*(A::Transpose{<:Any,<:RealHermSymComplexSym}, B::Transpose{<:Any,<:AbstractTriangular}) = A.parent * B

"""
    symmetrize!(X::AbstractMatrix; conjugate::Bool=true)

Make the matrix `X` symmetric in-place using the formula `(X + X') / 2`.
If `conjugate` is `true`, `X` is made Hermitian rather than symmetric.

See also: [`symmetrize`](@ref)

!!! compat "Julia 1.3"
    This function requires Julia 1.3 or later.
"""
symmetrize!(X::AbstractMatrix; conjugate::Bool=true) =
    conjugate ? _symmetrize!(conj, real, X) : _symmetrize!(identity, identity, X)

function _symmetrize!(f::Function, g::Function, X::AbstractMatrix)
    inds = axes(X, 1)
    inds == axes(X, 2) || throw(DimensionMismatch("matrix is not square"))
    r = first(inds)
    s = step(inds)
    @inbounds for j in inds
        X[j,j] = g(X[j,j])
        for i in r:s:j-s
            X[i,j] = (X[i,j] + f(X[j,i])) / 2
            X[j,i] = f(X[i,j])
        end
    end
    X
end

"""
    symmetrize(X::AbstractMatrix; conjugate::Bool=true)

Construct a symmetric matrix based on `X` using the formula `(X + X') / 2`.
If `conjugate` is `true`, the result is Hermitian rather than symmetric.

See also: [`symmetrize!`](@ref)

!!! compat "Julia 1.3"
    This function requires Julia 1.3 or later.
"""
symmetrize(X::AbstractMatrix{T}; conjugate::Bool=true) where {T<:Number} =
    symmetrize!(copyto!(similar(X, typeof(zero(T) / 2)), X), conjugate=conjugate)
symmetrize(X::Symmetric{<:Real}; conjugate::Bool=true) = X
symmetrize(X::Symmetric{<:Complex}; conjugate::Bool=true) = conjugate ? Hermitian(X) : X

Review comment (suggestion):
Is it possible to have `s != 1`?   `a[begin:end]` only works if the step is 1, so I thought we required consecutive ranges of indices.

Code after:
"""
hermitianpart(A::AbstractMatrix, uplo::Symbol=:U) = Hermitian(_hermitianpart(A), uplo)

"""
    hermitianpart!(A, uplo=:U) -> Hermitian

Overwrite the square matrix `A` in-place with its Hermitian part `(A + A') / 2`, and return
[`Hermitian(A, uplo)`](@ref). For real matrices `A`, this is also known as the symmetric
part of `A`.

See also [`hermitianpart`](@ref) for the corresponding out-of-place operation.

!!! compat "Julia 1.10"
    This function requires Julia 1.10 or later.
"""
hermitianpart!(A::AbstractMatrix, uplo::Symbol=:U) = Hermitian(_hermitianpart!(A), uplo)

_hermitianpart(A::AbstractMatrix) = _hermitianpart!(copy_similar(A, Base.promote_op(/, eltype(A), Int)))
_hermitianpart(a::Number) = real(a)

function _hermitianpart!(A::AbstractMatrix)
    require_one_based_indexing(A)
    n = checksquare(A)
    @inbounds for j in 1:n
        A[j, j] = _hermitianpart(A[j, j])
        for i in 1:j-1
            A[i, j] = val = (A[i, j] + adjoint(A[j, i])) / 2
            A[j, i] = adjoint(val)
        end
    end
    return A
end