Skip to content

Commit

Permalink
change dim arguments for diff and unique to keyword args (#26776)
Browse files Browse the repository at this point in the history
  • Loading branch information
JeffBezanson authored and mbauman committed Apr 12, 2018
1 parent 4131580 commit 1508942
Show file tree
Hide file tree
Showing 9 changed files with 50 additions and 37 deletions.
6 changes: 5 additions & 1 deletion base/deprecated.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1368,6 +1368,9 @@ export readandwrite

@deprecate squeeze(A, dims) squeeze(A, dims=dims)

@deprecate diff(A::AbstractMatrix, dim::Integer) diff(A, dims=dim)
@deprecate unique(A::AbstractArray, dim::Int) unique(A, dims=dim)

# PR #25196
@deprecate_binding ObjectIdDict IdDict{Any,Any}

Expand Down Expand Up @@ -1528,7 +1531,8 @@ end
@deprecate(matchall(r::Regex, s::AbstractString; overlap::Bool = false),
collect(m.match for m in eachmatch(r, s, overlap = overlap)))

@deprecate diff(A::AbstractMatrix) diff(A, 1)
# remove depwarn for `diff` in multidimensional.jl
# @deprecate diff(A::AbstractMatrix) diff(A, dims=1)

# PR 26194
export assert
Expand Down
32 changes: 20 additions & 12 deletions base/multidimensional.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1083,10 +1083,10 @@ diff(a::AbstractVector) = [ a[i+1] - a[i] for i=1:length(a)-1 ]

"""
diff(A::AbstractVector)
diff(A::AbstractMatrix, dim::Integer)
diff(A::AbstractMatrix; dims::Integer)
Finite difference operator of matrix or vector `A`. If `A` is a matrix,
specify the dimension over which to operate with the `dim` argument.
specify the dimension over which to operate with the `dims` keyword argument.
# Examples
```jldoctest
Expand All @@ -1095,7 +1095,7 @@ julia> a = [2 4; 6 16]
2 4
6 16
julia> diff(a,2)
julia> diff(a, dims=2)
2×1 Array{Int64,2}:
2
10
Expand All @@ -1107,13 +1107,17 @@ julia> diff(vec(a))
12
```
"""
function diff(A::AbstractMatrix, dim::Integer)
if dim == 1
function diff(A::AbstractMatrix; dims::Union{Integer,Nothing}=nothing)
if dims === nothing
depwarn("`diff(A::AbstractMatrix)` is deprecated, use `diff(A, dims=1)` instead.", :diff)
dims = 1
end
if dims == 1
[A[i+1,j] - A[i,j] for i=1:size(A,1)-1, j=1:size(A,2)]
elseif dim == 2
elseif dims == 2
[A[i,j+1] - A[i,j] for i=1:size(A,1), j=1:size(A,2)-1]
else
throw(ArgumentError("dimension dim must be 1 or 2, got $dim"))
throw(ArgumentError("dimension must be 1 or 2, got $dims"))
end
end

Expand Down Expand Up @@ -1724,9 +1728,9 @@ end
hash(x::Prehashed) = x.hash

"""
unique(A::AbstractArray, dim::Int)
unique(A::AbstractArray; dims::Int)
Return unique regions of `A` along dimension `dim`.
Return unique regions of `A` along dimension `dims`.
# Examples
```jldoctest
Expand All @@ -1745,7 +1749,7 @@ julia> unique(A)
true
false
julia> unique(A, 2)
julia> unique(A, dims=2)
2×1×2 Array{Bool,3}:
[:, :, 1] =
true
Expand All @@ -1755,14 +1759,18 @@ julia> unique(A, 2)
true
false
julia> unique(A, 3)
julia> unique(A, dims=3)
2×2×1 Array{Bool,3}:
[:, :, 1] =
true true
false false
```
"""
@generated function unique(A::AbstractArray{T,N}, dim::Int) where {T,N}
unique(A::AbstractArray; dims::Union{Colon,Integer} = :) = _unique_dims(A, dims)

_unique_dims(A::AbstractArray, dims::Colon) = invoke(unique, Tuple{Any}, A)

@generated function _unique_dims(A::AbstractArray{T,N}, dim::Integer) where {T,N}
inds = inds -> zeros(UInt, inds)
quote
1 <= dim <= $N || return copy(A)
Expand Down
2 changes: 1 addition & 1 deletion stdlib/LinearAlgebra/test/lapack.jl
Original file line number Diff line number Diff line change
Expand Up @@ -202,7 +202,7 @@ end
Bvs = eigvecs(B)
Avs = eigvecs(A)
Bvs = LAPACK.gebak!('S','R',ilo,ihi,scale,Bvs)
@test norm(diff(Avs ./ Bvs, 1)) < 100 * eps(abs(float(one(elty))))
@test norm(diff(Avs ./ Bvs, dims=1)) < 100 * eps(abs(float(one(elty))))
end
end

Expand Down
1 change: 1 addition & 0 deletions stdlib/SparseArrays/src/deprecated.jl
Original file line number Diff line number Diff line change
Expand Up @@ -233,6 +233,7 @@ Base.@deprecate_binding blkdiag blockdiag
@deprecate complex(x::AbstractVector{<:Real}, y::AbstractSparseVector{<:Real}) complex.(x, y)
@deprecate complex(x::AbstractSparseVector{<:Real}, y::AbstractVector{<:Real}) complex.(x, y)

@deprecate diff(a::SparseMatrixCSC, dim::Integer) diff(a, dims=dim)

# END 0.7 deprecations

Expand Down
2 changes: 1 addition & 1 deletion stdlib/SparseArrays/src/linalg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -526,7 +526,7 @@ function sparse_diff2(a::SparseMatrixCSC{Tv,Ti}) where {Tv,Ti}
return SparseMatrixCSC(m, n-1, colptr, rowval, nzval)
end

diff(a::SparseMatrixCSC, dim::Integer)= dim==1 ? sparse_diff1(a) : sparse_diff2(a)
diff(a::SparseMatrixCSC; dims::Integer) = dims==1 ? sparse_diff1(a) : sparse_diff2(a)

## norm and rank
vecnorm(A::SparseMatrixCSC, p::Real=2) = vecnorm(view(A.nzval, 1:nnz(A)), p)
Expand Down
4 changes: 2 additions & 2 deletions stdlib/SparseArrays/test/sparse.jl
Original file line number Diff line number Diff line change
Expand Up @@ -569,13 +569,13 @@ end

@testset "issue #5853, sparse diff" begin
for i=1:2, a=Any[[1 2 3], reshape([1, 2, 3],(3,1)), Matrix(1.0I, 3, 3)]
@test all(diff(sparse(a),i) == diff(a,i))
@test all(diff(sparse(a),dims=i) == diff(a,dims=i))
end
end

@testset "access to undefined error types that initially allocate elements as #undef" begin
@test all(sparse(1:2, 1:2, Number[1,2])^2 == sparse(1:2, 1:2, [1,4]))
sd1 = diff(sparse([1,1,1], [1,2,3], Number[1,2,3]), 1)
sd1 = diff(sparse([1,1,1], [1,2,3], Number[1,2,3]), dims=1)
end

@testset "issue #6036" begin
Expand Down
32 changes: 16 additions & 16 deletions test/arrayops.jl
Original file line number Diff line number Diff line change
Expand Up @@ -618,23 +618,23 @@ Base.hash(::HashCollision, h::UInt) = h
let A, B, C, D
A = fill(1., 10, 10)
A[diagind(A)] = shuffle!([1:10;])
@test unique(A, 1) == A
@test unique(A, 2) == A
@test unique(A, dims=1) == A
@test unique(A, dims=2) == A

# 10 repeats of each row
B = A[shuffle!(repeat(1:10, 10)), :]
C = unique(B, 1)
C = unique(B, dims=1)
@test sortrows(C) == sortrows(A)
@test unique(B, 2) == B
@test unique(B', 2)' == C
@test unique(B, dims=2) == B
@test unique(B', dims=2)' == C

# Along third dimension
D = cat(3, B, B)
@test unique(D, 1) == cat(3, C, C)
@test unique(D, 3) == cat(3, B)
@test unique(D, dims=1) == cat(3, C, C)
@test unique(D, dims=3) == cat(3, B)

# With hash collisions
@test map(x -> x.x, unique(map(HashCollision, B), 1)) == C
@test map(x -> x.x, unique(map(HashCollision, B), dims=1)) == C
end

@testset "large matrices transpose" begin
Expand Down Expand Up @@ -2216,14 +2216,14 @@ end
X = [3 9 5;
7 4 2;
2 1 10]
@test diff(X,1) == [4 -5 -3; -5 -3 8]
@test diff(X,2) == [6 -4; -3 -2; -1 9]
@test diff(view(X, 1:2, 1:2),1) == [4 -5]
@test diff(view(X, 1:2, 1:2),2) == reshape([6; -3], (2,1))
@test diff(view(X, 2:3, 2:3),1) == [-3 8]
@test diff(view(X, 2:3, 2:3),2) == reshape([-2; 9], (2,1))
@test_throws ArgumentError diff(X,3)
@test_throws ArgumentError diff(X,-1)
@test diff(X,dims=1) == [4 -5 -3; -5 -3 8]
@test diff(X,dims=2) == [6 -4; -3 -2; -1 9]
@test diff(view(X, 1:2, 1:2),dims=1) == [4 -5]
@test diff(view(X, 1:2, 1:2),dims=2) == reshape([6; -3], (2,1))
@test diff(view(X, 2:3, 2:3),dims=1) == [-3 8]
@test diff(view(X, 2:3, 2:3),dims=2) == reshape([-2; 9], (2,1))
@test_throws ArgumentError diff(X,dims=3)
@test_throws ArgumentError diff(X,dims=-1)
end

@testset "accumulate, accumulate!" begin
Expand Down
4 changes: 2 additions & 2 deletions test/bitarray.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1412,8 +1412,8 @@ timesofar("cat")
@check_bit_operation diff(b1) Vector{Int}

b1 = bitrand(n1, n2)
@check_bit_operation diff(b1, 1) Matrix{Int}
@check_bit_operation diff(b1, 2) Matrix{Int}
@check_bit_operation diff(b1, dims=1) Matrix{Int}
@check_bit_operation diff(b1, dims=2) Matrix{Int}

b1 = bitrand(n1, n1)
@check_bit_operation svd(b1)
Expand Down
4 changes: 2 additions & 2 deletions test/offsetarray.jl
Original file line number Diff line number Diff line change
Expand Up @@ -406,8 +406,8 @@ seek(io, 0)
amin, amax = extrema(parent(A))
@test clamp.(A, (amax+amin)/2, amax).parent == clamp.(parent(A), (amax+amin)/2, amax)

@test unique(A, 1) == parent(A)
@test unique(A, 2) == parent(A)
@test unique(A, dims=1) == parent(A)
@test unique(A, dims=2) == parent(A)
v = OffsetArray(rand(8), (-2,))
@test sort(v) == OffsetArray(sort(parent(v)), v.offsets)
@test sortrows(A) == OffsetArray(sortrows(parent(A)), A.offsets)
Expand Down

0 comments on commit 1508942

Please sign in to comment.