diff --git a/base/abstractarraymath.jl b/base/abstractarraymath.jl index a924a083eddd9..7f343b6365d37 100644 --- a/base/abstractarraymath.jl +++ b/base/abstractarraymath.jl @@ -163,6 +163,18 @@ function cumsum_kbn{T<:AbstractFloat}(A::AbstractArray{T}, axis::Integer=1) return B + C end +## Permute array dims ## + +function permutedims(B::AbstractArray, perm) + dimsB = size(B) + ndimsB = length(dimsB) + (ndimsB == length(perm) && isperm(perm)) || throw(ArgumentError("no valid permutation of dimensions")) + dimsP = ntuple(i->dimsB[perm[i]], ndimsB)::typeof(dimsB) + P = similar(B, dimsP) + permutedims!(P, B, perm) +end + + ## ipermutedims in terms of permutedims ## function ipermutedims(A::AbstractArray,perm) diff --git a/base/array.jl b/base/array.jl index 89596506a1206..8e7fa1200977e 100644 --- a/base/array.jl +++ b/base/array.jl @@ -582,7 +582,6 @@ function lexcmp(a::Array{UInt8,1}, b::Array{UInt8,1}) c < 0 ? -1 : c > 0 ? +1 : cmp(length(a),length(b)) end -# note: probably should be StridedVector or AbstractVector function reverse(A::AbstractVector, s=1, n=length(A)) B = similar(A) for i = 1:s-1 @@ -598,9 +597,7 @@ function reverse(A::AbstractVector, s=1, n=length(A)) end reverseind(a::AbstractVector, i::Integer) = length(a) + 1 - i -reverse(v::StridedVector) = (n=length(v); [ v[n-i+1] for i=1:n ]) -reverse(v::StridedVector, s, n=length(v)) = reverse!(copy(v), s, n) -function reverse!(v::StridedVector, s=1, n=length(v)) +function reverse!(v::AbstractVector, s=1, n=length(v)) if n <= s # empty case; ok elseif !(1 ≤ s ≤ endof(v)) throw(BoundsError(v, s)) @@ -724,7 +721,7 @@ function find(testf::Function, A::AbstractArray) I end -function find(A::StridedArray) +function find(A::AbstractArray) nnzA = countnz(A) I = similar(A, Int, nnzA) count = 1 @@ -742,7 +739,7 @@ find(testf::Function, x::Number) = !testf(x) ? Array(Int,0) : [1] findn(A::AbstractVector) = find(A) -function findn(A::StridedMatrix) +function findn(A::AbstractMatrix) nnzA = countnz(A) I = similar(A, Int, nnzA) J = similar(A, Int, nnzA) diff --git a/base/arraymath.jl b/base/arraymath.jl index cefffb73c9d65..26cd5cdcb5cc8 100644 --- a/base/arraymath.jl +++ b/base/arraymath.jl @@ -11,7 +11,7 @@ end for f in (:-, :~, :conj, :sign) @eval begin - function ($f)(A::StridedArray) + function ($f)(A::AbstractArray) F = similar(A) for i in eachindex(A) F[i] = ($f)(A[i]) @@ -21,14 +21,12 @@ for f in (:-, :~, :conj, :sign) end end -(-)(A::StridedArray{Bool}) = reshape([ -A[i] for i in eachindex(A) ], size(A)) +(-)(A::AbstractArray{Bool}) = reshape([ -A[i] for i in eachindex(A) ], size(A)) -real(A::StridedArray) = reshape([ real(x) for x in A ], size(A)) -imag(A::StridedArray) = reshape([ imag(x) for x in A ], size(A)) -real{T<:Real}(x::StridedArray{T}) = x -imag{T<:Real}(x::StridedArray{T}) = zero(x) +real(A::AbstractArray) = reshape([ real(x) for x in A ], size(A)) +imag(A::AbstractArray) = reshape([ imag(x) for x in A ], size(A)) -function !(A::StridedArray{Bool}) +function !(A::AbstractArray{Bool}) F = similar(A) for i in eachindex(A) F[i] = !A[i] @@ -42,6 +40,7 @@ promote_array_type{Scalar, Arry}(F, ::Type{Scalar}, ::Type{Arry}) = promote_op(F promote_array_type{S<:Real, A<:AbstractFloat}(F, ::Type{S}, ::Type{A}) = A promote_array_type{S<:Integer, A<:Integer}(F, ::Type{S}, ::Type{A}) = A promote_array_type{S<:Integer}(F, ::Type{S}, ::Type{Bool}) = S +promote_array_type(F, ::Type{Bool}, ::Type{Bool}) = promote_op(F, Bool, Bool) # Handle operations that return different types ./(x::Number, Y::AbstractArray) = @@ -142,39 +141,9 @@ end (-)(A::AbstractArray,x::Number) = A .- x (-)(x::Number,A::AbstractArray) = x .- A -# functions that should give an Int result for Bool arrays -for f in (:.+, :.-) - @eval begin - function ($f)(A::Bool, B::StridedArray{Bool}) - F = similar(B, Int, size(B)) - for i in eachindex(B) - @inbounds F[i] = ($f)(A, B[i]) - end - return F - end - function ($f)(A::StridedArray{Bool}, B::Bool) - F = similar(A, Int, size(A)) - for i in eachindex(A) - @inbounds F[i] = ($f)(A[i], B) - end - return F - end - end -end -for f in (:+, :-) - @eval begin - function ($f)(A::StridedArray{Bool}, B::StridedArray{Bool}) - F = similar(A, Int, promote_shape(size(A), size(B))) - for i in eachindex(A,B) - @inbounds F[i] = ($f)(A[i], B[i]) - end - return F - end - end -end - ## data movement ## +# TODO?: replace with slice? function slicedim(A::Array, d::Integer, i::Integer) if d < 1 throw(ArgumentError("dimension must be ≥ 1")) @@ -272,7 +241,7 @@ function flipdim{T}(A::Array{T}, d::Integer) return B end -function rotl90(A::StridedMatrix) +function rotl90(A::AbstractMatrix) m,n = size(A) B = similar(A,(n,m)) for i=1:m, j=1:n @@ -280,7 +249,7 @@ function rotl90(A::StridedMatrix) end return B end -function rotr90(A::StridedMatrix) +function rotr90(A::AbstractMatrix) m,n = size(A) B = similar(A,(n,m)) for i=1:m, j=1:n @@ -288,7 +257,7 @@ function rotr90(A::StridedMatrix) end return B end -function rot180(A::StridedMatrix) +function rot180(A::AbstractMatrix) m,n = size(A) B = similar(A) for i=1:m, j=1:n @@ -308,7 +277,7 @@ rot180(A::AbstractMatrix, k::Integer) = mod(k, 2) == 1 ? rot180(A) : copy(A) ## Transpose ## const transposebaselength=64 -function transpose!(B::StridedMatrix,A::StridedMatrix) +function transpose!(B::AbstractMatrix,A::AbstractMatrix) m, n = size(A) size(B,1) == n && size(B,2) == m || throw(DimensionMismatch("transpose")) @@ -325,15 +294,15 @@ function transpose!(B::StridedMatrix,A::StridedMatrix) end return B end -function transpose!(B::StridedVector, A::StridedMatrix) +function transpose!(B::AbstractVector, A::AbstractMatrix) length(B) == length(A) && size(A,1) == 1 || throw(DimensionMismatch("transpose")) copy!(B, A) end -function transpose!(B::StridedMatrix, A::StridedVector) +function transpose!(B::AbstractMatrix, A::AbstractVector) length(B) == length(A) && size(B,1) == 1 || throw(DimensionMismatch("transpose")) copy!(B, A) end -function transposeblock!(B::StridedMatrix,A::StridedMatrix,m::Int,n::Int,offseti::Int,offsetj::Int) +function transposeblock!(B::AbstractMatrix,A::AbstractMatrix,m::Int,n::Int,offseti::Int,offsetj::Int) if m*n<=transposebaselength @inbounds begin for j = offsetj+(1:n) @@ -353,7 +322,7 @@ function transposeblock!(B::StridedMatrix,A::StridedMatrix,m::Int,n::Int,offseti end return B end -function ctranspose!(B::StridedMatrix,A::StridedMatrix) +function ctranspose!(B::AbstractMatrix,A::AbstractMatrix) m, n = size(A) size(B,1) == n && size(B,2) == m || throw(DimensionMismatch("transpose")) @@ -370,15 +339,15 @@ function ctranspose!(B::StridedMatrix,A::StridedMatrix) end return B end -function ctranspose!(B::StridedVector, A::StridedMatrix) +function ctranspose!(B::AbstractVector, A::AbstractMatrix) length(B) == length(A) && size(A,1) == 1 || throw(DimensionMismatch("transpose")) ccopy!(B, A) end -function ctranspose!(B::StridedMatrix, A::StridedVector) +function ctranspose!(B::AbstractMatrix, A::AbstractVector) length(B) == length(A) && size(B,1) == 1 || throw(DimensionMismatch("transpose")) ccopy!(B, A) end -function ctransposeblock!(B::StridedMatrix,A::StridedMatrix,m::Int,n::Int,offseti::Int,offsetj::Int) +function ctransposeblock!(B::AbstractMatrix,A::AbstractMatrix,m::Int,n::Int,offseti::Int,offsetj::Int) if m*n<=transposebaselength @inbounds begin for j = offsetj+(1:n) @@ -404,18 +373,18 @@ function ccopy!(B, A) end end -function transpose(A::StridedMatrix) +function transpose(A::AbstractMatrix) B = similar(A, size(A, 2), size(A, 1)) transpose!(B, A) end -function ctranspose(A::StridedMatrix) +function ctranspose(A::AbstractMatrix) B = similar(A, size(A, 2), size(A, 1)) ctranspose!(B, A) end -ctranspose{T<:Real}(A::StridedVecOrMat{T}) = transpose(A) +ctranspose{T<:Real}(A::AbstractVecOrMat{T}) = transpose(A) -transpose(x::StridedVector) = [ transpose(x[j]) for i=1, j=1:size(x,1) ] -ctranspose{T}(x::StridedVector{T}) = T[ ctranspose(x[j]) for i=1, j=1:size(x,1) ] +transpose(x::AbstractVector) = [ transpose(x[j]) for i=1, j=1:size(x,1) ] +ctranspose{T}(x::AbstractVector{T}) = T[ ctranspose(x[j]) for i=1, j=1:size(x,1) ] _cumsum_type{T<:Number}(v::AbstractArray{T}) = typeof(+zero(T)) _cumsum_type(v) = typeof(v[1]+v[1]) @@ -452,45 +421,3 @@ for (f, f!, fp, op) = ((:cumsum, :cumsum!, :cumsum_pairwise!, :+), return ($f!)(c, v) end end - -for (f, op) = ((:cummin, :min), (:cummax, :max)) - @eval function ($f)(v::AbstractVector) - n = length(v) - cur_val = v[1] - res = similar(v, n) - res[1] = cur_val - for i in 2:n - cur_val = ($op)(v[i], cur_val) - res[i] = cur_val - end - return res - end - - @eval function ($f)(A::StridedArray, axis::Integer) - dimsA = size(A) - ndimsA = ndims(A) - axis_size = dimsA[axis] - axis_stride = 1 - for i = 1:(axis-1) - axis_stride *= size(A,i) - end - - if axis_size < 1 - return A - end - - B = similar(A) - - for i = 1:length(A) - if div(i-1, axis_stride) % axis_size == 0 - B[i] = A[i] - else - B[i] = ($op)(A[i], B[i-axis_stride]) - end - end - - return B - end - - @eval ($f)(A::AbstractArray) = ($f)(A, 1) -end diff --git a/base/bitarray.jl b/base/bitarray.jl index d5acc46e19310..57cab47eb34d1 100644 --- a/base/bitarray.jl +++ b/base/bitarray.jl @@ -1641,18 +1641,6 @@ end ctranspose(B::BitArray) = transpose(B) -## Permute array dims ## - -function permutedims(B::Union{BitArray,StridedArray}, perm) - dimsB = size(B) - ndimsB = length(dimsB) - (ndimsB == length(perm) && isperm(perm)) || throw(ArgumentError("no valid permutation of dimensions")) - dimsP = ntuple(i->dimsB[perm[i]], ndimsB)::typeof(dimsB) - P = similar(B, dimsP) - permutedims!(P, B, perm) -end - - ## Concatenation ## function hcat(B::BitVector...) diff --git a/base/coreimg.jl b/base/coreimg.jl index fcb6eec947442..47a2be523a02a 100644 --- a/base/coreimg.jl +++ b/base/coreimg.jl @@ -47,10 +47,6 @@ const checked_sub = - # core array operations include("abstractarray.jl") -typealias StridedArray{T,N,A<:DenseArray,I<:Tuple{Vararg{RangeIndex}}} DenseArray{T,N} -typealias StridedVector{T,A<:DenseArray,I<:Tuple{Vararg{RangeIndex}}} DenseArray{T,1} -typealias StridedMatrix{T,A<:DenseArray,I<:Tuple{Vararg{RangeIndex}}} DenseArray{T,2} -typealias StridedVecOrMat{T} Union{StridedVector{T}, StridedMatrix{T}} include("array.jl") #TODO: eliminate Dict from inference diff --git a/base/linalg/matmul.jl b/base/linalg/matmul.jl index 8b2f123e9066a..93f87f8ae51ac 100644 --- a/base/linalg/matmul.jl +++ b/base/linalg/matmul.jl @@ -96,35 +96,35 @@ for elty in (Float32,Float64) end end end -A_mul_B!(y::StridedVector, A::StridedVecOrMat, x::StridedVector) = generic_matvecmul!(y, 'N', A, x) +A_mul_B!(y::AbstractVector, A::AbstractVecOrMat, x::AbstractVector) = generic_matvecmul!(y, 'N', A, x) function At_mul_B{T<:BlasFloat,S}(A::StridedMatrix{T}, x::StridedVector{S}) TS = promote_op(MulFun(),arithtype(T),arithtype(S)) At_mul_B!(similar(x,TS,size(A,2)), A, convert(AbstractVector{TS}, x)) end -function At_mul_B{T,S}(A::StridedMatrix{T}, x::StridedVector{S}) +function At_mul_B{T,S}(A::AbstractMatrix{T}, x::AbstractVector{S}) TS = promote_op(MulFun(),arithtype(T),arithtype(S)) At_mul_B!(similar(x,TS,size(A,2)), A, x) end At_mul_B!{T<:BlasFloat}(y::StridedVector{T}, A::StridedVecOrMat{T}, x::StridedVector{T}) = gemv!(y, 'T', A, x) -At_mul_B!(y::StridedVector, A::StridedVecOrMat, x::StridedVector) = generic_matvecmul!(y, 'T', A, x) +At_mul_B!(y::AbstractVector, A::AbstractVecOrMat, x::AbstractVector) = generic_matvecmul!(y, 'T', A, x) function Ac_mul_B{T<:BlasFloat,S}(A::StridedMatrix{T}, x::StridedVector{S}) TS = promote_op(MulFun(),arithtype(T),arithtype(S)) Ac_mul_B!(similar(x,TS,size(A,2)),A,convert(AbstractVector{TS},x)) end -function Ac_mul_B{T,S}(A::StridedMatrix{T}, x::StridedVector{S}) +function Ac_mul_B{T,S}(A::AbstractMatrix{T}, x::AbstractVector{S}) TS = promote_op(MulFun(),arithtype(T),arithtype(S)) Ac_mul_B!(similar(x,TS,size(A,2)), A, x) end Ac_mul_B!{T<:BlasReal}(y::StridedVector{T}, A::StridedVecOrMat{T}, x::StridedVector{T}) = At_mul_B!(y, A, x) Ac_mul_B!{T<:BlasComplex}(y::StridedVector{T}, A::StridedVecOrMat{T}, x::StridedVector{T}) = gemv!(y, 'C', A, x) -Ac_mul_B!(y::StridedVector, A::StridedVecOrMat, x::StridedVector) = generic_matvecmul!(y, 'C', A, x) +Ac_mul_B!(y::AbstractVector, A::AbstractVecOrMat, x::AbstractVector) = generic_matvecmul!(y, 'C', A, x) # Matrix-matrix multiplication -function (*){T,S}(A::AbstractMatrix{T}, B::StridedMatrix{S}) +function (*){T,S}(A::AbstractMatrix{T}, B::AbstractMatrix{S}) TS = promote_op(MulFun(), arithtype(T), arithtype(S)) A_mul_B!(similar(B, TS, (size(A,1), size(B,2))), A, B) end @@ -139,16 +139,16 @@ for elty in (Float32,Float64) end end end -A_mul_B!(C::StridedMatrix, A::StridedVecOrMat, B::StridedVecOrMat) = generic_matmatmul!(C, 'N', 'N', A, B) +A_mul_B!(C::AbstractMatrix, A::AbstractVecOrMat, B::AbstractVecOrMat) = generic_matmatmul!(C, 'N', 'N', A, B) -function At_mul_B{T,S}(A::AbstractMatrix{T}, B::StridedMatrix{S}) +function At_mul_B{T,S}(A::AbstractMatrix{T}, B::AbstractMatrix{S}) TS = promote_op(MulFun(),arithtype(T), arithtype(S)) At_mul_B!(similar(B, TS, (size(A,2), size(B,2))), A, B) end At_mul_B!{T<:BlasFloat}(C::StridedMatrix{T}, A::StridedVecOrMat{T}, B::StridedVecOrMat{T}) = is(A,B) ? syrk_wrapper!(C, 'T', A) : gemm_wrapper!(C, 'T', 'N', A, B) -At_mul_B!(C::StridedMatrix, A::StridedVecOrMat, B::StridedVecOrMat) = generic_matmatmul!(C, 'T', 'N', A, B) +At_mul_B!(C::AbstractMatrix, A::AbstractVecOrMat, B::AbstractVecOrMat) = generic_matmatmul!(C, 'T', 'N', A, B) -function A_mul_Bt{T,S}(A::AbstractMatrix{T}, B::StridedMatrix{S}) +function A_mul_Bt{T,S}(A::AbstractMatrix{T}, B::AbstractMatrix{S}) TS = promote_op(MulFun(),arithtype(T), arithtype(S)) A_mul_Bt!(similar(B, TS, (size(A,1), size(B,1))), A, B) end @@ -163,42 +163,42 @@ for elty in (Float32,Float64) end end end -A_mul_Bt!(C::StridedVecOrMat, A::StridedVecOrMat, B::StridedVecOrMat) = generic_matmatmul!(C, 'N', 'T', A, B) +A_mul_Bt!(C::AbstractVecOrMat, A::AbstractVecOrMat, B::AbstractVecOrMat) = generic_matmatmul!(C, 'N', 'T', A, B) -function At_mul_Bt{T,S}(A::AbstractMatrix{T}, B::StridedVecOrMat{S}) +function At_mul_Bt{T,S}(A::AbstractMatrix{T}, B::AbstractVecOrMat{S}) TS = promote_op(MulFun(),arithtype(T), arithtype(S)) At_mul_Bt!(similar(B, TS, (size(A,2), size(B,1))), A, B) end At_mul_Bt!{T<:BlasFloat}(C::StridedMatrix{T}, A::StridedVecOrMat{T}, B::StridedVecOrMat{T}) = gemm_wrapper!(C, 'T', 'T', A, B) -At_mul_Bt!(C::StridedMatrix, A::StridedVecOrMat, B::StridedVecOrMat) = generic_matmatmul!(C, 'T', 'T', A, B) +At_mul_Bt!(C::AbstractMatrix, A::AbstractVecOrMat, B::AbstractVecOrMat) = generic_matmatmul!(C, 'T', 'T', A, B) Ac_mul_B{T<:BlasReal}(A::StridedMatrix{T}, B::StridedMatrix{T}) = At_mul_B(A, B) Ac_mul_B!{T<:BlasReal}(C::StridedMatrix{T}, A::StridedVecOrMat{T}, B::StridedVecOrMat{T}) = At_mul_B!(C, A, B) -function Ac_mul_B{T,S}(A::StridedMatrix{T}, B::StridedMatrix{S}) +function Ac_mul_B{T,S}(A::AbstractMatrix{T}, B::AbstractMatrix{S}) TS = promote_op(MulFun(),arithtype(T), arithtype(S)) Ac_mul_B!(similar(B, TS, (size(A,2), size(B,2))), A, B) end Ac_mul_B!{T<:BlasComplex}(C::StridedMatrix{T}, A::StridedVecOrMat{T}, B::StridedVecOrMat{T}) = is(A,B) ? herk_wrapper!(C,'C',A) : gemm_wrapper!(C,'C', 'N', A, B) -Ac_mul_B!(C::StridedMatrix, A::StridedVecOrMat, B::StridedVecOrMat) = generic_matmatmul!(C, 'C', 'N', A, B) +Ac_mul_B!(C::AbstractMatrix, A::AbstractVecOrMat, B::AbstractVecOrMat) = generic_matmatmul!(C, 'C', 'N', A, B) A_mul_Bc{T<:BlasFloat,S<:BlasReal}(A::StridedMatrix{T}, B::StridedMatrix{S}) = A_mul_Bt(A, B) A_mul_Bc!{T<:BlasFloat,S<:BlasReal}(C::StridedMatrix{T}, A::StridedVecOrMat{T}, B::StridedVecOrMat{S}) = A_mul_Bt!(C, A, B) -function A_mul_Bc{T,S}(A::StridedMatrix{T}, B::StridedMatrix{S}) +function A_mul_Bc{T,S}(A::AbstractMatrix{T}, B::AbstractMatrix{S}) TS = promote_op(MulFun(),arithtype(T),arithtype(S)) A_mul_Bc!(similar(B,TS,(size(A,1),size(B,1))),A,B) end A_mul_Bc!{T<:BlasComplex}(C::StridedMatrix{T}, A::StridedVecOrMat{T}, B::StridedVecOrMat{T}) = is(A,B) ? herk_wrapper!(C, 'N', A) : gemm_wrapper!(C, 'N', 'C', A, B) -A_mul_Bc!(C::StridedMatrix, A::StridedVecOrMat, B::StridedVecOrMat) = generic_matmatmul!(C, 'N', 'C', A, B) +A_mul_Bc!(C::AbstractMatrix, A::AbstractVecOrMat, B::AbstractVecOrMat) = generic_matmatmul!(C, 'N', 'C', A, B) -Ac_mul_Bc{T,S}(A::AbstractMatrix{T}, B::StridedMatrix{S}) = Ac_mul_Bc!(similar(B, promote_op(MulFun(),arithtype(T), arithtype(S)), (size(A,2), size(B,1))), A, B) +Ac_mul_Bc{T,S}(A::AbstractMatrix{T}, B::AbstractMatrix{S}) = Ac_mul_Bc!(similar(B, promote_op(MulFun(),arithtype(T), arithtype(S)), (size(A,2), size(B,1))), A, B) Ac_mul_Bc!{T<:BlasFloat}(C::StridedMatrix{T}, A::StridedVecOrMat{T}, B::StridedVecOrMat{T}) = gemm_wrapper!(C, 'C', 'C', A, B) -Ac_mul_Bc!(C::StridedMatrix, A::StridedVecOrMat, B::StridedVecOrMat) = generic_matmatmul!(C, 'C', 'C', A, B) -Ac_mul_Bt{T,S}(A::AbstractMatrix{T}, B::StridedMatrix{S}) = Ac_mul_Bt(similar(B, promote_op(MulFun(),arithtype(A), arithtype(B)), (size(A,2), size(B,1))), A, B) -Ac_mul_Bt!(C::StridedMatrix, A::StridedVecOrMat, B::StridedVecOrMat) = generic_matmatmul!(C, 'C', 'T', A, B) +Ac_mul_Bc!(C::AbstractMatrix, A::AbstractVecOrMat, B::AbstractVecOrMat) = generic_matmatmul!(C, 'C', 'C', A, B) +Ac_mul_Bt{T,S}(A::AbstractMatrix{T}, B::AbstractMatrix{S}) = Ac_mul_Bt(similar(B, promote_op(MulFun(),arithtype(A), arithtype(B)), (size(A,2), size(B,1))), A, B) +Ac_mul_Bt!(C::AbstractMatrix, A::AbstractVecOrMat, B::AbstractVecOrMat) = generic_matmatmul!(C, 'C', 'T', A, B) # Supporting functions for matrix multiplication -function copytri!(A::StridedMatrix, uplo::Char, conjugate::Bool=false) +function copytri!(A::AbstractMatrix, uplo::Char, conjugate::Bool=false) n = checksquare(A) if uplo == 'U' for i = 1:(n-1), j = (i+1):n diff --git a/base/linalg/qr.jl b/base/linalg/qr.jl index b4ec2b2582499..e98ec6e87c69c 100644 --- a/base/linalg/qr.jl +++ b/base/linalg/qr.jl @@ -390,7 +390,7 @@ function A_mul_Bc!{T}(A::AbstractMatrix{T},Q::QRPackedQ{T}) end A end -function A_mul_Bc{TA,TB}(A::AbstractArray{TA}, B::Union{QRCompactWYQ{TB},QRPackedQ{TB}}) +function A_mul_Bc{TA,TB}(A::AbstractMatrix{TA}, B::Union{QRCompactWYQ{TB},QRPackedQ{TB}}) TAB = promote_type(TA,TB) BB = convert(AbstractMatrix{TAB}, B) if size(A,2) == size(B.factors, 1) diff --git a/base/linalg/triangular.jl b/base/linalg/triangular.jl index 09688dd10e4aa..34248e24b0a67 100644 --- a/base/linalg/triangular.jl +++ b/base/linalg/triangular.jl @@ -376,8 +376,11 @@ scale!(c::Number, A::Union{UpperTriangular,LowerTriangular}) = scale!(A,c) ###################### A_mul_B!(A::Tridiagonal, B::AbstractTriangular) = A*full!(B) +A_mul_B!(C::AbstractVector, A::AbstractTriangular, B::AbstractVector) = A_mul_B!(A, copy!(C, B)) +A_mul_B!(C::AbstractMatrix, A::AbstractTriangular, B::AbstractVecOrMat) = A_mul_B!(A, copy!(C, B)) A_mul_B!(C::AbstractVecOrMat, A::AbstractTriangular, B::AbstractVecOrMat) = A_mul_B!(A, copy!(C, B)) A_mul_Bt!(C::AbstractVecOrMat, A::AbstractTriangular, B::AbstractVecOrMat) = A_mul_B!(A, transpose!(C, B)) +A_mul_Bc!(C::AbstractMatrix, A::AbstractTriangular, B::AbstractVecOrMat) = A_mul_B!(A, ctranspose!(C, B)) A_mul_Bc!(C::AbstractVecOrMat, A::AbstractTriangular, B::AbstractVecOrMat) = A_mul_B!(A, ctranspose!(C, B)) for (t, uploc, isunitc) in ((:LowerTriangular, 'L', 'N'), diff --git a/base/linalg/tridiag.jl b/base/linalg/tridiag.jl index 83288d3f48262..e2be1b5909233 100644 --- a/base/linalg/tridiag.jl +++ b/base/linalg/tridiag.jl @@ -498,7 +498,11 @@ function convert{T}(::Type{SymTridiagonal{T}}, M::Tridiagonal) end end -function A_mul_B!(C::AbstractVecOrMat, A::Tridiagonal, B::AbstractVecOrMat) +A_mul_B!(C::AbstractVector, A::Tridiagonal, B::AbstractVector) = A_mul_B_td!(C, A, B) +A_mul_B!(C::AbstractMatrix, A::Tridiagonal, B::AbstractVecOrMat) = A_mul_B_td!(C, A, B) +A_mul_B!(C::AbstractVecOrMat, A::Tridiagonal, B::AbstractVecOrMat) = A_mul_B_td!(C, A, B) + +function A_mul_B_td!(C::AbstractVecOrMat, A::Tridiagonal, B::AbstractVecOrMat) nA = size(A,1) nB = size(B,2) if !(size(C,1) == size(B,1) == nA) diff --git a/base/multidimensional.jl b/base/multidimensional.jl index 7ba32048f70b6..7c7a157f8aa4e 100644 --- a/base/multidimensional.jl +++ b/base/multidimensional.jl @@ -377,6 +377,46 @@ end end end +for (f, fmod, op) = ((:cummin, :_cummin!, :min), (:cummax, :_cummax!, :max)) + @eval function ($f)(v::AbstractVector) + n = length(v) + cur_val = v[1] + res = similar(v, n) + res[1] = cur_val + for i in 2:n + cur_val = ($op)(v[i], cur_val) + res[i] = cur_val + end + return res + end + + @eval function ($f)(A::AbstractArray, axis::Integer) + res = similar(A) + if size(A, axis) < 1 + return res + end + R1 = CartesianRange(size(A)[1:axis-1]) + R2 = CartesianRange(size(A)[axis+1:end]) + ($fmod)(res, A, R1, R2, axis) + end + + @eval @noinline function ($fmod)(res, A::AbstractArray, R1::CartesianRange, R2::CartesianRange, axis::Integer) + for I2 in R2 + for I1 in R1 + res[I1, 1, I2] = A[I1, 1, I2] + end + for i = 2:size(A, axis) + for I1 in R1 + res[I1, i, I2] = ($op)(A[I1, i, I2], res[I1, i-1, I2]) + end + end + end + res + end + + @eval ($f)(A::AbstractArray) = ($f)(A, 1) +end + ## SubArray index merging # A view created like V = A[2:3:8, 5:2:17] can later be indexed as V[2:7], # creating a new 1d view. diff --git a/base/sparse.jl b/base/sparse.jl index 8ea18173d8187..963abc1a946f3 100644 --- a/base/sparse.jl +++ b/base/sparse.jl @@ -7,7 +7,9 @@ using Base.Sort: Forward using Base.LinAlg: AbstractTriangular, PosDefException import Base: +, -, *, \, &, |, $, .+, .-, .*, ./, .\, .^, .<, .!=, == -import Base: A_mul_B!, Ac_mul_B, Ac_mul_B!, At_mul_B, At_mul_B!, At_ldiv_B, Ac_ldiv_B, A_ldiv_B! +import Base: A_mul_B!, Ac_mul_B, Ac_mul_B!, At_mul_B, At_mul_B! +import Base: A_mul_Bc, A_mul_Bt, Ac_mul_Bc, At_mul_Bt +import Base: At_ldiv_B, Ac_ldiv_B, A_ldiv_B! import Base.LinAlg: At_ldiv_B!, Ac_ldiv_B! import Base: @get!, acos, acosd, acot, acotd, acsch, asech, asin, asind, asinh, diff --git a/base/sparse/linalg.jl b/base/sparse/linalg.jl index ba38612573135..9a8612533dd27 100644 --- a/base/sparse/linalg.jl +++ b/base/sparse/linalg.jl @@ -26,11 +26,24 @@ increment{T<:Integer}(A::AbstractArray{T}) = increment!(copy(A)) ## sparse matrix multiplication function (*){TvA,TiA,TvB,TiB}(A::SparseMatrixCSC{TvA,TiA}, B::SparseMatrixCSC{TvB,TiB}) + (*)(sppromote(A, B)...) +end +for f in (:A_mul_Bt, :A_mul_Bc, + :At_mul_B, :Ac_mul_B, + :At_mul_Bt, :Ac_mul_Bc) + @eval begin + function ($f){TvA,TiA,TvB,TiB}(A::SparseMatrixCSC{TvA,TiA}, B::SparseMatrixCSC{TvB,TiB}) + ($f)(sppromote(A, B)...) + end + end +end + +function sppromote{TvA,TiA,TvB,TiB}(A::SparseMatrixCSC{TvA,TiA}, B::SparseMatrixCSC{TvB,TiB}) Tv = promote_type(TvA, TvB) Ti = promote_type(TiA, TiB) A = convert(SparseMatrixCSC{Tv,Ti}, A) B = convert(SparseMatrixCSC{Tv,Ti}, B) - A * B + A, B end # In matrix-vector multiplication, the correct orientation of the vector is assumed. @@ -109,6 +122,18 @@ end # http://dl.acm.org/citation.cfm?id=355796 (*){Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, B::SparseMatrixCSC{Tv,Ti}) = spmatmul(A,B) +for (f, opA, opB) in ((:A_mul_Bt, :identity, :transpose), + (:A_mul_Bc, :identity, :ctranspose), + (:At_mul_B, :transpose, :identity), + (:Ac_mul_B, :ctranspose, :identity), + (:At_mul_Bt, :transpose, :transpose), + (:Ac_mul_Bc, :ctranspose, :ctranspose)) + @eval begin + function ($f){Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, B::SparseMatrixCSC{Tv,Ti}) + spmatmul(($opA)(A), ($opB)(B)) + end + end +end function spmatmul{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, B::SparseMatrixCSC{Tv,Ti}; sortindices::Symbol = :sortcols) diff --git a/base/subarray.jl b/base/subarray.jl index 2ae359c29acf7..57d73015d807e 100644 --- a/base/subarray.jl +++ b/base/subarray.jl @@ -48,10 +48,6 @@ function getindex(N::NoSlice, r::Range{Int}) end abstract AbstractCartesianIndex{N} # This is a hacky forward declaration for CartesianIndex -typealias StridedArray{T,N,A<:DenseArray,I<:Tuple{Vararg{Union{RangeIndex, NoSlice, AbstractCartesianIndex}}}} Union{DenseArray{T,N}, SubArray{T,N,A,I}} -typealias StridedVector{T,A<:DenseArray,I<:Tuple{Vararg{Union{RangeIndex, NoSlice, AbstractCartesianIndex}}}} Union{DenseArray{T,1}, SubArray{T,1,A,I}} -typealias StridedMatrix{T,A<:DenseArray,I<:Tuple{Vararg{Union{RangeIndex, NoSlice, AbstractCartesianIndex}}}} Union{DenseArray{T,2}, SubArray{T,2,A,I}} -typealias StridedVecOrMat{T} Union{StridedVector{T}, StridedMatrix{T}} # This computes the linear indexing compatability for a given tuple of indices viewindexing() = LinearFast() diff --git a/base/sysimg.jl b/base/sysimg.jl index 0090bda16acf7..30ae9d4dab9dd 100644 --- a/base/sysimg.jl +++ b/base/sysimg.jl @@ -101,6 +101,12 @@ include("dict.jl") include("set.jl") include("iterator.jl") +# StridedArrays +typealias StridedArray{T,N,A<:DenseArray,I<:Tuple{Vararg{Union{RangeIndex, NoSlice, AbstractCartesianIndex}}}} Union{DenseArray{T,N}, SubArray{T,N,A,I}} +typealias StridedVector{T,A<:DenseArray,I<:Tuple{Vararg{Union{RangeIndex, NoSlice, AbstractCartesianIndex}}}} Union{DenseArray{T,1}, SubArray{T,1,A,I}} +typealias StridedMatrix{T,A<:DenseArray,I<:Tuple{Vararg{Union{RangeIndex, NoSlice, AbstractCartesianIndex}}}} Union{DenseArray{T,2}, SubArray{T,2,A,I}} +typealias StridedVecOrMat{T} Union{StridedVector{T}, StridedMatrix{T}} + # For OS specific stuff include(UTF8String(vcat(length(Core.ARGS)>=2?Core.ARGS[2].data:"".data, "build_h.jl".data))) # include($BUILDROOT/base/build_h.jl) include(UTF8String(vcat(length(Core.ARGS)>=2?Core.ARGS[2].data:"".data, "version_git.jl".data))) # include($BUILDROOT/base/version_git.jl) diff --git a/test/linalg/matmul.jl b/test/linalg/matmul.jl index 1b21f6d215afb..7885d850e0f85 100644 --- a/test/linalg/matmul.jl +++ b/test/linalg/matmul.jl @@ -65,6 +65,21 @@ v = [1,2] C = Array(Int, 2, 2) @test @inferred(A_mul_Bc!(C, v, v)) == [1 2; 2 4] +# Generic AbstractArrays +module MyArray15367 + using Base.Test + immutable MyArray{T,N} <: AbstractArray{T,N} + data::Array{T,N} + end + + Base.size(A::MyArray) = size(A.data) + Base.getindex(A::MyArray, indexes...) = A.data[indexes...] + + A = MyArray(rand(4,5)) + b = rand(5) + @test_approx_eq A*b A.data*b +end + # Preallocated C = Array(Int, size(A, 1), size(B, 2)) @test A_mul_B!(C, A, B) == A*B diff --git a/test/sparsedir/sparse.jl b/test/sparsedir/sparse.jl index 56f498d6484da..436bd549825b4 100644 --- a/test/sparsedir/sparse.jl +++ b/test/sparsedir/sparse.jl @@ -272,6 +272,11 @@ end cA = sprandn(5,5,0.2) + im*sprandn(5,5,0.2) @test full(conj(cA)) == conj(full(cA)) +# transpose of SubArrays +A = sub(sprandn(10, 10, 0.3), 1:4, 1:4) +@test transpose(full(A)) == full(transpose(A)) +@test ctranspose(full(A)) == full(ctranspose(A)) + # exp A = sprandn(5,5,0.2) @test e.^A ≈ e.^full(A) @@ -1187,9 +1192,6 @@ Ari = ceil(Int64,100*Ar) @test_throws ArgumentError Base.SparseArrays.normestinv(Ac,21) @test_throws DimensionMismatch Base.SparseArrays.normestinv(sprand(3,5,.9)) -@test_throws ErrorException transpose(sub(sprandn(10, 10, 0.3), 1:4, 1:4)) -@test_throws ErrorException ctranspose(sub(sprandn(10, 10, 0.3), 1:4, 1:4)) - # csc_permute A = sprand(10,10,0.2) p = randperm(10)