Deprecate +/- methods for array+scalar etc

[andyferris]

There has been much discussion around this lately. This PR deprecates array + scalar (and related) in favor of array .+ scalar for v0.7 (with the possibility of introducing behavior such that matrix + scalar == matrix + scalar*I in v1.0 in a separate PR).

The arguments for changing the behavior of array + scalar from elementwise are:

• Apart from the binary Number/AbstractArray methods deprecated in this PR, all the other methods of +, -, *, / involving AbstractArray are motivated by and consistent with linear algebra. For instance vector + vector is defined because vector belongs to a vector space (not because vector .+ vector is inconvenient) and a such, the entire container adds (not just the elements). However, the current element-wise definition seems to be incorrect for linear algebra. For instance, if I write a polynomial function such as f(a,b,c,x) = a*x^2 + b*x + c, I expect it to work for any x that belongs to a field/ring/whatever. However, this definition is broken if x is a square matrix - you expect behavior like f(a,b,c,x) = a*x^2 + b*x + c*I but currently the c is added to all parts of the matrix, not just the diagonal part. It is therefore difficult to write generic code that works for a wide range of types of x (though I would note a*x*x + b*x + c*one(x) currently works for x::Number and x::AbstractMatrix, though the extra one seems unnatural).
• It seems to violate a property of one that you might expect, whereby array + scalar != array + scalar*one(array)
• AFAIK these are the last auto-vectorized methods that haven't been deprecated as a part of the dot-call changes (e.g. sin(array) no longer works; rather sin.(array) is necessary). In this situation, I would actually argue that array .+ scalar for elementwise addition is less confusing... it's less of a special case.

The primary arguments against (that I am aware of) are:

• The current array + scalar behavior is convenient (though .+ is hardly much harder to type that +).
• We tried this once before and it was difficult (before the dot-call convention was introduced and well-known?)

[fredrikekre]

xref #22880

[antoine-levitt]

While +(matrix, scalar) is usually more a source of bugs than anything, +(vector, scalar), while wrong in the linear algebra sense, looks pretty harmless to me, and can sometimes be useful (not every array is used as a vector in linear algebra, after all). If +(vector, scalar) is too inconvenient to remove, a middle ground could be defining it only for rank-1 arrays.

[ararslan]

Should reference the PR number here like the other bullet points do

[StefanKarpinski]

For extra credit reference the commit that the NEWS change is made in.

[stevengj]

Previous attempt was: #5807, #5810, reverted in #7226.

The main difference, as @andyferris says, is that now we've eliminated most of the implicit broadcasting in favor of dot calls, and dot calls have some performance advantages (when they are combined) so it is not merely a matter of spelling.

[pabloferz]

Could you add a note here indicating that when these deprecations get removed, we should rewrite these methods throwing a MethodError that suggests using the dot version and explaining why it fails? (I imagine this could be regular pitfall for a lot of people)

[tkelman]

I think these would need deprecating too, GeneralPeriod isn't a subtype of Number right?

[andyferris]

Good idea.

[andyferris]

done

[tkelman]

not really broken so much as not done by broadcast, since it's tough in general to know whether an arbitrarily fused operation preserves spacing of a range

[andyferris]

Right. I've been playing with this to make the test better... I just realized there was no 1 .+ r specialization either, for ranges. (In fact, there's a couple other strange definitions there I will deal with, like 3 / (1:3) still exists whereas 3 / [1,2,3] has already been removed in a previous vectorized operations purge).

[andyferris]

Oh god, I'm in some kind of hell where (a) I can't use stacktraces to figure out why my sysimg won't build and (b) we still have the v0.6 deprecations and some of my new range broadcasts are somehow being taken over by the old .+ operator or something...

[andyferris]

Yep, I can confirm that 1 .+ (1:10) does not call broadcast(+, 1, 1:10) on v0.7.

Wow.

[andyferris]

If anyone is curious, this one was a pain to find...

[JeffBezanson]

Any reason not to use .+ here?

[andyferris]

I was scared of what I think turned out to the issue of literals "disguising" the broadcasted function. I really should go through this PR and revert calls like this to .+ everywhere literals aren't used.

[pabloferz]

Are all broadcast(f, ...) for f = {*, /, \} methods really needed? Either literals or broadcast fusion can defeat the purpose of this definition returning and array iterating over each element of the range. I'd say that is why we have *, / and \ special-cased.

[andyferris]

Good point - I'm open for debate - any other opinions?

