StatProfilerHTML.jl report
Generated on Thu, 26 Mar 2020 19:09:01
File source code
Line Exclusive Inclusive Code
1 # This file is a part of Julia. License is MIT: https://julialang.org/license
2
3 ## Basic functions ##
4
5 """
6 AbstractArray{T,N}
7
8 Supertype for `N`-dimensional arrays (or array-like types) with elements of type `T`.
9 [`Array`](@ref) and other types are subtypes of this. See the manual section on the
10 [`AbstractArray` interface](@ref man-interface-array).
11 """
12 AbstractArray
13
14 convert(::Type{T}, a::T) where {T<:AbstractArray} = a
15 convert(::Type{AbstractArray{T}}, a::AbstractArray) where {T} = AbstractArray{T}(a)
16 convert(::Type{AbstractArray{T,N}}, a::AbstractArray{<:Any,N}) where {T,N} = AbstractArray{T,N}(a)
17
18 """
19 size(A::AbstractArray, [dim])
20
21 Return a tuple containing the dimensions of `A`. Optionally you can specify a
22 dimension to just get the length of that dimension.
23
24 Note that `size` may not be defined for arrays with non-standard indices, in which case [`axes`](@ref)
25 may be useful. See the manual chapter on [arrays with custom indices](@ref man-custom-indices).
26
27 # Examples
28 ```jldoctest
29 julia> A = fill(1, (2,3,4));
30
31 julia> size(A)
32 (2, 3, 4)
33
34 julia> size(A, 2)
35 3
36 ```
37 """
38 size(t::AbstractArray{T,N}, d) where {T,N} = d::Integer <= N ? size(t)[d] : 1
39
40 """
41 axes(A, d)
42
43 Return the valid range of indices for array `A` along dimension `d`.
44
45 See also [`size`](@ref), and the manual chapter on [arrays with custom indices](@ref man-custom-indices).
46
47 # Examples
48 ```jldoctest
49 julia> A = fill(1, (5,6,7));
50
51 julia> axes(A, 2)
52 Base.OneTo(6)
53 ```
54 """
55 function axes(A::AbstractArray{T,N}, d) where {T,N}
56 @_inline_meta
57 d::Integer <= N ? axes(A)[d] : OneTo(1)
58 end
59
60 """
61 axes(A)
62
63 Return the tuple of valid indices for array `A`.
64
65 # Examples
66 ```jldoctest
67 julia> A = fill(1, (5,6,7));
68
69 julia> axes(A)
70 (Base.OneTo(5), Base.OneTo(6), Base.OneTo(7))
71 ```
72 """
73 function axes(A)
74 @_inline_meta
75 1 (0 %)
1 (0 %) samples spent in axes
1 (100 %) (incl.) when called from axes1 line 95
1 (100 %) samples spent calling size
map(OneTo, size(A))
76 end
77
78 """
79 has_offset_axes(A)
80 has_offset_axes(A, B, ...)
81
82 Return `true` if the indices of `A` start with something other than 1 along any axis.
83 If multiple arguments are passed, equivalent to `has_offset_axes(A) | has_offset_axes(B) | ...`.
84 """
85 has_offset_axes(A) = _tuple_any(x->first(x)!=1, axes(A))
86 has_offset_axes(A...) = _tuple_any(has_offset_axes, A)
87 has_offset_axes(::Colon) = false
88
89 require_one_based_indexing(A...) = !has_offset_axes(A...) || throw(ArgumentError("offset arrays are not supported but got an array with index other than 1"))
90
91 # Performance optimization: get rid of a branch on `d` in `axes(A, d)`
92 # for d=1. 1d arrays are heavily used, and the first dimension comes up
93 # in other applications.
94 axes1(A::AbstractArray{<:Any,0}) = OneTo(1)
95 1 (0 %)
1 (0 %) samples spent in axes1
1 (100 %) (incl.) when called from eachindex line 267
1 (100 %) samples spent calling axes
axes1(A::AbstractArray) = (@_inline_meta; axes(A)[1])
96 axes1(iter) = OneTo(length(iter))
97
98 unsafe_indices(A) = axes(A)
99 unsafe_indices(r::AbstractRange) = (OneTo(unsafe_length(r)),) # Ranges use checked_sub for size
100
101 keys(a::AbstractArray) = CartesianIndices(axes(a))
102 keys(a::AbstractVector) = LinearIndices(a)
103
104 """
105 keytype(T::Type{<:AbstractArray})
106 keytype(A::AbstractArray)
107
108 Return the key type of an array. This is equal to the
109 `eltype` of the result of `keys(...)`, and is provided
110 mainly for compatibility with the dictionary interface.
111
112 # Examples
113 ```jldoctest
114 julia> keytype([1, 2, 3]) == Int
115 true
116
117 julia> keytype([1 2; 3 4])
118 CartesianIndex{2}
119 ```
120
121 !!! compat "Julia 1.2"
122 For arrays, this function requires at least Julia 1.2.
123 """
124 keytype(a::AbstractArray) = keytype(typeof(a))
125
126 keytype(A::Type{<:AbstractArray}) = CartesianIndex{ndims(A)}
127 keytype(A::Type{<:AbstractVector}) = Int
128
129 valtype(a::AbstractArray) = valtype(typeof(a))
130
131 """
132 valtype(T::Type{<:AbstractArray})
133 valtype(A::AbstractArray)
134
135 Return the value type of an array. This is identical to `eltype` and is
136 provided mainly for compatibility with the dictionary interface.
137
138 # Examples
139 ```jldoctest
140 julia> valtype(["one", "two", "three"])
141 String
142 ```
143
144 !!! compat "Julia 1.2"
145 For arrays, this function requires at least Julia 1.2.
146 """
147 valtype(A::Type{<:AbstractArray}) = eltype(A)
148
149 prevind(::AbstractArray, i::Integer) = Int(i)-1
150 nextind(::AbstractArray, i::Integer) = Int(i)+1
151
152 eltype(::Type{<:AbstractArray{E}}) where {E} = @isdefined(E) ? E : Any
153 elsize(A::AbstractArray) = elsize(typeof(A))
154
155 """
156 ndims(A::AbstractArray) -> Integer
157
158 Return the number of dimensions of `A`.
159
160 # Examples
161 ```jldoctest
162 julia> A = fill(1, (3,4,5));
163
164 julia> ndims(A)
165 3
166 ```
167 """
168 ndims(::AbstractArray{T,N}) where {T,N} = N
169 ndims(::Type{<:AbstractArray{T,N}}) where {T,N} = N
170
171 """
172 length(collection) -> Integer
173
174 Return the number of elements in the collection.
175
176 Use [`lastindex`](@ref) to get the last valid index of an indexable collection.
177
178 # Examples
179 ```jldoctest
180 julia> length(1:5)
181 5
182
183 julia> length([1, 2, 3, 4])
184 4
185
186 julia> length([1 2; 3 4])
187 4
188 ```
189 """
190 length
191
192 """
193 length(A::AbstractArray)
194
195 Return the number of elements in the array, defaults to `prod(size(A))`.
196
197 # Examples
198 ```jldoctest
199 julia> length([1, 2, 3, 4])
200 4
201
202 julia> length([1 2; 3 4])
203 4
204 ```
205 """
206 length(t::AbstractArray) = (@_inline_meta; prod(size(t)))
207
208 # `eachindex` is mostly an optimization of `keys`
209 eachindex(itrs...) = keys(itrs...)
210
211 # eachindex iterates over all indices. IndexCartesian definitions are later.
212 eachindex(A::AbstractVector) = (@_inline_meta(); axes1(A))
213
214 @noinline function throw_eachindex_mismatch(::IndexLinear, A...)
215 throw(DimensionMismatch("all inputs to eachindex must have the same indices, got $(join(eachindex.(A), ", ", " and "))"))
216 end
217 @noinline function throw_eachindex_mismatch(::IndexCartesian, A...)
218 throw(DimensionMismatch("all inputs to eachindex must have the same axes, got $(join(axes.(A), ", ", " and "))"))
219 end
220
221 """
222 eachindex(A...)
223
224 Create an iterable object for visiting each index of an `AbstractArray` `A` in an efficient
225 manner. For array types that have opted into fast linear indexing (like `Array`), this is
226 simply the range `1:length(A)`. For other array types, return a specialized Cartesian
227 range to efficiently index into the array with indices specified for every dimension. For
228 other iterables, including strings and dictionaries, return an iterator object
229 supporting arbitrary index types (e.g. unevenly spaced or non-integer indices).
230
231 If you supply more than one `AbstractArray` argument, `eachindex` will create an
232 iterable object that is fast for all arguments (a [`UnitRange`](@ref)
233 if all inputs have fast linear indexing, a [`CartesianIndices`](@ref)
234 otherwise).
235 If the arrays have different sizes and/or dimensionalities, a DimensionMismatch exception
236 will be thrown.
237 # Examples
238 ```jldoctest
239 julia> A = [1 2; 3 4];
240
241 julia> for i in eachindex(A) # linear indexing
242 println(i)
243 end
244 1
245 2
246 3
247 4
248
249 julia> for i in eachindex(view(A, 1:2, 1:1)) # Cartesian indexing
250 println(i)
251 end
252 CartesianIndex(1, 1)
253 CartesianIndex(2, 1)
254 ```
255 """
256 eachindex(A::AbstractArray) = (@_inline_meta(); eachindex(IndexStyle(A), A))
257
258 function eachindex(A::AbstractArray, B::AbstractArray)
259 @_inline_meta
260 eachindex(IndexStyle(A,B), A, B)
261 end
262 function eachindex(A::AbstractArray, B::AbstractArray...)
263 @_inline_meta
264 eachindex(IndexStyle(A,B...), A, B...)
265 end
266 eachindex(::IndexLinear, A::AbstractArray) = (@_inline_meta; OneTo(length(A)))
267 1 (0 %)
1 (0 %) samples spent in eachindex
1 (100 %) (incl.) when called from lastindex line 302
1 (100 %) samples spent calling axes1
eachindex(::IndexLinear, A::AbstractVector) = (@_inline_meta; axes1(A))
268 function eachindex(::IndexLinear, A::AbstractArray, B::AbstractArray...)
269 @_inline_meta
270 indsA = eachindex(IndexLinear(), A)
271 _all_match_first(X->eachindex(IndexLinear(), X), indsA, B...) ||
272 throw_eachindex_mismatch(IndexLinear(), A, B...)
273 indsA
274 end
275 function _all_match_first(f::F, inds, A, B...) where F<:Function
276 @_inline_meta
277 (inds == f(A)) & _all_match_first(f, inds, B...)
278 end
279 _all_match_first(f::F, inds) where F<:Function = true
280
281 # keys with an IndexStyle
282 keys(s::IndexStyle, A::AbstractArray, B::AbstractArray...) = eachindex(s, A, B...)
283
284 """
285 lastindex(collection) -> Integer
286 lastindex(collection, d) -> Integer
287
288 Return the last index of `collection`. If `d` is given, return the last index of `collection` along dimension `d`.
289
290 The syntaxes `A[end]` and `A[end, end]` lower to `A[lastindex(A)]` and
291 `A[lastindex(A, 1), lastindex(A, 2)]`, respectively.
292
293 # Examples
294 ```jldoctest
295 julia> lastindex([1,2,4])
296 3
297
298 julia> lastindex(rand(3,4,5), 2)
299 4
300 ```
301 """
302 1 (0 %)
1 (0 %) samples spent in lastindex
1 (100 %) (incl.) when called from push! line 913
1 (100 %) samples spent calling eachindex
lastindex(a::AbstractArray) = (@_inline_meta; last(eachindex(IndexLinear(), a)))
303 lastindex(a::AbstractArray, d) = (@_inline_meta; last(axes(a, d)))
304
305 """
306 firstindex(collection) -> Integer
307 firstindex(collection, d) -> Integer
308
309 Return the first index of `collection`. If `d` is given, return the first index of `collection` along dimension `d`.
310
311 # Examples
312 ```jldoctest
313 julia> firstindex([1,2,4])
314 1
315
316 julia> firstindex(rand(3,4,5), 2)
317 1
318 ```
319 """
320 firstindex(a::AbstractArray) = (@_inline_meta; first(eachindex(IndexLinear(), a)))
321 firstindex(a::AbstractArray, d) = (@_inline_meta; first(axes(a, d)))
322
323 first(a::AbstractArray) = a[first(eachindex(a))]
324
325 """
326 first(coll)
327
328 Get the first element of an iterable collection. Return the start point of an
329 [`AbstractRange`](@ref) even if it is empty.
330
331 # Examples
332 ```jldoctest
333 julia> first(2:2:10)
334 2
335
336 julia> first([1; 2; 3; 4])
337 1
338 ```
339 """
340 function first(itr)
341 x = iterate(itr)
342 x === nothing && throw(ArgumentError("collection must be non-empty"))
343 x[1]
344 end
345
346 """
347 last(coll)
348
349 Get the last element of an ordered collection, if it can be computed in O(1) time. This is
350 accomplished by calling [`lastindex`](@ref) to get the last index. Return the end
351 point of an [`AbstractRange`](@ref) even if it is empty.
352
353 # Examples
354 ```jldoctest
355 julia> last(1:2:10)
356 9
357
358 julia> last([1; 2; 3; 4])
359 4
360 ```
361 """
362 last(a) = a[end]
363
364 """
365 strides(A)
366
367 Return a tuple of the memory strides in each dimension.
368
369 # Examples
370 ```jldoctest
371 julia> A = fill(1, (3,4,5));
372
373 julia> strides(A)
374 (1, 3, 12)
375 ```
376 """
377 function strides end
378
379 """
380 stride(A, k::Integer)
381
382 Return the distance in memory (in number of elements) between adjacent elements in dimension `k`.
383
384 # Examples
385 ```jldoctest
386 julia> A = fill(1, (3,4,5));
387
388 julia> stride(A,2)
389 3
390
391 julia> stride(A,3)
392 12
393 ```
394 """
395 stride(A::AbstractArray, k::Integer) = strides(A)[k]
396
397 @inline size_to_strides(s, d, sz...) = (s, size_to_strides(s * d, sz...)...)
398 size_to_strides(s, d) = (s,)
399 size_to_strides(s) = ()
400
401
402 function isassigned(a::AbstractArray, i::Integer...)
403 try
404 a[i...]
405 true
406 catch e
407 if isa(e, BoundsError) || isa(e, UndefRefError)
408 return false
409 else
410 rethrow()
411 end
412 end
413 end
414
415 # used to compute "end" for last index
416 function trailingsize(A, n)
417 s = 1
418 for i=n:ndims(A)
419 s *= size(A,i)
420 end
421 return s
422 end
423 function trailingsize(inds::Indices, n)
424 s = 1
425 for i=n:length(inds)
426 s *= unsafe_length(inds[i])
427 end
428 return s
429 end
430 # This version is type-stable even if inds is heterogeneous
431 function trailingsize(inds::Indices)
432 @_inline_meta
433 prod(map(unsafe_length, inds))
434 end
435
436 ## Bounds checking ##
437
438 # The overall hierarchy is
439 # `checkbounds(A, I...)` ->
440 # `checkbounds(Bool, A, I...)` ->
441 # `checkbounds_indices(Bool, IA, I)`, which recursively calls
442 # `checkindex` for each dimension
443 #
444 # See the "boundscheck" devdocs for more information.
445 #
446 # Note this hierarchy has been designed to reduce the likelihood of
447 # method ambiguities. We try to make `checkbounds` the place to
448 # specialize on array type, and try to avoid specializations on index
449 # types; conversely, `checkindex` is intended to be specialized only
450 # on index type (especially, its last argument).
451
452 """
453 checkbounds(Bool, A, I...)
454
455 Return `true` if the specified indices `I` are in bounds for the given
456 array `A`. Subtypes of `AbstractArray` should specialize this method
457 if they need to provide custom bounds checking behaviors; however, in
458 many cases one can rely on `A`'s indices and [`checkindex`](@ref).
459
460 See also [`checkindex`](@ref).
461
462 # Examples
463 ```jldoctest
464 julia> A = rand(3, 3);
465
466 julia> checkbounds(Bool, A, 2)
467 true
468
469 julia> checkbounds(Bool, A, 3, 4)
470 false
471
472 julia> checkbounds(Bool, A, 1:3)
473 true
474
475 julia> checkbounds(Bool, A, 1:3, 2:4)
476 false
477 ```
478 """
479 function checkbounds(::Type{Bool}, A::AbstractArray, I...)
480 @_inline_meta
481 checkbounds_indices(Bool, axes(A), I)
482 end
483
484 # Linear indexing is explicitly allowed when there is only one (non-cartesian) index
485 function checkbounds(::Type{Bool}, A::AbstractArray, i)
486 @_inline_meta
487 checkindex(Bool, eachindex(IndexLinear(), A), i)
488 end
489 # As a special extension, allow using logical arrays that match the source array exactly
490 function checkbounds(::Type{Bool}, A::AbstractArray{<:Any,N}, I::AbstractArray{Bool,N}) where N
491 @_inline_meta
492 axes(A) == axes(I)
493 end
494
495 """
496 checkbounds(A, I...)
497
498 Throw an error if the specified indices `I` are not in bounds for the given array `A`.
499 """
500 function checkbounds(A::AbstractArray, I...)
501 @_inline_meta
502 checkbounds(Bool, A, I...) || throw_boundserror(A, I)
503 nothing
504 end
505
506 """
507 checkbounds_indices(Bool, IA, I)
508
509 Return `true` if the "requested" indices in the tuple `I` fall within
510 the bounds of the "permitted" indices specified by the tuple
511 `IA`. This function recursively consumes elements of these tuples,
512 usually in a 1-for-1 fashion,
513
514 checkbounds_indices(Bool, (IA1, IA...), (I1, I...)) = checkindex(Bool, IA1, I1) &
515 checkbounds_indices(Bool, IA, I)
516
517 Note that [`checkindex`](@ref) is being used to perform the actual
518 bounds-check for a single dimension of the array.
519
520 There are two important exceptions to the 1-1 rule: linear indexing and
521 CartesianIndex{N}, both of which may "consume" more than one element
522 of `IA`.
523
524 See also [`checkbounds`](@ref).
525 """
526 function checkbounds_indices(::Type{Bool}, IA::Tuple, I::Tuple)
527 @_inline_meta
528 checkindex(Bool, IA[1], I[1]) & checkbounds_indices(Bool, tail(IA), tail(I))
529 end
530 function checkbounds_indices(::Type{Bool}, ::Tuple{}, I::Tuple)
531 @_inline_meta
532 checkindex(Bool, OneTo(1), I[1]) & checkbounds_indices(Bool, (), tail(I))
533 end
534 checkbounds_indices(::Type{Bool}, IA::Tuple, ::Tuple{}) = (@_inline_meta; all(x->unsafe_length(x)==1, IA))
535 checkbounds_indices(::Type{Bool}, ::Tuple{}, ::Tuple{}) = true
536
537 throw_boundserror(A, I) = (@_noinline_meta; throw(BoundsError(A, I)))
538
539 # check along a single dimension
540 """
541 checkindex(Bool, inds::AbstractUnitRange, index)
542
543 Return `true` if the given `index` is within the bounds of
544 `inds`. Custom types that would like to behave as indices for all
545 arrays can extend this method in order to provide a specialized bounds
546 checking implementation.
547
548 # Examples
549 ```jldoctest
550 julia> checkindex(Bool, 1:20, 8)
551 true
552
553 julia> checkindex(Bool, 1:20, 21)
554 false
555 ```
556 """
557 checkindex(::Type{Bool}, inds::AbstractUnitRange, i) =
558 throw(ArgumentError("unable to check bounds for indices of type $(typeof(i))"))
559 checkindex(::Type{Bool}, inds::AbstractUnitRange, i::Real) = (first(inds) <= i) & (i <= last(inds))
560 checkindex(::Type{Bool}, inds::AbstractUnitRange, ::Colon) = true
561 checkindex(::Type{Bool}, inds::AbstractUnitRange, ::Slice) = true
562 function checkindex(::Type{Bool}, inds::AbstractUnitRange, r::AbstractRange)
563 @_propagate_inbounds_meta
564 isempty(r) | (checkindex(Bool, inds, first(r)) & checkindex(Bool, inds, last(r)))
565 end
566 checkindex(::Type{Bool}, indx::AbstractUnitRange, I::AbstractVector{Bool}) = indx == axes1(I)
567 checkindex(::Type{Bool}, indx::AbstractUnitRange, I::AbstractArray{Bool}) = false
568 function checkindex(::Type{Bool}, inds::AbstractUnitRange, I::AbstractArray)
569 @_inline_meta
570 b = true
571 for i in I
572 b &= checkindex(Bool, inds, i)
573 end
574 b
575 end
576
577 # See also specializations in multidimensional
578
579 ## Constructors ##
580
581 # default arguments to similar()
582 """
583 similar(array, [element_type=eltype(array)], [dims=size(array)])
584
585 Create an uninitialized mutable array with the given element type and size, based upon the
586 given source array. The second and third arguments are both optional, defaulting to the
587 given array's `eltype` and `size`. The dimensions may be specified either as a single tuple
588 argument or as a series of integer arguments.
589
590 Custom AbstractArray subtypes may choose which specific array type is best-suited to return
591 for the given element type and dimensionality. If they do not specialize this method, the
592 default is an `Array{element_type}(undef, dims...)`.
593
594 For example, `similar(1:10, 1, 4)` returns an uninitialized `Array{Int,2}` since ranges are
595 neither mutable nor support 2 dimensions:
596
597 ```julia-repl
598 julia> similar(1:10, 1, 4)
599 1×4 Array{Int64,2}:
600 4419743872 4374413872 4419743888 0
601 ```
602
603 Conversely, `similar(trues(10,10), 2)` returns an uninitialized `BitVector` with two
604 elements since `BitArray`s are both mutable and can support 1-dimensional arrays:
605
606 ```julia-repl
607 julia> similar(trues(10,10), 2)
608 2-element BitArray{1}:
609 0
610 0
611 ```
612
613 Since `BitArray`s can only store elements of type [`Bool`](@ref), however, if you request a
614 different element type it will create a regular `Array` instead:
615
616 ```julia-repl
617 julia> similar(falses(10), Float64, 2, 4)
618 2×4 Array{Float64,2}:
619 2.18425e-314 2.18425e-314 2.18425e-314 2.18425e-314
620 2.18425e-314 2.18425e-314 2.18425e-314 2.18425e-314
621 ```
622
623 """
624 similar(a::AbstractArray{T}) where {T} = similar(a, T)
625 similar(a::AbstractArray, ::Type{T}) where {T} = similar(a, T, to_shape(axes(a)))
626 similar(a::AbstractArray{T}, dims::Tuple) where {T} = similar(a, T, to_shape(dims))
627 similar(a::AbstractArray{T}, dims::DimOrInd...) where {T} = similar(a, T, to_shape(dims))
628 similar(a::AbstractArray, ::Type{T}, dims::DimOrInd...) where {T} = similar(a, T, to_shape(dims))
629 # Similar supports specifying dims as either Integers or AbstractUnitRanges or any mixed combination
630 # thereof. Ideally, we'd just convert Integers to OneTos and then call a canonical method with the axes,
631 # but we don't want to require all AbstractArray subtypes to dispatch on Base.OneTo. So instead we
632 # define this method to convert supported axes to Ints, with the expectation that an offset array
633 # package will define a method with dims::Tuple{Union{Integer, UnitRange}, Vararg{Union{Integer, UnitRange}}}
634 similar(a::AbstractArray, ::Type{T}, dims::Tuple{Union{Integer, OneTo}, Vararg{Union{Integer, OneTo}}}) where {T} = similar(a, T, to_shape(dims))
635 # similar creates an Array by default
636 similar(a::AbstractArray, ::Type{T}, dims::Dims{N}) where {T,N} = Array{T,N}(undef, dims)
637
638 to_shape(::Tuple{}) = ()
639 to_shape(dims::Dims) = dims
640 to_shape(dims::DimsOrInds) = map(to_shape, dims)::DimsOrInds
641 # each dimension
642 to_shape(i::Int) = i
643 to_shape(i::Integer) = Int(i)
644 to_shape(r::OneTo) = Int(last(r))
645 to_shape(r::AbstractUnitRange) = r
646
647 """
648 similar(storagetype, axes)
649
650 Create an uninitialized mutable array analogous to that specified by
651 `storagetype`, but with `axes` specified by the last
652 argument. `storagetype` might be a type or a function.
653
654 **Examples**:
655
656 similar(Array{Int}, axes(A))
657
658 creates an array that "acts like" an `Array{Int}` (and might indeed be
659 backed by one), but which is indexed identically to `A`. If `A` has
660 conventional indexing, this will be identical to
661 `Array{Int}(undef, size(A))`, but if `A` has unconventional indexing then the
662 indices of the result will match `A`.
663
664 similar(BitArray, (axes(A, 2),))
665
666 would create a 1-dimensional logical array whose indices match those
667 of the columns of `A`.
668 """
669 similar(::Type{T}, dims::DimOrInd...) where {T<:AbstractArray} = similar(T, dims)
670 similar(::Type{T}, shape::Tuple{Union{Integer, OneTo}, Vararg{Union{Integer, OneTo}}}) where {T<:AbstractArray} = similar(T, to_shape(shape))
671 similar(::Type{T}, dims::Dims) where {T<:AbstractArray} = T(undef, dims)
672
673 """
674 empty(v::AbstractVector, [eltype])
675
676 Create an empty vector similar to `v`, optionally changing the `eltype`.
677
678 # Examples
679
680 ```jldoctest
681 julia> empty([1.0, 2.0, 3.0])
682 0-element Array{Float64,1}
683
684 julia> empty([1.0, 2.0, 3.0], String)
685 0-element Array{String,1}
686 ```
687 """
688 empty(a::AbstractVector{T}, ::Type{U}=T) where {T,U} = Vector{U}()
689
690 # like empty, but should return a mutable collection, a Vector by default
691 emptymutable(a::AbstractVector{T}, ::Type{U}=T) where {T,U} = Vector{U}()
692 emptymutable(itr, ::Type{U}) where {U} = Vector{U}()
693
694 """
695 copy!(dst, src) -> dst
696
697 In-place [`copy`](@ref) of `src` into `dst`, discarding any pre-existing
698 elements in `dst`.
699 If `dst` and `src` are of the same type, `dst == src` should hold after
700 the call. If `dst` and `src` are multidimensional arrays, they must have
701 equal [`axes`](@ref).
702 See also [`copyto!`](@ref).
703
704 !!! compat "Julia 1.1"
705 This method requires at least Julia 1.1. In Julia 1.0 this method
706 is available from the `Future` standard library as `Future.copy!`.
707 """
708 copy!(dst::AbstractVector, src::AbstractVector) = append!(empty!(dst), src)
709
710 function copy!(dst::AbstractArray, src::AbstractArray)
711 axes(dst) == axes(src) || throw(ArgumentError(
712 "arrays must have the same axes for copy! (consider using `copyto!`)"))
713 copyto!(dst, src)
714 end
715
716 ## from general iterable to any array
717
718 function copyto!(dest::AbstractArray, src)
719 destiter = eachindex(dest)
720 y = iterate(destiter)
721 for x in src
722 y === nothing &&
723 throw(ArgumentError("destination has fewer elements than required"))
724 dest[y[1]] = x
725 y = iterate(destiter, y[2])
726 end
727 return dest
728 end
729
730 function copyto!(dest::AbstractArray, dstart::Integer, src)
731 i = Int(dstart)
732 for x in src
733 dest[i] = x
734 i += 1
735 end
736 return dest
737 end
738
739 # copy from an some iterable object into an AbstractArray
740 function copyto!(dest::AbstractArray, dstart::Integer, src, sstart::Integer)
741 if (sstart < 1)
742 throw(ArgumentError(string("source start offset (",sstart,") is < 1")))
743 end
744 y = iterate(src)
745 for j = 1:(sstart-1)
746 if y === nothing
747 throw(ArgumentError(string("source has fewer elements than required, ",
748 "expected at least ",sstart,", got ",j-1)))
749 end
750 y = iterate(src, y[2])
751 end
752 if y === nothing
753 throw(ArgumentError(string("source has fewer elements than required, ",
754 "expected at least ",sstart,", got ",sstart-1)))
755 end
756 i = Int(dstart)
757 while y !== nothing
758 val, st = y
759 dest[i] = val
760 i += 1
761 y = iterate(src, st)
762 end
763 return dest
764 end
765
766 # this method must be separate from the above since src might not have a length
767 function copyto!(dest::AbstractArray, dstart::Integer, src, sstart::Integer, n::Integer)
768 n < 0 && throw(ArgumentError(string("tried to copy n=", n, " elements, but n should be nonnegative")))
769 n == 0 && return dest
770 dmax = dstart + n - 1
771 inds = LinearIndices(dest)
772 if (dstart ∉ inds || dmax ∉ inds) | (sstart < 1)
773 sstart < 1 && throw(ArgumentError(string("source start offset (",sstart,") is < 1")))
774 throw(BoundsError(dest, dstart:dmax))
775 end
776 y = iterate(src)
777 for j = 1:(sstart-1)
778 if y === nothing
779 throw(ArgumentError(string("source has fewer elements than required, ",
780 "expected at least ",sstart,", got ",j-1)))
781 end
782 y = iterate(src, y[2])
783 end
784 i = Int(dstart)
785 while i <= dmax && y !== nothing
786 val, st = y
787 @inbounds dest[i] = val
788 y = iterate(src, st)
789 i += 1
790 end
791 i <= dmax && throw(BoundsError(dest, i))
792 return dest
793 end
794
795 ## copy between abstract arrays - generally more efficient
796 ## since a single index variable can be used.
797
798 copyto!(dest::AbstractArray, src::AbstractArray) =
799 copyto!(IndexStyle(dest), dest, IndexStyle(src), src)
800
801 function copyto!(::IndexStyle, dest::AbstractArray, ::IndexStyle, src::AbstractArray)
802 destinds, srcinds = LinearIndices(dest), LinearIndices(src)
803 isempty(srcinds) || (checkbounds(Bool, destinds, first(srcinds)) && checkbounds(Bool, destinds, last(srcinds))) ||
804 throw(BoundsError(dest, srcinds))
805 @inbounds for i in srcinds
806 dest[i] = src[i]
807 end
808 return dest
809 end
810
811 function copyto!(::IndexStyle, dest::AbstractArray, ::IndexCartesian, src::AbstractArray)
812 destinds, srcinds = LinearIndices(dest), LinearIndices(src)
813 isempty(srcinds) || (checkbounds(Bool, destinds, first(srcinds)) && checkbounds(Bool, destinds, last(srcinds))) ||
814 throw(BoundsError(dest, srcinds))
815 i = 0
816 @inbounds for a in src
817 dest[i+=1] = a
818 end
819 return dest
820 end
821
822 function copyto!(dest::AbstractArray, dstart::Integer, src::AbstractArray)
823 copyto!(dest, dstart, src, first(LinearIndices(src)), length(src))
824 end
825
826 function copyto!(dest::AbstractArray, dstart::Integer, src::AbstractArray, sstart::Integer)
827 srcinds = LinearIndices(src)
828 checkbounds(Bool, srcinds, sstart) || throw(BoundsError(src, sstart))
829 copyto!(dest, dstart, src, sstart, last(srcinds)-sstart+1)
830 end
831
832 function copyto!(dest::AbstractArray, dstart::Integer,
833 src::AbstractArray, sstart::Integer,
834 n::Integer)
835 n == 0 && return dest
836 n < 0 && throw(ArgumentError(string("tried to copy n=", n, " elements, but n should be nonnegative")))
837 destinds, srcinds = LinearIndices(dest), LinearIndices(src)
838 (checkbounds(Bool, destinds, dstart) && checkbounds(Bool, destinds, dstart+n-1)) || throw(BoundsError(dest, dstart:dstart+n-1))
839 (checkbounds(Bool, srcinds, sstart) && checkbounds(Bool, srcinds, sstart+n-1)) || throw(BoundsError(src, sstart:sstart+n-1))
840 @inbounds for i = 0:(n-1)
841 dest[dstart+i] = src[sstart+i]
842 end
843 return dest
844 end
845
846 function copy(a::AbstractArray)
847 @_propagate_inbounds_meta
848 copymutable(a)
849 end
850
851 function copyto!(B::AbstractVecOrMat{R}, ir_dest::AbstractRange{Int}, jr_dest::AbstractRange{Int},
852 A::AbstractVecOrMat{S}, ir_src::AbstractRange{Int}, jr_src::AbstractRange{Int}) where {R,S}
853 if length(ir_dest) != length(ir_src)
854 throw(ArgumentError(string("source and destination must have same size (got ",
855 length(ir_src)," and ",length(ir_dest),")")))
856 end
857 if length(jr_dest) != length(jr_src)
858 throw(ArgumentError(string("source and destination must have same size (got ",
859 length(jr_src)," and ",length(jr_dest),")")))
860 end
861 @boundscheck checkbounds(B, ir_dest, jr_dest)
862 @boundscheck checkbounds(A, ir_src, jr_src)
863 jdest = first(jr_dest)
864 for jsrc in jr_src
865 idest = first(ir_dest)
866 for isrc in ir_src
867 @inbounds B[idest,jdest] = A[isrc,jsrc]
868 idest += step(ir_dest)
869 end
870 jdest += step(jr_dest)
871 end
872 return B
873 end
874
875
876 """
877 copymutable(a)
878
879 Make a mutable copy of an array or iterable `a`. For `a::Array`,
880 this is equivalent to `copy(a)`, but for other array types it may
881 differ depending on the type of `similar(a)`. For generic iterables
882 this is equivalent to `collect(a)`.
883
884 # Examples
885 ```jldoctest
886 julia> tup = (1, 2, 3)
887 (1, 2, 3)
888
889 julia> Base.copymutable(tup)
890 3-element Array{Int64,1}:
891 1
892 2
893 3
894 ```
895 """
896 function copymutable(a::AbstractArray)
897 @_propagate_inbounds_meta
898 copyto!(similar(a), a)
899 end
900 copymutable(itr) = collect(itr)
901
902 zero(x::AbstractArray{T}) where {T} = fill!(similar(x), zero(T))
903
904 ## iteration support for arrays by iterating over `eachindex` in the array ##
905 # Allows fast iteration by default for both IndexLinear and IndexCartesian arrays
906
907 # While the definitions for IndexLinear are all simple enough to inline on their
908 # own, IndexCartesian's CartesianIndices is more complicated and requires explicit
909 # inlining.
910 function iterate(A::AbstractArray, state=(eachindex(A),))
911 y = iterate(state...)
912 y === nothing && return nothing
913 A[y[1]], (state[1], tail(y)...)
914 end
915
916 isempty(a::AbstractArray) = (length(a) == 0)
917
918
919 ## range conversions ##
920
921 map(::Type{T}, r::StepRange) where {T<:Real} = T(r.start):T(r.step):T(last(r))
922 map(::Type{T}, r::UnitRange) where {T<:Real} = T(r.start):T(last(r))
923 map(::Type{T}, r::StepRangeLen) where {T<:AbstractFloat} = convert(StepRangeLen{T}, r)
924 function map(::Type{T}, r::LinRange) where T<:AbstractFloat
925 LinRange(T(r.start), T(r.stop), length(r))
926 end
927
928 ## unsafe/pointer conversions ##
929
930 # note: the following type definitions don't mean any AbstractArray is convertible to
931 # a data Ref. they just map the array element type to the pointer type for
932 # convenience in cases that work.
933 pointer(x::AbstractArray{T}) where {T} = unsafe_convert(Ptr{T}, x)
934 function pointer(x::AbstractArray{T}, i::Integer) where T
935 @_inline_meta
936 unsafe_convert(Ptr{T}, x) + (i - first(LinearIndices(x)))*elsize(x)
937 end
938
939 ## Approach:
940 # We only define one fallback method on getindex for all argument types.
941 # That dispatches to an (inlined) internal _getindex function, where the goal is
942 # to transform the indices such that we can call the only getindex method that
943 # we require the type A{T,N} <: AbstractArray{T,N} to define; either:
944 # getindex(::A, ::Int) # if IndexStyle(A) == IndexLinear() OR
945 # getindex(::A{T,N}, ::Vararg{Int, N}) where {T,N} # if IndexCartesian()
946 # If the subtype hasn't defined the required method, it falls back to the
947 # _getindex function again where an error is thrown to prevent stack overflows.
948 """
949 getindex(A, inds...)
950
951 Return a subset of array `A` as specified by `inds`, where each `ind` may be an
952 `Int`, an [`AbstractRange`](@ref), or a [`Vector`](@ref). See the manual section on
953 [array indexing](@ref man-array-indexing) for details.
954
955 # Examples
956 ```jldoctest
957 julia> A = [1 2; 3 4]
958 2×2 Array{Int64,2}:
959 1 2
960 3 4
961
962 julia> getindex(A, 1)
963 1
964
965 julia> getindex(A, [2, 1])
966 2-element Array{Int64,1}:
967 3
968 1
969
970 julia> getindex(A, 2:4)
971 3-element Array{Int64,1}:
972 3
973 2
974 4
975 ```
976 """
977 function getindex(A::AbstractArray, I...)
978 @_propagate_inbounds_meta
979 error_if_canonical_getindex(IndexStyle(A), A, I...)
980 _getindex(IndexStyle(A), A, to_indices(A, I)...)
981 end
982 function unsafe_getindex(A::AbstractArray, I...)
983 @_inline_meta
984 @inbounds r = getindex(A, I...)
985 r
986 end
987
988 error_if_canonical_getindex(::IndexLinear, A::AbstractArray, ::Int) =
989 error("getindex not defined for ", typeof(A))
990 error_if_canonical_getindex(::IndexCartesian, A::AbstractArray{T,N}, ::Vararg{Int,N}) where {T,N} =
991 error("getindex not defined for ", typeof(A))
992 error_if_canonical_getindex(::IndexStyle, ::AbstractArray, ::Any...) = nothing
993
994 ## Internal definitions
995 _getindex(::IndexStyle, A::AbstractArray, I...) =
996 error("getindex for $(typeof(A)) with types $(typeof(I)) is not supported")
997
998 ## IndexLinear Scalar indexing: canonical method is one Int
999 _getindex(::IndexLinear, A::AbstractArray, i::Int) = (@_propagate_inbounds_meta; getindex(A, i))
1000 function _getindex(::IndexLinear, A::AbstractArray, I::Vararg{Int,M}) where M
1001 @_inline_meta
1002 @boundscheck checkbounds(A, I...) # generally _to_linear_index requires bounds checking
1003 @inbounds r = getindex(A, _to_linear_index(A, I...))
1004 r
1005 end
1006 _to_linear_index(A::AbstractArray, i::Int) = i
1007 _to_linear_index(A::AbstractVector, i::Int, I::Int...) = i
1008 _to_linear_index(A::AbstractArray) = 1
1009 _to_linear_index(A::AbstractArray, I::Int...) = (@_inline_meta; _sub2ind(A, I...))
1010
1011 ## IndexCartesian Scalar indexing: Canonical method is full dimensionality of Ints
1012 function _getindex(::IndexCartesian, A::AbstractArray, I::Vararg{Int,M}) where M
1013 @_inline_meta
1014 @boundscheck checkbounds(A, I...) # generally _to_subscript_indices requires bounds checking
1015 @inbounds r = getindex(A, _to_subscript_indices(A, I...)...)
1016 r
1017 end
1018 function _getindex(::IndexCartesian, A::AbstractArray{T,N}, I::Vararg{Int, N}) where {T,N}
1019 @_propagate_inbounds_meta
1020 getindex(A, I...)
1021 end
1022 _to_subscript_indices(A::AbstractArray, i::Int) = (@_inline_meta; _unsafe_ind2sub(A, i))
1023 _to_subscript_indices(A::AbstractArray{T,N}) where {T,N} = (@_inline_meta; fill_to_length((), 1, Val(N)))
1024 _to_subscript_indices(A::AbstractArray{T,0}) where {T} = ()
1025 _to_subscript_indices(A::AbstractArray{T,0}, i::Int) where {T} = ()
1026 _to_subscript_indices(A::AbstractArray{T,0}, I::Int...) where {T} = ()
1027 function _to_subscript_indices(A::AbstractArray{T,N}, I::Int...) where {T,N}
1028 @_inline_meta
1029 J, Jrem = IteratorsMD.split(I, Val(N))
1030 _to_subscript_indices(A, J, Jrem)
1031 end
1032 _to_subscript_indices(A::AbstractArray, J::Tuple, Jrem::Tuple{}) =
1033 __to_subscript_indices(A, axes(A), J, Jrem)
1034 function __to_subscript_indices(A::AbstractArray,
1035 ::Tuple{AbstractUnitRange,Vararg{AbstractUnitRange}}, J::Tuple, Jrem::Tuple{})
1036 @_inline_meta
1037 (J..., map(first, tail(_remaining_size(J, axes(A))))...)
1038 end
1039 _to_subscript_indices(A, J::Tuple, Jrem::Tuple) = J # already bounds-checked, safe to drop
1040 _to_subscript_indices(A::AbstractArray{T,N}, I::Vararg{Int,N}) where {T,N} = I
1041 _remaining_size(::Tuple{Any}, t::Tuple) = t
1042 _remaining_size(h::Tuple, t::Tuple) = (@_inline_meta; _remaining_size(tail(h), tail(t)))
1043 _unsafe_ind2sub(::Tuple{}, i) = () # _ind2sub may throw(BoundsError()) in this case
1044 _unsafe_ind2sub(sz, i) = (@_inline_meta; _ind2sub(sz, i))
1045
1046 ## Setindex! is defined similarly. We first dispatch to an internal _setindex!
1047 # function that allows dispatch on array storage
1048
1049 """
1050 setindex!(A, X, inds...)
1051 A[inds...] = X
1052
1053 Store values from array `X` within some subset of `A` as specified by `inds`.
1054 The syntax `A[inds...] = X` is equivalent to `setindex!(A, X, inds...)`.
1055
1056 # Examples
1057 ```jldoctest
1058 julia> A = zeros(2,2);
1059
1060 julia> setindex!(A, [10, 20], [1, 2]);
1061
1062 julia> A[[3, 4]] = [30, 40];
1063
1064 julia> A
1065 2×2 Array{Float64,2}:
1066 10.0 30.0
1067 20.0 40.0
1068 ```
1069 """
1070 function setindex!(A::AbstractArray, v, I...)
1071 @_propagate_inbounds_meta
1072 error_if_canonical_setindex(IndexStyle(A), A, I...)
1073 _setindex!(IndexStyle(A), A, v, to_indices(A, I)...)
1074 end
1075 function unsafe_setindex!(A::AbstractArray, v, I...)
1076 @_inline_meta
1077 @inbounds r = setindex!(A, v, I...)
1078 r
1079 end
1080
1081 error_if_canonical_setindex(::IndexLinear, A::AbstractArray, ::Int) =
1082 error("setindex! not defined for ", typeof(A))
1083 error_if_canonical_setindex(::IndexCartesian, A::AbstractArray{T,N}, ::Vararg{Int,N}) where {T,N} =
1084 error("setindex! not defined for ", typeof(A))
1085 error_if_canonical_setindex(::IndexStyle, ::AbstractArray, ::Any...) = nothing
1086
1087 ## Internal definitions
1088 _setindex!(::IndexStyle, A::AbstractArray, v, I...) =
1089 error("setindex! for $(typeof(A)) with types $(typeof(I)) is not supported")
1090
1091 ## IndexLinear Scalar indexing
1092 _setindex!(::IndexLinear, A::AbstractArray, v, i::Int) = (@_propagate_inbounds_meta; setindex!(A, v, i))
1093 function _setindex!(::IndexLinear, A::AbstractArray, v, I::Vararg{Int,M}) where M
1094 @_inline_meta
1095 @boundscheck checkbounds(A, I...)
1096 @inbounds r = setindex!(A, v, _to_linear_index(A, I...))
1097 r
1098 end
1099
1100 # IndexCartesian Scalar indexing
1101 function _setindex!(::IndexCartesian, A::AbstractArray{T,N}, v, I::Vararg{Int, N}) where {T,N}
1102 @_propagate_inbounds_meta
1103 setindex!(A, v, I...)
1104 end
1105 function _setindex!(::IndexCartesian, A::AbstractArray, v, I::Vararg{Int,M}) where M
1106 @_inline_meta
1107 @boundscheck checkbounds(A, I...)
1108 @inbounds r = setindex!(A, v, _to_subscript_indices(A, I...)...)
1109 r
1110 end
1111
1112 """
1113 parent(A)
1114
1115 Returns the "parent array" of an array view type (e.g., `SubArray`), or the array itself if
1116 it is not a view.
1117
1118 # Examples
1119 ```jldoctest
1120 julia> A = [1 2; 3 4]
1121 2×2 Array{Int64,2}:
1122 1 2
1123 3 4
1124
1125 julia> V = view(A, 1:2, :)
1126 2×2 view(::Array{Int64,2}, 1:2, :) with eltype Int64:
1127 1 2
1128 3 4
1129
1130 julia> parent(V)
1131 2×2 Array{Int64,2}:
1132 1 2
1133 3 4
1134 ```
1135 """
1136 parent(a::AbstractArray) = a
1137
1138 ## rudimentary aliasing detection ##
1139 """
1140 Base.unalias(dest, A)
1141
1142 Return either `A` or a copy of `A` in a rough effort to prevent modifications to `dest` from
1143 affecting the returned object. No guarantees are provided.
1144
1145 Custom arrays that wrap or use fields containing arrays that might alias against other
1146 external objects should provide a [`Base.dataids`](@ref) implementation.
1147
1148 This function must return an object of exactly the same type as `A` for performance and type
1149 stability. Mutable custom arrays for which [`copy(A)`](@ref) is not `typeof(A)` should
1150 provide a [`Base.unaliascopy`](@ref) implementation.
1151
1152 See also [`Base.mightalias`](@ref).
1153 """
1154 unalias(dest, A::AbstractArray) = mightalias(dest, A) ? unaliascopy(A) : A
1155 unalias(dest, A::AbstractRange) = A
1156 unalias(dest, A) = A
1157
1158 """
1159 Base.unaliascopy(A)
1160
1161 Make a preventative copy of `A` in an operation where `A` [`Base.mightalias`](@ref) against
1162 another array in order to preserve consistent semantics as that other array is mutated.
1163
1164 This must return an object of the same type as `A` to preserve optimal performance in the
1165 much more common case where aliasing does not occur. By default,
1166 `unaliascopy(A::AbstractArray)` will attempt to use [`copy(A)`](@ref), but in cases where
1167 `copy(A)` is not a `typeof(A)`, then the array should provide a custom implementation of
1168 `Base.unaliascopy(A)`.
1169 """
1170 unaliascopy(A::Array) = copy(A)
1171 unaliascopy(A::AbstractArray)::typeof(A) = (@_noinline_meta; _unaliascopy(A, copy(A)))
1172 _unaliascopy(A::T, C::T) where {T} = C
1173 _unaliascopy(A, C) = throw(ArgumentError("""
1174 an array of type `$(typeof(A).name)` shares memory with another argument and must
1175 make a preventative copy of itself in order to maintain consistent semantics,
1176 but `copy(A)` returns a new array of type `$(typeof(C))`. To fix, implement:
1177 `Base.unaliascopy(A::$(typeof(A).name))::typeof(A)`"""))
1178 unaliascopy(A) = A
1179
1180 """
1181 Base.mightalias(A::AbstractArray, B::AbstractArray)
1182
1183 Perform a conservative test to check if arrays `A` and `B` might share the same memory.
1184
1185 By default, this simply checks if either of the arrays reference the same memory
1186 regions, as identified by their [`Base.dataids`](@ref).
1187 """
1188 mightalias(A::AbstractArray, B::AbstractArray) = !isbits(A) && !isbits(B) && !_isdisjoint(dataids(A), dataids(B))
1189 mightalias(x, y) = false
1190
1191 _isdisjoint(as::Tuple{}, bs::Tuple{}) = true
1192 _isdisjoint(as::Tuple{}, bs::Tuple{UInt}) = true
1193 _isdisjoint(as::Tuple{}, bs::Tuple) = true
1194 _isdisjoint(as::Tuple{UInt}, bs::Tuple{}) = true
1195 _isdisjoint(as::Tuple{UInt}, bs::Tuple{UInt}) = as[1] != bs[1]
1196 _isdisjoint(as::Tuple{UInt}, bs::Tuple) = !(as[1] in bs)
1197 _isdisjoint(as::Tuple, bs::Tuple{}) = true
1198 _isdisjoint(as::Tuple, bs::Tuple{UInt}) = !(bs[1] in as)
1199 _isdisjoint(as::Tuple, bs::Tuple) = !(as[1] in bs) && _isdisjoint(tail(as), bs)
1200
1201 """
1202 Base.dataids(A::AbstractArray)
1203
1204 Return a tuple of `UInt`s that represent the mutable data segments of an array.
1205
1206 Custom arrays that would like to opt-in to aliasing detection of their component
1207 parts can specialize this method to return the concatenation of the `dataids` of
1208 their component parts. A typical definition for an array that wraps a parent is
1209 `Base.dataids(C::CustomArray) = dataids(C.parent)`.
1210 """
1211 dataids(A::AbstractArray) = (UInt(objectid(A)),)
1212 dataids(A::Array) = (UInt(pointer(A)),)
1213 dataids(::AbstractRange) = ()
1214 dataids(x) = ()
1215
1216 ## get (getindex with a default value) ##
1217
1218 RangeVecIntList{A<:AbstractVector{Int}} = Union{Tuple{Vararg{Union{AbstractRange, AbstractVector{Int}}}},
1219 AbstractVector{UnitRange{Int}}, AbstractVector{AbstractRange{Int}}, AbstractVector{A}}
1220
1221 get(A::AbstractArray, i::Integer, default) = checkbounds(Bool, A, i) ? A[i] : default
1222 get(A::AbstractArray, I::Tuple{}, default) = checkbounds(Bool, A) ? A[] : default
1223 get(A::AbstractArray, I::Dims, default) = checkbounds(Bool, A, I...) ? A[I...] : default
1224
1225 function get!(X::AbstractVector{T}, A::AbstractVector, I::Union{AbstractRange,AbstractVector{Int}}, default::T) where T
1226 # 1d is not linear indexing
1227 ind = findall(in(axes1(A)), I)
1228 X[ind] = A[I[ind]]
1229 Xind = axes1(X)
1230 X[first(Xind):first(ind)-1] = default
1231 X[last(ind)+1:last(Xind)] = default
1232 X
1233 end
1234 function get!(X::AbstractArray{T}, A::AbstractArray, I::Union{AbstractRange,AbstractVector{Int}}, default::T) where T
1235 # Linear indexing
1236 ind = findall(in(1:length(A)), I)
1237 X[ind] = A[I[ind]]
1238 fill!(view(X, 1:first(ind)-1), default)
1239 fill!(view(X, last(ind)+1:length(X)), default)
1240 X
1241 end
1242
1243 get(A::AbstractArray, I::AbstractRange, default) = get!(similar(A, typeof(default), index_shape(I)), A, I, default)
1244
1245 function get!(X::AbstractArray{T}, A::AbstractArray, I::RangeVecIntList, default::T) where T
1246 fill!(X, default)
1247 dst, src = indcopy(size(A), I)
1248 X[dst...] = A[src...]
1249 X
1250 end
1251
1252 get(A::AbstractArray, I::RangeVecIntList, default) =
1253 get!(similar(A, typeof(default), index_shape(I...)), A, I, default)
1254
1255 ## structured matrix methods ##
1256 replace_in_print_matrix(A::AbstractMatrix,i::Integer,j::Integer,s::AbstractString) = s
1257 replace_in_print_matrix(A::AbstractVector,i::Integer,j::Integer,s::AbstractString) = s
1258
1259 ## Concatenation ##
1260 eltypeof(x) = typeof(x)
1261 eltypeof(x::AbstractArray) = eltype(x)
1262
1263 promote_eltypeof() = Bottom
1264 promote_eltypeof(v1, vs...) = promote_type(eltypeof(v1), promote_eltypeof(vs...))
1265
1266 promote_eltype() = Bottom
1267 promote_eltype(v1, vs...) = promote_type(eltype(v1), promote_eltype(vs...))
1268
1269 #TODO: ERROR CHECK
1270 _cat(catdim::Integer) = Vector{Any}()
1271
1272 typed_vcat(::Type{T}) where {T} = Vector{T}()
1273 typed_hcat(::Type{T}) where {T} = Vector{T}()
1274
1275 ## cat: special cases
1276 vcat(X::T...) where {T} = T[ X[i] for i=1:length(X) ]
1277 vcat(X::T...) where {T<:Number} = T[ X[i] for i=1:length(X) ]
1278 hcat(X::T...) where {T} = T[ X[j] for i=1:1, j=1:length(X) ]
1279 hcat(X::T...) where {T<:Number} = T[ X[j] for i=1:1, j=1:length(X) ]
1280
1281 vcat(X::Number...) = hvcat_fill(Vector{promote_typeof(X...)}(undef, length(X)), X)
1282 hcat(X::Number...) = hvcat_fill(Matrix{promote_typeof(X...)}(undef, 1,length(X)), X)
1283 typed_vcat(::Type{T}, X::Number...) where {T} = hvcat_fill(Vector{T}(undef, length(X)), X)
1284 typed_hcat(::Type{T}, X::Number...) where {T} = hvcat_fill(Matrix{T}(undef, 1,length(X)), X)
1285
1286 vcat(V::AbstractVector...) = typed_vcat(promote_eltype(V...), V...)
1287 vcat(V::AbstractVector{T}...) where {T} = typed_vcat(T, V...)
1288
1289 # FIXME: this alias would better be Union{AbstractVector{T}, Tuple{Vararg{T}}}
1290 # and method signatures should do AbstractVecOrTuple{<:T} when they want covariance,
1291 # but that solution currently fails (see #27188 and #27224)
1292 AbstractVecOrTuple{T} = Union{AbstractVector{<:T}, Tuple{Vararg{T}}}
1293
1294 function _typed_vcat(::Type{T}, V::AbstractVecOrTuple{AbstractVector}) where T
1295 n::Int = 0
1296 for Vk in V
1297 n += length(Vk)
1298 end
1299 a = similar(V[1], T, n)
1300 pos = 1
1301 for k=1:length(V)
1302 Vk = V[k]
1303 p1 = pos+length(Vk)-1
1304 a[pos:p1] = Vk
1305 pos = p1+1
1306 end
1307 a
1308 end
1309
1310 typed_hcat(::Type{T}, A::AbstractVecOrMat...) where {T} = _typed_hcat(T, A)
1311
1312 hcat(A::AbstractVecOrMat...) = typed_hcat(promote_eltype(A...), A...)
1313 hcat(A::AbstractVecOrMat{T}...) where {T} = typed_hcat(T, A...)
1314
1315 function _typed_hcat(::Type{T}, A::AbstractVecOrTuple{AbstractVecOrMat}) where T
1316 nargs = length(A)
1317 nrows = size(A[1], 1)
1318 ncols = 0
1319 dense = true
1320 for j = 1:nargs
1321 Aj = A[j]
1322 if size(Aj, 1) != nrows
1323 throw(ArgumentError("number of rows of each array must match (got $(map(x->size(x,1), A)))"))
1324 end
1325 dense &= isa(Aj,Array)
1326 nd = ndims(Aj)
1327 ncols += (nd==2 ? size(Aj,2) : 1)
1328 end
1329 B = similar(A[1], T, nrows, ncols)
1330 pos = 1
1331 if dense
1332 for k=1:nargs
1333 Ak = A[k]
1334 n = length(Ak)
1335 copyto!(B, pos, Ak, 1, n)
1336 pos += n
1337 end
1338 else
1339 for k=1:nargs
1340 Ak = A[k]
1341 p1 = pos+(isa(Ak,AbstractMatrix) ? size(Ak, 2) : 1)-1
1342 B[:, pos:p1] = Ak
1343 pos = p1+1
1344 end
1345 end
1346 return B
1347 end
1348
1349 vcat(A::AbstractVecOrMat...) = typed_vcat(promote_eltype(A...), A...)
1350 vcat(A::AbstractVecOrMat{T}...) where {T} = typed_vcat(T, A...)
1351
1352 function _typed_vcat(::Type{T}, A::AbstractVecOrTuple{AbstractVecOrMat}) where T
1353 nargs = length(A)
1354 nrows = sum(a->size(a, 1), A)::Int
1355 ncols = size(A[1], 2)
1356 for j = 2:nargs
1357 if size(A[j], 2) != ncols
1358 throw(ArgumentError("number of columns of each array must match (got $(map(x->size(x,2), A)))"))
1359 end
1360 end
1361 B = similar(A[1], T, nrows, ncols)
1362 pos = 1
1363 for k=1:nargs
1364 Ak = A[k]
1365 p1 = pos+size(Ak,1)-1
1366 B[pos:p1, :] = Ak
1367 pos = p1+1
1368 end
1369 return B
1370 end
1371
1372 typed_vcat(::Type{T}, A::AbstractVecOrMat...) where {T} = _typed_vcat(T, A)
1373
1374 reduce(::typeof(vcat), A::AbstractVector{<:AbstractVecOrMat}) =
1375 _typed_vcat(mapreduce(eltype, promote_type, A), A)
1376
1377 reduce(::typeof(hcat), A::AbstractVector{<:AbstractVecOrMat}) =
1378 _typed_hcat(mapreduce(eltype, promote_type, A), A)
1379
1380 ## cat: general case
1381
1382 # helper functions
1383 cat_size(A) = (1,)
1384 cat_size(A::AbstractArray) = size(A)
1385 cat_size(A, d) = 1
1386 cat_size(A::AbstractArray, d) = size(A, d)
1387
1388 cat_indices(A, d) = OneTo(1)
1389 cat_indices(A::AbstractArray, d) = axes(A, d)
1390
1391 cat_similar(A, T, shape) = Array{T}(undef, shape)
1392 cat_similar(A::AbstractArray, T, shape) = similar(A, T, shape)
1393
1394 cat_shape(dims, shape::Tuple) = shape
1395 @inline cat_shape(dims, shape::Tuple, nshape::Tuple, shapes::Tuple...) =
1396 cat_shape(dims, _cshp(1, dims, shape, nshape), shapes...)
1397
1398 _cshp(ndim::Int, ::Tuple{}, ::Tuple{}, ::Tuple{}) = ()
1399 _cshp(ndim::Int, ::Tuple{}, ::Tuple{}, nshape) = nshape
1400 _cshp(ndim::Int, dims, ::Tuple{}, ::Tuple{}) = ntuple(b -> 1, Val(length(dims)))
1401 @inline _cshp(ndim::Int, dims, shape, ::Tuple{}) =
1402 (shape[1] + dims[1], _cshp(ndim + 1, tail(dims), tail(shape), ())...)
1403 @inline _cshp(ndim::Int, dims, ::Tuple{}, nshape) =
1404 (nshape[1], _cshp(ndim + 1, tail(dims), (), tail(nshape))...)
1405 @inline function _cshp(ndim::Int, ::Tuple{}, shape, ::Tuple{})
1406 _cs(ndim, shape[1], 1)
1407 (1, _cshp(ndim + 1, (), tail(shape), ())...)
1408 end
1409 @inline function _cshp(ndim::Int, ::Tuple{}, shape, nshape)
1410 next = _cs(ndim, shape[1], nshape[1])
1411 (next, _cshp(ndim + 1, (), tail(shape), tail(nshape))...)
1412 end
1413 @inline function _cshp(ndim::Int, dims, shape, nshape)
1414 a = shape[1]
1415 b = nshape[1]
1416 next = dims[1] ? a + b : _cs(ndim, a, b)
1417 (next, _cshp(ndim + 1, tail(dims), tail(shape), tail(nshape))...)
1418 end
1419
1420 _cs(d, a, b) = (a == b ? a : throw(DimensionMismatch(
1421 "mismatch in dimension $d (expected $a got $b)")))
1422
1423 function dims2cat(::Val{n}) where {n}
1424 n <= 0 && throw(ArgumentError("cat dimension must be a positive integer, but got $n"))
1425 ntuple(i -> (i == n), Val(n))
1426 end
1427
1428 function dims2cat(dims)
1429 if any(dims .<= 0)
1430 throw(ArgumentError("All cat dimensions must be positive integers, but got $dims"))
1431 end
1432 ntuple(in(dims), maximum(dims))
1433 end
1434
1435 _cat(dims, X...) = cat_t(promote_eltypeof(X...), X...; dims=dims)
1436
1437 @inline cat_t(::Type{T}, X...; dims) where {T} = _cat_t(dims, T, X...)
1438 @inline function _cat_t(dims, T::Type, X...)
1439 catdims = dims2cat(dims)
1440 shape = cat_shape(catdims, (), map(cat_size, X)...)
1441 A = cat_similar(X[1], T, shape)
1442 if T <: Number && count(!iszero, catdims) > 1
1443 fill!(A, zero(T))
1444 end
1445 return __cat(A, shape, catdims, X...)
1446 end
1447
1448 function __cat(A, shape::NTuple{N}, catdims, X...) where N
1449 offsets = zeros(Int, N)
1450 inds = Vector{UnitRange{Int}}(undef, N)
1451 concat = copyto!(zeros(Bool, N), catdims)
1452 for x in X
1453 for i = 1:N
1454 if concat[i]
1455 inds[i] = offsets[i] .+ cat_indices(x, i)
1456 offsets[i] += cat_size(x, i)
1457 else
1458 inds[i] = 1:shape[i]
1459 end
1460 end
1461 I::NTuple{N, UnitRange{Int}} = (inds...,)
1462 if x isa AbstractArray
1463 A[I...] = x
1464 else
1465 fill!(view(A, I...), x)
1466 end
1467 end
1468 return A
1469 end
1470
1471 """
1472 vcat(A...)
1473
1474 Concatenate along dimension 1.
1475
1476 # Examples
1477 ```jldoctest
1478 julia> a = [1 2 3 4 5]
1479 1×5 Array{Int64,2}:
1480 1 2 3 4 5
1481
1482 julia> b = [6 7 8 9 10; 11 12 13 14 15]
1483 2×5 Array{Int64,2}:
1484 6 7 8 9 10
1485 11 12 13 14 15
1486
1487 julia> vcat(a,b)
1488 3×5 Array{Int64,2}:
1489 1 2 3 4 5
1490 6 7 8 9 10
1491 11 12 13 14 15
1492
1493 julia> c = ([1 2 3], [4 5 6])
1494 ([1 2 3], [4 5 6])
1495
1496 julia> vcat(c...)
1497 2×3 Array{Int64,2}:
1498 1 2 3
1499 4 5 6
1500 ```
1501 """
1502 vcat(X...) = cat(X...; dims=Val(1))
1503 """
1504 hcat(A...)
1505
1506 Concatenate along dimension 2.
1507
1508 # Examples
1509 ```jldoctest
1510 julia> a = [1; 2; 3; 4; 5]
1511 5-element Array{Int64,1}:
1512 1
1513 2
1514 3
1515 4
1516 5
1517
1518 julia> b = [6 7; 8 9; 10 11; 12 13; 14 15]
1519 5×2 Array{Int64,2}:
1520 6 7
1521 8 9
1522 10 11
1523 12 13
1524 14 15
1525
1526 julia> hcat(a,b)
1527 5×3 Array{Int64,2}:
1528 1 6 7
1529 2 8 9
1530 3 10 11
1531 4 12 13
1532 5 14 15
1533
1534 julia> c = ([1; 2; 3], [4; 5; 6])
1535 ([1, 2, 3], [4, 5, 6])
1536
1537 julia> hcat(c...)
1538 3×2 Array{Int64,2}:
1539 1 4
1540 2 5
1541 3 6
1542
1543 julia> x = Matrix(undef, 3, 0) # x = [] would have created an Array{Any, 1}, but need an Array{Any, 2}
1544 3×0 Array{Any,2}
1545
1546 julia> hcat(x, [1; 2; 3])
1547 3×1 Array{Any,2}:
1548 1
1549 2
1550 3
1551 ```
1552 """
1553 hcat(X...) = cat(X...; dims=Val(2))
1554
1555 typed_vcat(T::Type, X...) = cat_t(T, X...; dims=Val(1))
1556 typed_hcat(T::Type, X...) = cat_t(T, X...; dims=Val(2))
1557
1558 """
1559 cat(A...; dims=dims)
1560
1561 Concatenate the input arrays along the specified dimensions in the iterable `dims`. For
1562 dimensions not in `dims`, all input arrays should have the same size, which will also be the
1563 size of the output array along that dimension. For dimensions in `dims`, the size of the
1564 output array is the sum of the sizes of the input arrays along that dimension. If `dims` is
1565 a single number, the different arrays are tightly stacked along that dimension. If `dims` is
1566 an iterable containing several dimensions, this allows one to construct block diagonal
1567 matrices and their higher-dimensional analogues by simultaneously increasing several
1568 dimensions for every new input array and putting zero blocks elsewhere. For example,
1569 `cat(matrices...; dims=(1,2))` builds a block diagonal matrix, i.e. a block matrix with
1570 `matrices[1]`, `matrices[2]`, ... as diagonal blocks and matching zero blocks away from the
1571 diagonal.
1572 """
1573 @inline cat(A...; dims) = _cat(dims, A...)
1574 _cat(catdims, A::AbstractArray{T}...) where {T} = cat_t(T, A...; dims=catdims)
1575
1576 # The specializations for 1 and 2 inputs are important
1577 # especially when running with --inline=no, see #11158
1578 vcat(A::AbstractArray) = cat(A; dims=Val(1))
1579 vcat(A::AbstractArray, B::AbstractArray) = cat(A, B; dims=Val(1))
1580 vcat(A::AbstractArray...) = cat(A...; dims=Val(1))
1581 hcat(A::AbstractArray) = cat(A; dims=Val(2))
1582 hcat(A::AbstractArray, B::AbstractArray) = cat(A, B; dims=Val(2))
1583 hcat(A::AbstractArray...) = cat(A...; dims=Val(2))
1584
1585 typed_vcat(T::Type, A::AbstractArray) = cat_t(T, A; dims=Val(1))
1586 typed_vcat(T::Type, A::AbstractArray, B::AbstractArray) = cat_t(T, A, B; dims=Val(1))
1587 typed_vcat(T::Type, A::AbstractArray...) = cat_t(T, A...; dims=Val(1))
1588 typed_hcat(T::Type, A::AbstractArray) = cat_t(T, A; dims=Val(2))
1589 typed_hcat(T::Type, A::AbstractArray, B::AbstractArray) = cat_t(T, A, B; dims=Val(2))
1590 typed_hcat(T::Type, A::AbstractArray...) = cat_t(T, A...; dims=Val(2))
1591
1592 # 2d horizontal and vertical concatenation
1593
1594 function hvcat(nbc::Integer, as...)
1595 # nbc = # of block columns
1596 n = length(as)
1597 mod(n,nbc) != 0 &&
1598 throw(ArgumentError("number of arrays $n is not a multiple of the requested number of block columns $nbc"))
1599 nbr = div(n,nbc)
1600 hvcat(ntuple(i->nbc, nbr), as...)
1601 end
1602
1603 """
1604 hvcat(rows::Tuple{Vararg{Int}}, values...)
1605
1606 Horizontal and vertical concatenation in one call. This function is called for block matrix
1607 syntax. The first argument specifies the number of arguments to concatenate in each block
1608 row.
1609
1610 # Examples
1611 ```jldoctest
1612 julia> a, b, c, d, e, f = 1, 2, 3, 4, 5, 6
1613 (1, 2, 3, 4, 5, 6)
1614
1615 julia> [a b c; d e f]
1616 2×3 Array{Int64,2}:
1617 1 2 3
1618 4 5 6
1619
1620 julia> hvcat((3,3), a,b,c,d,e,f)
1621 2×3 Array{Int64,2}:
1622 1 2 3
1623 4 5 6
1624
1625 julia> [a b;c d; e f]
1626 3×2 Array{Int64,2}:
1627 1 2
1628 3 4
1629 5 6
1630
1631 julia> hvcat((2,2,2), a,b,c,d,e,f)
1632 3×2 Array{Int64,2}:
1633 1 2
1634 3 4
1635 5 6
1636 ```
1637
1638 If the first argument is a single integer `n`, then all block rows are assumed to have `n`
1639 block columns.
1640 """
1641 hvcat(rows::Tuple{Vararg{Int}}, xs::AbstractVecOrMat...) = typed_hvcat(promote_eltype(xs...), rows, xs...)
1642 hvcat(rows::Tuple{Vararg{Int}}, xs::AbstractVecOrMat{T}...) where {T} = typed_hvcat(T, rows, xs...)
1643
1644 function typed_hvcat(::Type{T}, rows::Tuple{Vararg{Int}}, as::AbstractVecOrMat...) where T
1645 nbr = length(rows) # number of block rows
1646
1647 nc = 0
1648 for i=1:rows[1]
1649 nc += size(as[i],2)
1650 end
1651
1652 nr = 0
1653 a = 1
1654 for i = 1:nbr
1655 nr += size(as[a],1)
1656 a += rows[i]
1657 end
1658
1659 out = similar(as[1], T, nr, nc)
1660
1661 a = 1
1662 r = 1
1663 for i = 1:nbr
1664 c = 1
1665 szi = size(as[a],1)
1666 for j = 1:rows[i]
1667 Aj = as[a+j-1]
1668 szj = size(Aj,2)
1669 if size(Aj,1) != szi
1670 throw(ArgumentError("mismatched height in block row $(i) (expected $szi, got $(size(Aj,1)))"))
1671 end
1672 if c-1+szj > nc
1673 throw(ArgumentError("block row $(i) has mismatched number of columns (expected $nc, got $(c-1+szj))"))
1674 end
1675 out[r:r-1+szi, c:c-1+szj] = Aj
1676 c += szj
1677 end
1678 if c != nc+1
1679 throw(ArgumentError("block row $(i) has mismatched number of columns (expected $nc, got $(c-1))"))
1680 end
1681 r += szi
1682 a += rows[i]
1683 end
1684 out
1685 end
1686
1687 hvcat(rows::Tuple{Vararg{Int}}) = []
1688 typed_hvcat(::Type{T}, rows::Tuple{Vararg{Int}}) where {T} = Vector{T}()
1689
1690 function hvcat(rows::Tuple{Vararg{Int}}, xs::T...) where T<:Number
1691 nr = length(rows)
1692 nc = rows[1]
1693
1694 a = Matrix{T}(undef, nr, nc)
1695 if length(a) != length(xs)
1696 throw(ArgumentError("argument count does not match specified shape (expected $(length(a)), got $(length(xs)))"))
1697 end
1698 k = 1
1699 @inbounds for i=1:nr
1700 if nc != rows[i]
1701 throw(ArgumentError("row $(i) has mismatched number of columns (expected $nc, got $(rows[i]))"))
1702 end
1703 for j=1:nc
1704 a[i,j] = xs[k]
1705 k += 1
1706 end
1707 end
1708 a
1709 end
1710
1711 function hvcat_fill(a::Array, xs::Tuple)
1712 k = 1
1713 nr, nc = size(a,1), size(a,2)
1714 for i=1:nr
1715 @inbounds for j=1:nc
1716 a[i,j] = xs[k]
1717 k += 1
1718 end
1719 end
1720 a
1721 end
1722
1723 hvcat(rows::Tuple{Vararg{Int}}, xs::Number...) = typed_hvcat(promote_typeof(xs...), rows, xs...)
1724 hvcat(rows::Tuple{Vararg{Int}}, xs...) = typed_hvcat(promote_eltypeof(xs...), rows, xs...)
1725
1726 function typed_hvcat(::Type{T}, rows::Tuple{Vararg{Int}}, xs::Number...) where T
1727 nr = length(rows)
1728 nc = rows[1]
1729 for i = 2:nr
1730 if nc != rows[i]
1731 throw(ArgumentError("row $(i) has mismatched number of columns (expected $nc, got $(rows[i]))"))
1732 end
1733 end
1734 len = length(xs)
1735 if nr*nc != len
1736 throw(ArgumentError("argument count $(len) does not match specified shape $((nr,nc))"))
1737 end
1738 hvcat_fill(Matrix{T}(undef, nr, nc), xs)
1739 end
1740
1741 function typed_hvcat(::Type{T}, rows::Tuple{Vararg{Int}}, as...) where T
1742 nbr = length(rows) # number of block rows
1743 rs = Vector{Any}(undef, nbr)
1744 a = 1
1745 for i = 1:nbr
1746 rs[i] = typed_hcat(T, as[a:a-1+rows[i]]...)
1747 a += rows[i]
1748 end
1749 T[rs...;]
1750 end
1751
1752 ## Reductions and accumulates ##
1753
1754 function isequal(A::AbstractArray, B::AbstractArray)
1755 if A === B return true end
1756 if axes(A) != axes(B)
1757 return false
1758 end
1759 for (a, b) in zip(A, B)
1760 if !isequal(a, b)
1761 return false
1762 end
1763 end
1764 return true
1765 end
1766
1767 function cmp(A::AbstractVector, B::AbstractVector)
1768 for (a, b) in zip(A, B)
1769 if !isequal(a, b)
1770 return isless(a, b) ? -1 : 1
1771 end
1772 end
1773 return cmp(length(A), length(B))
1774 end
1775
1776 isless(A::AbstractVector, B::AbstractVector) = cmp(A, B) < 0
1777
1778 function (==)(A::AbstractArray, B::AbstractArray)
1779 if axes(A) != axes(B)
1780 return false
1781 end
1782 anymissing = false
1783 for (a, b) in zip(A, B)
1784 eq = (a == b)
1785 if ismissing(eq)
1786 anymissing = true
1787 elseif !eq
1788 return false
1789 end
1790 end
1791 return anymissing ? missing : true
1792 end
1793
1794 # _sub2ind and _ind2sub
1795 # fallbacks
1796 function _sub2ind(A::AbstractArray, I...)
1797 @_inline_meta
1798 _sub2ind(axes(A), I...)
1799 end
1800
1801 function _ind2sub(A::AbstractArray, ind)
1802 @_inline_meta
1803 _ind2sub(axes(A), ind)
1804 end
1805
1806 # 0-dimensional arrays and indexing with []
1807 _sub2ind(::Tuple{}) = 1
1808 _sub2ind(::DimsInteger) = 1
1809 _sub2ind(::Indices) = 1
1810 _sub2ind(::Tuple{}, I::Integer...) = (@_inline_meta; _sub2ind_recurse((), 1, 1, I...))
1811
1812 # Generic cases
1813 _sub2ind(dims::DimsInteger, I::Integer...) = (@_inline_meta; _sub2ind_recurse(dims, 1, 1, I...))
1814 _sub2ind(inds::Indices, I::Integer...) = (@_inline_meta; _sub2ind_recurse(inds, 1, 1, I...))
1815 # In 1d, there's a question of whether we're doing cartesian indexing
1816 # or linear indexing. Support only the former.
1817 _sub2ind(inds::Indices{1}, I::Integer...) =
1818 throw(ArgumentError("Linear indexing is not defined for one-dimensional arrays"))
1819 _sub2ind(inds::Tuple{OneTo}, I::Integer...) = (@_inline_meta; _sub2ind_recurse(inds, 1, 1, I...)) # only OneTo is safe
1820 _sub2ind(inds::Tuple{OneTo}, i::Integer) = i
1821
1822 _sub2ind_recurse(::Any, L, ind) = ind
1823 function _sub2ind_recurse(::Tuple{}, L, ind, i::Integer, I::Integer...)
1824 @_inline_meta
1825 _sub2ind_recurse((), L, ind+(i-1)*L, I...)
1826 end
1827 function _sub2ind_recurse(inds, L, ind, i::Integer, I::Integer...)
1828 @_inline_meta
1829 r1 = inds[1]
1830 _sub2ind_recurse(tail(inds), nextL(L, r1), ind+offsetin(i, r1)*L, I...)
1831 end
1832
1833 nextL(L, l::Integer) = L*l
1834 nextL(L, r::AbstractUnitRange) = L*unsafe_length(r)
1835 nextL(L, r::Slice) = L*unsafe_length(r.indices)
1836 offsetin(i, l::Integer) = i-1
1837 offsetin(i, r::AbstractUnitRange) = i-first(r)
1838
1839 _ind2sub(::Tuple{}, ind::Integer) = (@_inline_meta; ind == 1 ? () : throw(BoundsError()))
1840 _ind2sub(dims::DimsInteger, ind::Integer) = (@_inline_meta; _ind2sub_recurse(dims, ind-1))
1841 _ind2sub(inds::Indices, ind::Integer) = (@_inline_meta; _ind2sub_recurse(inds, ind-1))
1842 _ind2sub(inds::Indices{1}, ind::Integer) =
1843 throw(ArgumentError("Linear indexing is not defined for one-dimensional arrays"))
1844 _ind2sub(inds::Tuple{OneTo}, ind::Integer) = (ind,)
1845
1846 _ind2sub_recurse(::Tuple{}, ind) = (ind+1,)
1847 function _ind2sub_recurse(indslast::NTuple{1}, ind)
1848 @_inline_meta
1849 (_lookup(ind, indslast[1]),)
1850 end
1851 function _ind2sub_recurse(inds, ind)
1852 @_inline_meta
1853 r1 = inds[1]
1854 indnext, f, l = _div(ind, r1)
1855 (ind-l*indnext+f, _ind2sub_recurse(tail(inds), indnext)...)
1856 end
1857
1858 _lookup(ind, d::Integer) = ind+1
1859 _lookup(ind, r::AbstractUnitRange) = ind+first(r)
1860 _div(ind, d::Integer) = div(ind, d), 1, d
1861 _div(ind, r::AbstractUnitRange) = (d = unsafe_length(r); (div(ind, d), first(r), d))
1862
1863 # Vectorized forms
1864 function _sub2ind(inds::Indices{1}, I1::AbstractVector{T}, I::AbstractVector{T}...) where T<:Integer
1865 throw(ArgumentError("Linear indexing is not defined for one-dimensional arrays"))
1866 end
1867 _sub2ind(inds::Tuple{OneTo}, I1::AbstractVector{T}, I::AbstractVector{T}...) where {T<:Integer} =
1868 _sub2ind_vecs(inds, I1, I...)
1869 _sub2ind(inds::Union{DimsInteger,Indices}, I1::AbstractVector{T}, I::AbstractVector{T}...) where {T<:Integer} =
1870 _sub2ind_vecs(inds, I1, I...)
1871 function _sub2ind_vecs(inds, I::AbstractVector...)
1872 I1 = I[1]
1873 Iinds = axes1(I1)
1874 for j = 2:length(I)
1875 axes1(I[j]) == Iinds || throw(DimensionMismatch("indices of I[1] ($(Iinds)) does not match indices of I[$j] ($(axes1(I[j])))"))
1876 end
1877 Iout = similar(I1)
1878 _sub2ind!(Iout, inds, Iinds, I)
1879 Iout
1880 end
1881
1882 function _sub2ind!(Iout, inds, Iinds, I)
1883 @_noinline_meta
1884 for i in Iinds
1885 # Iout[i] = _sub2ind(inds, map(Ij -> Ij[i], I)...)
1886 Iout[i] = sub2ind_vec(inds, i, I)
1887 end
1888 Iout
1889 end
1890
1891 sub2ind_vec(inds, i, I) = (@_inline_meta; _sub2ind(inds, _sub2ind_vec(i, I...)...))
1892 _sub2ind_vec(i, I1, I...) = (@_inline_meta; (I1[i], _sub2ind_vec(i, I...)...))
1893 _sub2ind_vec(i) = ()
1894
1895 function _ind2sub(inds::Union{DimsInteger{N},Indices{N}}, ind::AbstractVector{<:Integer}) where N
1896 M = length(ind)
1897 t = ntuple(n->similar(ind),Val(N))
1898 for (i,idx) in pairs(IndexLinear(), ind)
1899 sub = _ind2sub(inds, idx)
1900 for j = 1:N
1901 t[j][i] = sub[j]
1902 end
1903 end
1904 t
1905 end
1906
1907 ## iteration utilities ##
1908
1909 """
1910 foreach(f, c...) -> Nothing
1911
1912 Call function `f` on each element of iterable `c`.
1913 For multiple iterable arguments, `f` is called elementwise.
1914 `foreach` should be used instead of `map` when the results of `f` are not
1915 needed, for example in `foreach(println, array)`.
1916
1917 # Examples
1918 ```jldoctest
1919 julia> a = 1:3:7;
1920
1921 julia> foreach(x -> println(x^2), a)
1922 1
1923 16
1924 49
1925 ```
1926 """
1927 foreach(f) = (f(); nothing)
1928 foreach(f, itr) = (for x in itr; f(x); end; nothing)
1929 foreach(f, itrs...) = (for z in zip(itrs...); f(z...); end; nothing)
1930
1931 ## map over arrays ##
1932
1933 ## transform any set of dimensions
1934 ## dims specifies which dimensions will be transformed. for example
1935 ## dims==1:2 will call f on all slices A[:,:,...]
1936 """
1937 mapslices(f, A; dims)
1938
1939 Transform the given dimensions of array `A` using function `f`. `f` is called on each slice
1940 of `A` of the form `A[...,:,...,:,...]`. `dims` is an integer vector specifying where the
1941 colons go in this expression. The results are concatenated along the remaining dimensions.
1942 For example, if `dims` is `[1,2]` and `A` is 4-dimensional, `f` is called on `A[:,:,i,j]`
1943 for all `i` and `j`.
1944
1945 # Examples
1946 ```jldoctest
1947 julia> a = reshape(Vector(1:16),(2,2,2,2))
1948 2×2×2×2 Array{Int64,4}:
1949 [:, :, 1, 1] =
1950 1 3
1951 2 4
1952
1953 [:, :, 2, 1] =
1954 5 7
1955 6 8
1956
1957 [:, :, 1, 2] =
1958 9 11
1959 10 12
1960
1961 [:, :, 2, 2] =
1962 13 15
1963 14 16
1964
1965 julia> mapslices(sum, a, dims = [1,2])
1966 1×1×2×2 Array{Int64,4}:
1967 [:, :, 1, 1] =
1968 10
1969
1970 [:, :, 2, 1] =
1971 26
1972
1973 [:, :, 1, 2] =
1974 42
1975
1976 [:, :, 2, 2] =
1977 58
1978 ```
1979 """
1980 function mapslices(f, A::AbstractArray; dims)
1981 if isempty(dims)
1982 return map(f,A)
1983 end
1984 if !isa(dims, AbstractVector)
1985 dims = [dims...]
1986 end
1987
1988 dimsA = [axes(A)...]
1989 ndimsA = ndims(A)
1990 alldims = [1:ndimsA;]
1991
1992 otherdims = setdiff(alldims, dims)
1993
1994 idx = Any[first(ind) for ind in axes(A)]
1995 itershape = tuple(dimsA[otherdims]...)
1996 for d in dims
1997 idx[d] = Slice(axes(A, d))
1998 end
1999
2000 # Apply the function to the first slice in order to determine the next steps
2001 Aslice = A[idx...]
2002 r1 = f(Aslice)
2003 # In some cases, we can re-use the first slice for a dramatic performance
2004 # increase. The slice itself must be mutable and the result cannot contain
2005 # any mutable containers. The following errs on the side of being overly
2006 # strict (#18570 & #21123).
2007 safe_for_reuse = isa(Aslice, StridedArray) &&
2008 (isa(r1, Number) || (isa(r1, AbstractArray) && eltype(r1) <: Number))
2009
2010 # determine result size and allocate
2011 Rsize = copy(dimsA)
2012 # TODO: maybe support removing dimensions
2013 if !isa(r1, AbstractArray) || ndims(r1) == 0
2014 # If the result of f on a single slice is a scalar then we add singleton
2015 # dimensions. When adding the dimensions, we have to respect the
2016 # index type of the input array (e.g. in the case of OffsetArrays)
2017 tmp = similar(Aslice, typeof(r1), reduced_indices(Aslice, 1:ndims(Aslice)))
2018 tmp[firstindex(tmp)] = r1
2019 r1 = tmp
2020 end
2021 nextra = max(0, length(dims)-ndims(r1))
2022 if eltype(Rsize) == Int
2023 Rsize[dims] = [size(r1)..., ntuple(d->1, nextra)...]
2024 else
2025 Rsize[dims] = [axes(r1)..., ntuple(d->OneTo(1), nextra)...]
2026 end
2027 R = similar(r1, tuple(Rsize...,))
2028
2029 ridx = Any[map(first, axes(R))...]
2030 for d in dims
2031 ridx[d] = axes(R,d)
2032 end
2033
2034 concatenate_setindex!(R, r1, ridx...)
2035
2036 nidx = length(otherdims)
2037 indices = Iterators.drop(CartesianIndices(itershape), 1) # skip the first element, we already handled it
2038 inner_mapslices!(safe_for_reuse, indices, nidx, idx, otherdims, ridx, Aslice, A, f, R)
2039 end
2040
2041 @noinline function inner_mapslices!(safe_for_reuse, indices, nidx, idx, otherdims, ridx, Aslice, A, f, R)
2042 if safe_for_reuse
2043 # when f returns an array, R[ridx...] = f(Aslice) line copies elements,
2044 # so we can reuse Aslice
2045 for I in indices
2046 replace_tuples!(nidx, idx, ridx, otherdims, I)
2047 _unsafe_getindex!(Aslice, A, idx...)
2048 concatenate_setindex!(R, f(Aslice), ridx...)
2049 end
2050 else
2051 # we can't guarantee safety (#18524), so allocate new storage for each slice
2052 for I in indices
2053 replace_tuples!(nidx, idx, ridx, otherdims, I)
2054 concatenate_setindex!(R, f(A[idx...]), ridx...)
2055 end
2056 end
2057
2058 return R
2059 end
2060
2061 function replace_tuples!(nidx, idx, ridx, otherdims, I)
2062 for i in 1:nidx
2063 idx[otherdims[i]] = ridx[otherdims[i]] = I.I[i]
2064 end
2065 end
2066
2067 concatenate_setindex!(R, v, I...) = (R[I...] .= (v,); R)
2068 concatenate_setindex!(R, X::AbstractArray, I...) = (R[I...] = X)
2069
2070 ## 1 argument
2071
2072 function map!(f::F, dest::AbstractArray, A::AbstractArray) where F
2073 for (i,j) in zip(eachindex(dest),eachindex(A))
2074 val = f(@inbounds A[j])
2075 @inbounds dest[i] = val
2076 end
2077 return dest
2078 end
2079
2080 # map on collections
2081 map(f, A::AbstractArray) = collect_similar(A, Generator(f,A))
2082
2083 # default to returning an Array for `map` on general iterators
2084 """
2085 map(f, c...) -> collection
2086
2087 Transform collection `c` by applying `f` to each element. For multiple collection arguments,
2088 apply `f` elementwise.
2089
2090 See also: [`mapslices`](@ref)
2091
2092 # Examples
2093 ```jldoctest
2094 julia> map(x -> x * 2, [1, 2, 3])
2095 3-element Array{Int64,1}:
2096 2
2097 4
2098 6
2099
2100 julia> map(+, [1, 2, 3], [10, 20, 30])
2101 3-element Array{Int64,1}:
2102 11
2103 22
2104 33
2105 ```
2106 """
2107 map(f, A) = collect(Generator(f,A))
2108
2109 map(f, ::AbstractDict) = error("map is not defined on dictionaries")
2110 map(f, ::AbstractSet) = error("map is not defined on sets")
2111
2112 ## 2 argument
2113 function map!(f::F, dest::AbstractArray, A::AbstractArray, B::AbstractArray) where F
2114 for (i, j, k) in zip(eachindex(dest), eachindex(A), eachindex(B))
2115 @inbounds a, b = A[j], B[k]
2116 val = f(a, b)
2117 @inbounds dest[i] = val
2118 end
2119 return dest
2120 end
2121
2122 ## N argument
2123
2124 @inline ith_all(i, ::Tuple{}) = ()
2125 function ith_all(i, as)
2126 @_propagate_inbounds_meta
2127 return (as[1][i], ith_all(i, tail(as))...)
2128 end
2129
2130 function map_n!(f::F, dest::AbstractArray, As) where F
2131 idxs1 = LinearIndices(As[1])
2132 @boundscheck LinearIndices(dest) == idxs1 && all(x -> LinearIndices(x) == idxs1, As)
2133 for i = idxs1
2134 @inbounds I = ith_all(i, As)
2135 val = f(I...)
2136 @inbounds dest[i] = val
2137 end
2138 return dest
2139 end
2140
2141 """
2142 map!(function, destination, collection...)
2143
2144 Like [`map`](@ref), but stores the result in `destination` rather than a new
2145 collection. `destination` must be at least as large as the first collection.
2146
2147 # Examples
2148 ```jldoctest
2149 julia> a = zeros(3);
2150
2151 julia> map!(x -> x * 2, a, [1, 2, 3]);
2152
2153 julia> a
2154 3-element Array{Float64,1}:
2155 2.0
2156 4.0
2157 6.0
2158 ```
2159 """
2160 map!(f::F, dest::AbstractArray, As::AbstractArray...) where {F} = map_n!(f, dest, As)
2161
2162 map(f) = f()
2163 map(f, iters...) = collect(Generator(f, iters...))
2164
2165 # multi-item push!, pushfirst! (built on top of type-specific 1-item version)
2166 # (note: must not cause a dispatch loop when 1-item case is not defined)
2167 push!(A, a, b) = push!(push!(A, a), b)
2168 push!(A, a, b, c...) = push!(push!(A, a, b), c...)
2169 pushfirst!(A, a, b) = pushfirst!(pushfirst!(A, b), a)
2170 pushfirst!(A, a, b, c...) = pushfirst!(pushfirst!(A, c...), a, b)
2171
2172 ## hashing AbstractArray ##
2173
2174 function hash(A::AbstractArray, h::UInt)
2175 h = hash(AbstractArray, h)
2176 # Axes are themselves AbstractArrays, so hashing them directly would stack overflow
2177 # Instead hash the tuple of firsts and lasts along each dimension
2178 h = hash(map(first, axes(A)), h)
2179 h = hash(map(last, axes(A)), h)
2180 isempty(A) && return h
2181
2182 # Goal: Hash approximately log(N) entries with a higher density of hashed elements
2183 # weighted towards the end and special consideration for repeated values. Colliding
2184 # hashes will often subsequently be compared by equality -- and equality between arrays
2185 # works elementwise forwards and is short-circuiting. This means that a collision
2186 # between arrays that differ by elements at the beginning is cheaper than one where the
2187 # difference is towards the end. Furthermore, blindly choosing log(N) entries from a
2188 # sparse array will likely only choose the same element repeatedly (zero in this case).
2189
2190 # To achieve this, we work backwards, starting by hashing the last element of the
2191 # array. After hashing each element, we skip `fibskip` elements, where `fibskip`
2192 # is pulled from the Fibonacci sequence -- Fibonacci was chosen as a simple
2193 # ~O(log(N)) algorithm that ensures we don't hit a common divisor of a dimension
2194 # and only end up hashing one slice of the array (as might happen with powers of
2195 # two). Finally, we find the next distinct value from the one we just hashed.
2196
2197 # This is a little tricky since skipping an integer number of values inherently works
2198 # with linear indices, but `findprev` uses `keys`. Hoist out the conversion "maps":
2199 ks = keys(A)
2200 key_to_linear = LinearIndices(ks) # Index into this map to compute the linear index
2201 linear_to_key = vec(ks) # And vice-versa
2202
2203 # Start at the last index
2204 keyidx = last(ks)
2205 linidx = key_to_linear[keyidx]
2206 fibskip = prevfibskip = oneunit(linidx)
2207 n = 0
2208 while true
2209 n += 1
2210 # Hash the current key-index and its element
2211 elt = A[keyidx]
2212 h = hash(keyidx=>elt, h)
2213
2214 # Skip backwards a Fibonacci number of indices -- this is a linear index operation
2215 linidx = key_to_linear[keyidx]
2216 linidx <= fibskip && break
2217 linidx -= fibskip
2218 keyidx = linear_to_key[linidx]
2219
2220 # Only increase the Fibonacci skip once every N iterations. This was chosen
2221 # to be big enough that all elements of small arrays get hashed while
2222 # obscenely large arrays are still tractable. With a choice of N=4096, an
2223 # entirely-distinct 8000-element array will have ~75% of its elements hashed,
2224 # with every other element hashed in the first half of the array. At the same
2225 # time, hashing a `typemax(Int64)`-length Float64 range takes about a second.
2226 if rem(n, 4096) == 0
2227 fibskip, prevfibskip = fibskip + prevfibskip, fibskip
2228 end
2229
2230 # Find a key index with a value distinct from `elt` -- might be `keyidx` itself
2231 keyidx = findprev(!isequal(elt), A, keyidx)
2232 keyidx === nothing && break
2233 end
2234
2235 return h
2236 end