Merge pull request #25364 from Sacha0/lowercase

adjoint/transpose transcend material concerns

Author: Sacha0

base/boot.jl

@@ -212,8 +212,6 @@ macro _noinline_meta()
     Expr(:meta, :noinline)
 end
 
-function postfixapostrophize end
-
 struct BoundsError <: Exception
     a::Any
     i::Any

base/deprecated.jl

@@ -51,22 +51,40 @@ Transpose(x::Number) = transpose(x)
 # unwrapping constructors
 Adjoint(A::Adjoint) = A.parent
 Transpose(A::Transpose) = A.parent
-# normalizing unwrapping constructors
-# technically suspect, but at least fine for now
-Adjoint(A::Transpose) = conj(A.parent)
-Transpose(A::Adjoint) = conj(A.parent)
 
-# eager lowercase quasi-constructors, unwrapping
-adjoint(A::Adjoint) = copy(A.parent)
-transpose(A::Transpose) = copy(A.parent)
-# eager lowercase quasi-constructors, normalizing
-# technically suspect, but at least fine for now
-adjoint(A::Transpose) = conj!(copy(A.parent))
-transpose(A::Adjoint) = conj!(copy(A.parent))
-
-# lowercase quasi-constructors for vectors, TODO: deprecate
-adjoint(sv::AbstractVector) = Adjoint(sv)
-transpose(sv::AbstractVector) = Transpose(sv)
+# wrapping lowercase quasi-constructors
+adjoint(A::AbstractVecOrMat) = Adjoint(A)
+"""
+    transpose(A::AbstractMatrix)
+
+Lazy matrix transpose. Mutating the returned object should appropriately mutate `A`. Often,
+but not always, yields `Transpose(A)`, where `Transpose` is a lazy transpose wrapper. Note
+that this operation is recursive.
+
+This operation is intended for linear algebra usage - for general data manipulation see
+[`permutedims`](@ref), which is non-recursive.
+
+# Examples
+```jldoctest
+julia> A = [1 2 3; 4 5 6; 7 8 9]
+3×3 Array{Int64,2}:
+ 1  2  3
+ 4  5  6
+ 7  8  9
+
+julia> transpose(A)
+3×3 Transpose{Int64,Array{Int64,2}}:
+ 1  4  7
+ 2  5  8
+ 3  6  9
+```
+"""
+transpose(A::AbstractVecOrMat) = Transpose(A)
+# unwrapping lowercase quasi-constructors
+adjoint(A::Adjoint) = A.parent
+transpose(A::Transpose) = A.parent
+adjoint(A::Transpose{<:Real}) = A.parent
+transpose(A::Adjoint{<:Real}) = A.parent
 
 
 # some aliases for internal convenience use
@@ -187,13 +205,14 @@ pinv(v::TransposeAbsVec, tol::Real = 0) = pinv(conj(v.parent)).parent
 ## right-division \
 /(u::AdjointAbsVec, A::AbstractMatrix) = Adjoint(Adjoint(A) \ u.parent)
 /(u::TransposeAbsVec, A::AbstractMatrix) = Transpose(Transpose(A) \ u.parent)
-
+/(u::AdjointAbsVec, A::Transpose{<:Any,<:AbstractMatrix}) = Adjoint(conj(A.parent) \ u.parent) # technically should be Adjoint(copy(Adjoint(copy(A))) \ u.parent)
+/(u::TransposeAbsVec, A::Adjoint{<:Any,<:AbstractMatrix}) = Transpose(conj(A.parent) \ u.parent) # technically should be Transpose(copy(Transpose(copy(A))) \ u.parent)
 
 # dismabiguation methods