[stevengj]

I would add a deprecation for these methods in favor of the non-dot functions, to help anyone who is expecting them to produce a range.

[andyferris]

So, what do we want here? We can't deprecate methods that don't exist on v0.6, but I could just remove these definitions. But as a user I'd be a bit confused why I need .+ but can't use .*, however the dot-call lowering is pretty complex also. (The ranges stuff interacting with broadcast was clearly the worst part of doing this PR; it's kind of awkward to preserve Rangeness now).

[andreasnoack]

Triage is for this so could you finish it @andyferris?

[StefanKarpinski]

Bumbity, bumbity.

[andyferris]

Yes, sorry, real world pressures ATM... I will be finishing this when I can.

[fredrikekre]

Is there anything else than rebasing to be done?

[andyferris]

@fredrikekre I was worried I messed up a force push at some point and lost some final fixes I had done... and I haven't gotten back to this to fix it up.

[StefanKarpinski]

Let us know if there's anything else we can do to help!

[morningkyle]

@andyferris

It seems to violate a property of one that you might expect, whereby array + scalar != array + scalar*one(array)

Is array conceptually the same thing as matrix in above context and in Julia? I have this question because array in general could represent more concept than matrix from linear algebra. For example, it's common to use a 2D array as a pixel container for a 2D bitmap image. Then "array + scalar" means operation on each pixel in the container, just nothing to do with the matrix math.

Edit: For arrays that are higher than 3D(eg. m x n x k), matrix multiplication has no definition and not supported. Also there is no "one"(either Identity or UniformScaling) defined for these arrays.

Sorry to ask this question here if it is not appropriate for a non-contributor.

[andyferris]

@morningkyle Julia has a standard and dedicated way of operating on each element of a container, such as the pixels of an image, using broadcast and its short . form like array .+ scalar. This last form makes it obvious to everyone that you intend to do an elementwise operation and is basically just as easy to write.

We've generally heavily overloaded arrays with linear algebra operations like * for matrix multiplication (opposed to elementwise .*) so I felt it was confusing having * specifically for linear algebra and + specifically inconsistent with linear algebra... Other moves in this direction include exp(matrix) vs exp.(matrix) and so-on.

[morningkyle]

@andyferris Yes, deprecating "matrix + scalar" makes sense to me. I am not sure if deprecating "array + scalar" also drops other conveniences. Let me give one more code example below:
a = [1 2 3]
f(x) = 2x + 1
Now is it valid to call function f in below format? If it is valid, it seems "array + scalar" is used in this case.
f(a)
Or should I always define f(x) in this format in order to make it more general?
f(x) = 2 .* x .+ 1

Edit myself: We should use the dot-call syntax f.(a) to accept the array input. There is no inconvenience if the dot-call syntax is well accepted in above example.

[andyferris]

Edit myself: We should use the dot-call syntax f.(a) to accept the array input. There is no inconvenience if the dot-call syntax is well accepted in above example.

You're right - exactly. As a general rule of thumb, genericism can only go so far. One thing we've tended to do in Julia is separate out functions which act on "containers of objects" from functions acting on "objects" (scalars, or whatever) themselves. Think max(a,b) vs maximum(a), where the former finds the largest of some objects, and the latter does something similar but searching through the elements of a container. I recommend that even for very generic code you should generally be aware whether your bindings are to collections of things you are interested in (e.g. an array of numbers) or if it is a binding to one of the things themselves (e.g. a number). We've removed most functions that tend to act on both - except for the important case where a "container" (like a AbstractVector or AbstractMatrix) can be considered an interesting object in its own right (in this case w.r.t. linear algebra, they model vector spaces and linear operators). My experience is that this approach allows you write clean and generic code, while still retaining flexibility.

deprecating "matrix + scalar" makes sense to me. I am not sure if deprecating "array + scalar" also drops other conveniences. Let me give one more code example below:

I feel it would be really confusing if this worked differently for 2D arrays, than all the other dimensions.

[StefanKarpinski]

Or should I always define f(x) in this format in order to make it more general?

It depends on what you mean. 2x and 2.*x mean the same thing, so that's not an issue, but 2x .+ 1 is not the only possible matrix analogue of 2x + 1 for scalars. Sometimes it may be, but often a more valid analogue is 2x + I, which adds ones only to the diagonal of 2x.

[morningkyle]

Array: a continuous block of data of the same element type. With element type be scalar, array supports arbitrary arithmetic on each single element, thus a very general and fundamental data structure.

Matrix: linear algebra concept, an m × n array of scalars from a given field F. So matrix is just a subset form of array and supports restricted mathematical operations. That's the main point I think deprecating "matrix + scalar" makes sense.

Especially, arrays that are equal or higher than 3D(m x n x k) are not like matrix at all(eg. matrix like multiplication or * operator has no definition for 3D array).

Now let's try to define "array + scalar".
array + scalar: modify all elements of array with the same operation(plus by scalar). Another way to do this is array .+ scalar. And the later expression is more explicit.

Now come back to these points:

but 2x .+ 1 is not the only possible matrix analogue of 2x + 1 for scalars. Sometimes it may be, but often a more valid analogue is 2x + I

deprecating "matrix + scalar" makes sense to me. I am not sure if deprecating "array + scalar" also drops other conveniences.

If x is a matrix, then we could either make "2x + 1" an invalid operation or make it a convenient form for of 2x + I. The current choice is the former. And it makes sense. I have no problem with this choice.

But if x is an array, based on above "array + scalar" definition, 2x + 1 is well defined and convenient. 2x .+ 1 is also well defined and not inconvenient. But why should we stop people using the 2x + 1 format? If x is an array, we should not assume it as a matrix unless in a predefined context that we know we are working with linear algebra. Without such specific context, it's easy and natural to just think x as a block of data. Because array is a common and widely used concept in other programing languages(C/C++/C#/Java/Python) before Julia. Taking "array + scalar" as a convenient form of block modification, then keeping the form is not confusing and makes sense to me.

Python array(from numpy) supports "array + scalar" and Matlab matrix supports "matrix + scalar". I don't argue that Julia should do the same, just meant I did not invent the "array + scalar" concept.

[oscardssmith]

I think that the problem is that matrix is a type of array, so defining array+scaler in a way that doesn't make sense for matrix will lead to a definition of matrix+ array that doesn't make sense

[TotalVerb]

@morningkyle The reason it is considered confusing to support vector + scalar is that this is not consistent with how other operations work. It is better for a user migrating from Python or Matlab to realize that things like sin(A) are not "pointwise" operations. If we allow some operations to be implicitly pointwise, like A + x, this can cause confusion.

Because as you mention the convention A .+ x is not inconvenient, this kind of "explicit pointwiseness" using . is preferred in all situations in Julia.

[andyferris]

@morningkyle I feel that you are using definitions here which aren't true in Julia. In particular, given:

With element type be scalar, array supports arbitrary arithmetic on each single element, thus a very general and fundamental data structure.

I think that you are taking our array to be similar to those in MATLAB, numpy, etc, where f(array) and op(scalar, array) generally does do some elementwise computation, which is generally not true for Julia. There are two important special cases which are introduced because of linear algebra; these are *(scalar, array) and +(array, array) (i.e., we say that arrays live in a vector space, just like vectors, matrices and tensors in math). But apart from "vector space" operations (here I also include things like -, /, real(array), etc) there are (after this PR, and as far as I'm aware) precisely zero operations which do arithmetic on each element. Even when @StefanKarpinski says 2x and 2.*x mean the same thing, this basically assumes that x is an element of a vector space and (if x is a collection) that getindex (or iteration) returns basis element coefficients (true for Number and Array).

[morningkyle]

@andyferris

@morningkyle I feel that you are using definitions here which aren't true in Julia.
With element type be scalar, array supports arbitrary arithmetic on each single element, thus a very general and fundamental data structure.

It's true in following sense(code). Julia array itself is still very similar to those in other languages. See the first statement in https://docs.julialang.org/en/stable/manual/arrays/. It just that Julia built-in operators/functions do not simulate other languages(eg. Python and Matlab). Anyway, I got your point(built-in arithmetic functions generally do not accept array input without dot-call syntax) and thanks very much!

a = [2, 1, 3]
a_new = [a[1]/2, a[2]+1, sin(a[3])]  #arbitrary arithmetic on each element
a_sorted = sort!(a)  #sort accepts array and updates elements

[StefanKarpinski]

It's not mentioned here, but the reason we need to get rid of A + 1 is somewhat involved:

1. We want to allow x + I to work when x is a square matrix or a scalar – this allows e.g. writing generic polynomials like x^2 + 2x + I that work on both scalars and square matrices.

2. Matrix + scalar and scalar + uniform scaling together badly violate associativity of +:

julia> (zeros(3,3) + 1) + I
3×3 Array{Float64,2}:
 2.0  1.0  1.0
 1.0  2.0  1.0
 1.0  1.0  2.0

julia> zeros(3,3) + (1 + I)
3×3 Array{Float64,2}:
 2.0  2.0  2.0
 2.0  2.0  2.0
 2.0  2.0  2.0

If the first + has to be written as .+ then there's no issue since .+ and + are not the same operator and not required or expected to be associative. We could potentially make zeros(3,3) + 1 do what zeros(3,3) + I currently does and no violate associativity. That would allow writing generic polynomials (and other arithmetic expressions of matrices or scalar) more naturally as x^2 + 2x + 1, but I suspect that behavior might be a little too surprising compared to other languages.

[andyferris]

Similarly, it would great if the distributive law worked, for example A*v - λ*v == (A - λ)*v for the Eigenvalue problem (matrix A, vector v, scalar λ).

[morningkyle]

@StefanKarpinski

(zeros(3,3) + 1) + I

Both "matrix + scalar" and "array + I(identity)" are odd operations without brain type conversion first.

Now your example mixed the two:
A + scalar + I
I wouldn't be surprised that if anything like associativity breaks here, especially if you presume A is a matrix(matrix operation has strict rules). So, I feel it is not a fair example:)

BTW, zeros(m, n) returns an literal array for all documentations from Julia/Python/Matlab.
Julia: https://docs.julialang.org/en/latest/stdlib/arrays/
Python: https://docs.scipy.org/doc/numpy/reference/generated/numpy.zeros.html
Matlab: https://www.mathworks.com/help/matlab/matrices-and-arrays.html?s_tid=gn_loc_drop

As a conclusion to myself after these discussion, it is subtle to replace "array + scalar" with "array .+ scalar", but I am getting used to it.

[StefanKarpinski]

I wouldn't be surprised that if anything like associativity breaks here, especially if you presume A is a matrix (matrix operation has strict rules). So, I feel it is not a fair example:)

The point is that when writing generic code, you don't know what specific implementations of + are going to be called. That's why it's crucial that + always obey a semantically strict meaning and that general laws like associativity be followed, regardless of what kinds of arguments are involved.

[morningkyle]

That's why it's crucial that + always obey a semantically strict meaning and that general laws like associativity be followed, regardless of what kinds of arguments are involved.

This is a good point. What I think is that your example is confusing. It does not support your (this agreed) point very well, then does not support the logic to "get rid of A + 1" very well.

A + scalar + I

(1) If I is an identity matrix here, then associativity is satisfied.
(2) If I is an UniformScaling, associativity is broken. But it is mainly because (matrix + UniformScaling) and (scalar + UniformScaling) are broken definitions first. I guess both of them should use the multiplication syntax.

julia> 2 + I

WARNING: x::Number + J::UniformScaling is deprecated, use x + J.λ instead.
Stacktrace:
 [1] depwarn at ./deprecated.jl:70 [inlined]
 [2] +(::Int64, ::UniformScaling{Int64}) at ./deprecated.jl:57
 [3] include_string(::Module, ::String, ::String) at ./loading.jl:431
 [4] execute_request(::ZMQ.Socket, ::IJulia.Msg) at /opt/julia_packages/.julia/v0.7/IJulia/src/execute_request.jl:154
 [5] (::Base.#inner#4{Array{Any,1},IJulia.#execute_request,Tuple{ZMQ.Socket,IJulia.Msg}})() at ./essentials.jl:421
 [6] eventloop(::ZMQ.Socket) at /opt/julia_packages/.julia/v0.7/IJulia/src/eventloop.jl:8
 [7] (::IJulia.##14#17)() at ./task.jl:335
while loading In[14], in expression starting on line 1

3

julia> 2 .+ I

WARNING: x::Number + J::UniformScaling is deprecated, use x + J.λ instead.
Stacktrace:
 [1] depwarn at ./deprecated.jl:70 [inlined]
 [2] +(::Int64, ::UniformScaling{Int64}) at ./deprecated.jl:57
 [3] (::##7#8)(::UniformScaling{Int64}) at ./<missing>:0
 [4] broadcast(::Function, ::UniformScaling{Int64}) at ./broadcast.jl:434
 [5] include_string(::Module, ::String, ::String) at ./loading.jl:431
 [6] execute_request(::ZMQ.Socket, ::IJulia.Msg) at /opt/julia_packages/.julia/v0.7/IJulia/src/execute_request.jl:154
 [7] (::Base.#inner#4{Array{Any,1},IJulia.#execute_request,Tuple{ZMQ.Socket,IJulia.Msg}})() at ./essentials.jl:421
 [8] eventloop(::ZMQ.Socket) at /opt/julia_packages/.julia/v0.7/IJulia/src/eventloop.jl:8
 [9] (::IJulia.##14#17)() at ./task.jl:335
while loading In[15], in expression starting on line 1

3

[StefanKarpinski]

If I is an UniformScaling, associativity is broken.

Yes, that's why we're changing it.

[timholy]

The arguments behind this are very convincing, and I'm a supporter. But I just noticed a little issue that would be hilarious if it weren't so...well, you'll see. Watch the julia version numbers and return types:

  | | |_| | | | (_| |  |  Version 0.5.2 (2017-05-06 16:34 UTC)
 _/ |\__'_|_|_|\__'_|  |  
|__/                   |  x86_64-linux-gnu

julia> (1:5) + 1
2:6

julia> (1:5) .+ 1
2:6

julia> broadcast(+, 1:5, 1)
5-element Array{Int64,1}:
 2
 3
 4
 5
 6

\  | | |_| | | | (_| |  |  Version 0.6.1-pre.0 (2017-06-19 13:06 UTC)
 _/ |\__'_|_|_|\__'_|  |  Commit dcf39a1 (92 days old release-0.6)
|__/                   |  x86_64-linux-gnu

julia> (1:5) + 1
2:6

julia> (1:5) .+ 1
5-element Array{Int64,1}:
 2
 3
 4
 5
 6

julia> broadcast(+, 1:5, 1)
5-element Array{Int64,1}:
 2
 3
 4
 5
 6

  | | |_| | | | (_| |  |  Version 0.7.0-DEV.1846 (2017-09-19 12:18 UTC)
 _/ |\__'_|_|_|\__'_|  |  Commit 77cbd18* (0 days old master)
|__/                   |  x86_64-linux-gnu

julia> (1:5) + 1
WARNING: a::AbstractArray + b::Number is deprecated, use broadcast(+, a, b) instead.
Stacktrace:
 [1] depwarn(::String, ::Symbol) at ./deprecated.jl:68
 [2] +(::UnitRange{Int64}, ::Int64) at ./deprecated.jl:56
 [3] eval(::Module, ::Expr) at ./repl/REPL.jl:3
 [4] eval_user_input(::Any, ::Base.REPL.REPLBackend) at ./repl/REPL.jl:69
 [5] macro expansion at /home/tim/.julia/v0.7/Revise/src/Revise.jl:748 [inlined]
 [6] (::getfield(Revise, Symbol("##11#12")){Base.REPL.REPLBackend})() at ./event.jl:96
while loading no file, in expression starting on line 0
2:6

julia> (1:5) .+ 1
5-element Array{Int64,1}:
 2
 3
 4
 5
 6

julia> broadcast(+, 1:5, 1)
2:6

Ha ha. Ugh.

[andyferris]

Haha. Yeah... sorry for contributing to that :)

Making.broadcast fusion more "transparent" to code may solve this neatly.

[StefanKarpinski]

I don't think that's any sort of show-stopper, I think that's really just a series of accidental behaviors flip-flopping. Correct behaviors IMO:

(1:5) + 1            => error
(1:5) .+ 1           => 2:5
broadcast(+, 1:5, 1) => 2:5

Fortunately, we're quite close to that already and just need to fix broadcasting over ranges. Even if we don't it's not really all that significant since 2:5 is just an efficient representation of [2:5;].

[timholy]

It's definitely not a show-stopper. It does make it very likely that some packages (especially those that define new types of AbstractArrays) will have a fairly substantial pile of if VERSION statements to support basic indexing operations. Still, the obvious thing to do is to get (1:5) .+ 1 returning a UnitRange, and all will be better than it's ever been before.

EDIT: Probably best would be to have a Compat.broadcast that exhibits the 0.7 behavior.

[martinholters]

This removed /(x::Number, r::Range) without going through a deprecation, ref #23067 (comment).

[andyferris]

Hmm. Maybe I thought that method was a bug, not a feature?

A deprecation seems fine.

[andreasnoack]

Maybe I thought that method was a bug, not a feature?

It was a bug. It must have been hiding the range code when the method for vectors was removed. Still fine to deprecate, though.