-*(A::AdjointAbsVec, B::Transpose{<:Any,<:AbstractMatrix}) = A * transpose(B.parent)
-*(A::TransposeAbsVec, B::Adjoint{<:Any,<:AbstractMatrix}) = A * adjoint(B.parent)
-*(A::Transpose{<:Any,<:AbstractMatrix}, B::Adjoint{<:Any,<:AbstractMatrix}) = transpose(A.parent) * B
-*(A::Adjoint{<:Any,<:AbstractMatrix}, B::Transpose{<:Any,<:AbstractMatrix}) = A * transpose(B.parent)
+*(A::AdjointAbsVec, B::Transpose{<:Any,<:AbstractMatrix}) = A * copy(B)
+*(A::TransposeAbsVec, B::Adjoint{<:Any,<:AbstractMatrix}) = A * copy(B)
+*(A::Transpose{<:Any,<:AbstractMatrix}, B::Adjoint{<:Any,<:AbstractMatrix}) = copy(A) * B
+*(A::Adjoint{<:Any,<:AbstractMatrix}, B::Transpose{<:Any,<:AbstractMatrix}) = A * copy(B)
 # Adj/Trans-vector * Trans/Adj-vector, shouldn't exist, here for ambiguity resolution? TODO: test removal
 *(A::Adjoint{<:Any,<:AbstractVector}, B::Transpose{<:Any,<:AbstractVector}) = throw(MethodError(*, (A, B)))
 *(A::Transpose{<:Any,<:AbstractVector}, B::Adjoint{<:Any,<:AbstractVector}) = throw(MethodError(*, (A, B)))

base/linalg/adjtrans.jl

@@ -247,8 +247,14 @@ broadcast(::typeof(trunc), ::Type{T}, M::Bidiagonal) where {T<:Integer} = Bidiag
 broadcast(::typeof(floor), ::Type{T}, M::Bidiagonal) where {T<:Integer} = Bidiagonal(floor.(T, M.dv), floor.(T, M.ev), M.uplo)
 broadcast(::typeof(ceil), ::Type{T}, M::Bidiagonal) where {T<:Integer} = Bidiagonal(ceil.(T, M.dv), ceil.(T, M.ev), M.uplo)
 
-transpose(M::Bidiagonal) = Bidiagonal(M.dv, M.ev, M.uplo == 'U' ? :L : :U)
-adjoint(M::Bidiagonal) = Bidiagonal(conj(M.dv), conj(M.ev), M.uplo == 'U' ? :L : :U)
+adjoint(B::Bidiagonal) = Adjoint(B)
+transpose(B::Bidiagonal) = Transpose(B)
+adjoint(B::Bidiagonal{<:Real}) = Bidiagonal(B.dv, B.ev, B.uplo == 'U' ? :L : :U)
+transpose(B::Bidiagonal{<:Number}) = Bidiagonal(B.dv, B.ev, B.uplo == 'U' ? :L : :U)
+Base.copy(aB::Adjoint{<:Any,<:Bidiagonal}) =
+    (B = aB.parent; Bidiagonal(map(x -> copy.(Adjoint.(x)), (B.dv, B.ev))..., B.uplo == 'U' ? :L : :U))
+Base.copy(tB::Transpose{<:Any,<:Bidiagonal}) =
+    (B = tB.parent; Bidiagonal(map(x -> copy.(Transpose.(x)), (B.dv, B.ev))..., B.uplo == 'U' ? :L : :U))
 
 istriu(M::Bidiagonal) = M.uplo == 'U' || iszero(M.ev)
 istril(M::Bidiagonal) = M.uplo == 'L' || iszero(M.ev)
@@ -494,15 +500,15 @@ const SpecialMatrix = Union{Bidiagonal,SymTridiagonal,Tridiagonal}
 
 #Generic multiplication
 *(A::Bidiagonal{T}, B::AbstractVector{T}) where {T} = *(Array(A), B)
-*(adjA::Adjoint{<:Any,<:Bidiagonal{T}}, B::AbstractVector{T}) where {T} = *(Adjoint(Array(adjA.parent)), B)
-*(A::Bidiagonal{T}, adjB::Adjoint{<:Any,<:AbstractVector{T}}) where {T} = *(Array(A), Adjoint(adjB.parent))
+*(adjA::Adjoint{<:Any,<:Bidiagonal{T}}, B::AbstractVector{T}) where {T} = *(adjoint(Array(adjA.parent)), B)
+*(A::Bidiagonal{T}, adjB::Adjoint{<:Any,<:AbstractVector{T}}) where {T} = *(Array(A), adjoint(adjB.parent))
 /(A::Bidiagonal{T}, B::AbstractVector{T}) where {T} = /(Array(A), B)
-/(A::Bidiagonal{T}, adjB::Adjoint{<:Any,<:AbstractVector{T}}) where {T} = /(Array(A), Adjoint(adjB.parent))
+/(A::Bidiagonal{T}, adjB::Adjoint{<:Any,<:AbstractVector{T}}) where {T} = /(Array(A), adjoint(adjB.parent))
 
 #Linear solvers
 ldiv!(A::Union{Bidiagonal, AbstractTriangular}, b::AbstractVector) = naivesub!(A, b)
-ldiv!(transA::Transpose{<:Any,<:Bidiagonal}, b::AbstractVector) = ldiv!(transpose(transA.parent), b)
-ldiv!(adjA::Adjoint{<:Any,<:Bidiagonal}, b::AbstractVector) = ldiv!(adjoint(adjA.parent), b)
+ldiv!(A::Transpose{<:Any,<:Bidiagonal}, b::AbstractVector) = ldiv!(copy(A), b)
+ldiv!(A::Adjoint{<:Any,<:Bidiagonal}, b::AbstractVector) = ldiv!(copy(A), b)
 function ldiv!(A::Union{Bidiagonal,AbstractTriangular}, B::AbstractMatrix)
     nA,mA = size(A)
     tmp = similar(B,size(B,1))
@@ -527,7 +533,7 @@ function ldiv!(adjA::Adjoint{<:Any,<:Union{Bidiagonal,AbstractTriangular}}, B::A
     end
     for i = 1:size(B,2)
         copyto!(tmp, 1, B, (i - 1)*n + 1, n)
-        ldiv!(Adjoint(A), tmp)
+        ldiv!(adjoint(A), tmp)
         copyto!(B, (i - 1)*n + 1, tmp, 1, n) # Modify this when array view are implemented.
     end
     B
@@ -542,7 +548,7 @@ function ldiv!(transA::Transpose{<:Any,<:Union{Bidiagonal,AbstractTriangular}},
     end
     for i = 1:size(B,2)
         copyto!(tmp, 1, B, (i - 1)*n + 1, n)
-        ldiv!(Transpose(A), tmp)
+        ldiv!(transpose(A), tmp)
         copyto!(B, (i - 1)*n + 1, tmp, 1, n) # Modify this when array view are implemented.
     end
     B
@@ -580,16 +586,16 @@ function \(transA::Transpose{<:Number,<:Bidiagonal{<:Number}}, B::AbstractVecOrM
     A = transA.parent
     TA, TB = eltype(A), eltype(B)
     TAB = typeof((zero(TA)*zero(TB) + zero(TA)*zero(TB))/one(TA))
-    ldiv!(Transpose(convert(AbstractArray{TAB}, A)), copy_oftype(B, TAB))
+    ldiv!(transpose(convert(AbstractArray{TAB}, A)), copy_oftype(B, TAB))
 end
-\(transA::Transpose{<:Any,<:Bidiagonal}, B::AbstractVecOrMat) = ldiv!(Transpose(transA.parent), copy(B))
+\(transA::Transpose{<:Any,<:Bidiagonal}, B::AbstractVecOrMat) = ldiv!(transpose(transA.parent), copy(B))
 function \(adjA::Adjoint{<:Number,<:Bidiagonal{<:Number}}, B::AbstractVecOrMat{<:Number})
     A = adjA.parent
     TA, TB = eltype(A), eltype(B)
     TAB = typeof((zero(TA)*zero(TB) + zero(TA)*zero(TB))/one(TA))
-    ldiv!(Adjoint(convert(AbstractArray{TAB}, A)), copy_oftype(B, TAB))
+    ldiv!(adjoint(convert(AbstractArray{TAB}, A)), copy_oftype(B, TAB))
 end
-\(adjA::Adjoint{<:Any,<:Bidiagonal}, B::AbstractVecOrMat) = ldiv!(Adjoint(adjA.parent), copy(B))
+\(adjA::Adjoint{<:Any,<:Bidiagonal}, B::AbstractVecOrMat) = ldiv!(adjoint(adjA.parent), copy(B))
 
 factorize(A::Bidiagonal) = A
 

base/linalg/bidiag.jl

@@ -123,7 +123,7 @@ end
 
 ## Structure query functions
 
-issymmetric(A::BitMatrix) = size(A, 1)==size(A, 2) && count(!iszero, A - transpose(A))==0
+issymmetric(A::BitMatrix) = size(A, 1)==size(A, 2) && count(!iszero, A - copy(A'))==0
 ishermitian(A::BitMatrix) = issymmetric(A)
 
 function nonzero_chunks(chunks::Vector{UInt64}, pos0::Int, pos1::Int)
@@ -250,16 +250,19 @@ function put_8x8_chunk(Bc::Vector{UInt64}, i1::Int, i2::Int, x::UInt64, m::Int,
     return
 end
 
-function transpose(B::BitMatrix)
-    l1 = size(B, 1)
-    l2 = size(B, 2)
-    Bt = falses(l2, l1)
+adjoint(B::Union{BitVector,BitMatrix}) = Adjoint(B)
+transpose(B::Union{BitVector,BitMatrix}) = Transpose(B)
+Base.copy(B::Adjoint{Bool,BitMatrix}) = transpose!(falses(size(B)), B.parent)
+Base.copy(B::Transpose{Bool,BitMatrix}) = transpose!(falses(size(B)), B.parent)
+function transpose!(C::BitMatrix, B::BitMatrix)
+    @boundscheck size(C) == reverse(size(B)) || throw(DimensionMismatch())
+    l1, l2 = size(B)
 
     cgap1, cinc1 = Base._div64(l1), Base._mod64(l1)
     cgap2, cinc2 = Base._div64(l2), Base._mod64(l2)
 
     Bc = B.chunks
-    Btc = Bt.chunks
+    Cc = C.chunks
 
     nc = length(Bc)
 
@@ -278,10 +281,8 @@ function transpose(B::BitMatrix)
                 msk8_2 >>>= j + 7 - l2
             end
 
-            put_8x8_chunk(Btc, j, i, x, l2, cgap2, cinc2, nc, msk8_2)
+            put_8x8_chunk(Cc, j, i, x, l2, cgap2, cinc2, nc, msk8_2)
         end
     end
-    return Bt
+    return C
 end
-
-adjoint(B::Union{BitMatrix,BitVector}) = transpose(B)

base/linalg/bitarray.jl

@@ -1000,7 +1000,7 @@ end
 """
     syr!(uplo, alpha, x, A)
 
-Rank-1 update of the symmetric matrix `A` with vector `x` as `alpha*x*Transpose(x) + A`.
+Rank-1 update of the symmetric matrix `A` with vector `x` as `alpha*x*transpose(x) + A`.
 [`uplo`](@ref stdlib-blas-uplo) controls which triangle of `A` is updated. Returns `A`.
 """
 function syr! end
@@ -1243,7 +1243,7 @@ end
 """
     syrk!(uplo, trans, alpha, A, beta, C)
 
-Rank-k update of the symmetric matrix `C` as `alpha*A*Transpose(A) + beta*C` or `alpha*Transpose(A)*A +
+Rank-k update of the symmetric matrix `C` as `alpha*A*transpose(A) + beta*C` or `alpha*transpose(A)*A +
 beta*C` according to [`trans`](@ref stdlib-blas-trans).
 Only the [`uplo`](@ref stdlib-blas-uplo) triangle of `C` is used. Returns `C`.
 """
@@ -1254,7 +1254,7 @@ function syrk! end
 
 Returns either the upper triangle or the lower triangle of `A`,
 according to [`uplo`](@ref stdlib-blas-uplo),
-of `alpha*A*Transpose(A)` or `alpha*Transpose(A)*A`,
+of `alpha*A*transpose(A)` or `alpha*transpose(A)*A`,
 according to [`trans`](@ref stdlib-blas-trans).
 """
 function syrk end

base/linalg/blas.jl

@@ -43,7 +43,7 @@ end
 """
     bkfact(A, rook::Bool=false) -> BunchKaufman
 
-Compute the Bunch-Kaufman [^Bunch1977] factorization of a symmetric or Hermitian matrix `A` as ``P'*U*D*U'*P`` or ``P'*L*D*L'*P``, depending on which triangle is stored in `A`, and return a `BunchKaufman` object. Note that if `A` is complex symmetric then `U'` and `L'` denote the unconjugated transposes, i.e. `Transpose(U)` and `Transpose(L)`.
+Compute the Bunch-Kaufman [^Bunch1977] factorization of a symmetric or Hermitian matrix `A` as ``P'*U*D*U'*P`` or ``P'*L*D*L'*P``, depending on which triangle is stored in `A`, and return a `BunchKaufman` object. Note that if `A` is complex symmetric then `U'` and `L'` denote the unconjugated transposes, i.e. `transpose(U)` and `transpose(L)`.
 
 If `rook` is `true`, rook pivoting is used. If `rook` is false, rook pivoting is not used.
 
@@ -115,7 +115,7 @@ end
     getproperty(B::BunchKaufman, d::Symbol)
 
 Extract the factors of the Bunch-Kaufman factorization `B`. The factorization can take the
-two forms `P'*L*D*L'*P` or `P'*U*D*U'*P` (or `L*D*Transpose(L)` in the complex symmetric case)
+two forms `P'*L*D*L'*P` or `P'*U*D*U'*P` (or `L*D*transpose(L)` in the complex symmetric case)
 where `P` is a (symmetric) permutation matrix, `L` is a `UnitLowerTriangular` matrix, `U` is a
 `UnitUpperTriangular`, and `D` is a block diagonal symmetric or Hermitian matrix with
 1x1 or 2x2 blocks. The argument `d` can be
@@ -191,13 +191,13 @@ function getproperty(B::BunchKaufman{T}, d::Symbol) where {T<:BlasFloat}
             if getfield(B, :uplo) == 'L'
                 return UnitLowerTriangular(LUD)
             else
-                throw(ArgumentError("factorization is U*D*Transpose(U) but you requested L"))
+                throw(ArgumentError("factorization is U*D*transpose(U) but you requested L"))
             end
         else # :U
             if B.uplo == 'U'
                 return UnitUpperTriangular(LUD)
             else
-                throw(ArgumentError("factorization is L*D*Transpose(L) but you requested U"))
+                throw(ArgumentError("factorization is L*D*transpose(L) but you requested U"))
             end
         end
     else

base/linalg/bunchkaufman.jl

@@ -87,7 +87,7 @@ function _chol!(A::AbstractMatrix, ::Type{UpperTriangular})
                 return UpperTriangular(A), info
             end
             A[k,k] = Akk
-            AkkInv = inv(adjoint(Akk))
+            AkkInv = inv(copy(Akk'))
             for j = k + 1:n
                 for i = 1:k - 1
                     A[k,j] -= A[i,k]'A[i,j]
@@ -384,9 +384,9 @@ function getproperty(C::Cholesky, d::Symbol)
     Cfactors = getfield(C, :factors)
     Cuplo    = getfield(C, :uplo)
     if d == :U
-        return UpperTriangular(Symbol(Cuplo) == d ? Cfactors : adjoint(Cfactors))
+        return UpperTriangular(Symbol(Cuplo) == d ? Cfactors : copy(Cfactors'))
     elseif d == :L
-        return LowerTriangular(Symbol(Cuplo) == d ? Cfactors : adjoint(Cfactors))
+        return LowerTriangular(Symbol(Cuplo) == d ? Cfactors : copy(Cfactors'))
     elseif d == :UL
         return Symbol(Cuplo) == :U ? UpperTriangular(Cfactors) : LowerTriangular(Cfactors)
     else
@@ -397,9 +397,9 @@ function getproperty(C::CholeskyPivoted{T}, d::Symbol) where T<:BlasFloat
     Cfactors = getfield(C, :factors)
     Cuplo    = getfield(C, :uplo)
     if d == :U
-        return UpperTriangular(Symbol(Cuplo) == d ? Cfactors : adjoint(Cfactors))
+        return UpperTriangular(Symbol(Cuplo) == d ? Cfactors : copy(Cfactors'))
     elseif d == :L
-        return LowerTriangular(Symbol(Cuplo) == d ? Cfactors : adjoint(Cfactors))
+        return LowerTriangular(Symbol(Cuplo) == d ? Cfactors : copy(Cfactors'))
     elseif d == :p
         return getfield(C, :piv)
     elseif d == :P
@@ -437,9 +437,9 @@ ldiv!(C::Cholesky{T,<:AbstractMatrix}, B::StridedVecOrMat{T}) where {T<:BlasFloa
 
 function ldiv!(C::Cholesky{<:Any,<:AbstractMatrix}, B::StridedVecOrMat)
     if C.uplo == 'L'
-        return ldiv!(Adjoint(LowerTriangular(C.factors)), ldiv!(LowerTriangular(C.factors), B))
+        return ldiv!(adjoint(LowerTriangular(C.factors)), ldiv!(LowerTriangular(C.factors), B))
     else
-        return ldiv!(UpperTriangular(C.factors), ldiv!(Adjoint(UpperTriangular(C.factors)), B))
+        return ldiv!(UpperTriangular(C.factors), ldiv!(adjoint(UpperTriangular(C.factors)), B))
     end
 end
 
@@ -462,21 +462,21 @@ end
 
 function ldiv!(C::CholeskyPivoted, B::StridedVector)
     if C.uplo == 'L'
-        ldiv!(Adjoint(LowerTriangular(C.factors)),
+        ldiv!(adjoint(LowerTriangular(C.factors)),
             ldiv!(LowerTriangular(C.factors), B[C.piv]))[invperm(C.piv)]
     else
         ldiv!(UpperTriangular(C.factors),
-            ldiv!(Adjoint(UpperTriangular(C.factors)), B[C.piv]))[invperm(C.piv)]
+            ldiv!(adjoint(UpperTriangular(C.factors)), B[C.piv]))[invperm(C.piv)]
     end
 end
 
 function ldiv!(C::CholeskyPivoted, B::StridedMatrix)
     if C.uplo == 'L'
-        ldiv!(Adjoint(LowerTriangular(C.factors)),
+        ldiv!(adjoint(LowerTriangular(C.factors)),
             ldiv!(LowerTriangular(C.factors), B[C.piv,:]))[invperm(C.piv),:]
     else
         ldiv!(UpperTriangular(C.factors),
-            ldiv!(Adjoint(UpperTriangular(C.factors)), B[C.piv,:]))[invperm(C.piv),:]
+            ldiv!(adjoint(UpperTriangular(C.factors)), B[C.piv,:]))[invperm(C.piv),:]
     end
 end
 

base/linalg/cholesky.jl

@@ -1336,7 +1336,7 @@ function nullspace(A::StridedMatrix{T}) where T
     (m == 0 || n == 0) && return Matrix{T}(I, n, n)
     SVD = svdfact(A, full = true)
     indstart = sum(SVD.S .> max(m,n)*maximum(SVD.S)*eps(eltype(SVD.S))) + 1
-    return adjoint(SVD.Vt[indstart:end,:])
+    return copy(SVD.Vt[indstart:end,:]')
 end
 nullspace(a::StridedVector) = nullspace(reshape(a, length(a), 1))
 
@@ -1402,9 +1402,9 @@ function sylvester(A::StridedMatrix{T},B::StridedMatrix{T},C::StridedMatrix{T})
     RA, QA = schur(A)
     RB, QB = schur(B)
 
-    D = -(Adjoint(QA) * (C*QB))
+    D = -(adjoint(QA) * (C*QB))
     Y, scale = LAPACK.trsyl!('N','N', RA, RB, D)
-    scale!(QA*(Y * Adjoint(QB)), inv(scale))
+    scale!(QA*(Y * adjoint(QB)), inv(scale))
 end
 sylvester(A::StridedMatrix{T}, B::StridedMatrix{T}, C::StridedMatrix{T}) where {T<:Integer} = sylvester(float(A), float(B), float(C))
 
@@ -1445,9 +1445,9 @@ julia> A*X + X*A' + B
 function lyap(A::StridedMatrix{T}, C::StridedMatrix{T}) where {T<:BlasFloat}
     R, Q = schur(A)
 
-    D = -(Adjoint(Q) * (C*Q))
+    D = -(adjoint(Q) * (C*Q))
     Y, scale = LAPACK.trsyl!('N', T <: Complex ? 'C' : 'T', R, R, D)
-    scale!(Q*(Y * Adjoint(Q)), inv(scale))
+    scale!(Q*(Y * adjoint(Q)), inv(scale))
 end
 lyap(A::StridedMatrix{T}, C::StridedMatrix{T}) where {T<:Integer} = lyap(float(A), float(C))
 lyap(a::T, c::T) where {T<:Number} = -c/(2a)