Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

linspace: try to "lift" linspace the float ranges are [close #9637] #9666

Merged
merged 13 commits into from
Apr 13, 2015

Conversation

StefanKarpinski
Copy link
Sponsor Member

No description provided.

@JeffBezanson
Copy link
Sponsor Member

Is there any way to have [linrange(...)] do this, so we don't need as many functions?

@StefanKarpinski
Copy link
Sponsor Member Author

Not quite – the lifted cases are the same, but the unlifted cases are different: ranges use start + k*step whereas linspace uses start*(n-k)/n + stop*k/n.

@StefanKarpinski
Copy link
Sponsor Member Author

This now seems to be uniformly better than what we had before – the tests are pretty comprehensive. There are still some corner cases like this:

julia> linspace(realmin(),realmax(),maxintfloat())
linspace(0.0,1.7976931348623157e308,9007199254740992)

Note that the start of the resulting LinSpace object is exactly zero, rather than realmin(). I don't actually think that can be fixed, however, without adding significant complexity to the implementation. I'm not sure if we should throw an error upon construction of such a linear space object or not.

@tkelman
Copy link
Contributor

tkelman commented Apr 3, 2015

looks like there's a 32 bit failure

@StefanKarpinski StefanKarpinski force-pushed the sk/linspace branch 2 times, most recently from aa30ada to 03c2607 Compare April 13, 2015 16:34
The length zero and one cases for LinSpace are really weird – arguably
they simply shouldn't be allowed at all. However, I think for usability
we really need to support these. This change uses -1 and 0 as special
values of the r.len field, encoding, somewhat counterintuitively, real
lengths of 0 and 1, respectively. The reason for this is so that the
element access computations for the length one case has r.len-1 == -1
and only flips sign, avoiding overflow when r.start is very large. In
the length zero case, r.len-1 == -2, so we need to divide by -2 to get
correct first and last values back. This is susceptible to overflow,
but since these values are kind of nonsense in the zero case, this
seems to be a preferrable situation.

There are still cases where the (r.len-i)*r.start + (i-1)*r.stop can
overflow for longer ranges, and it's unclear to me whether you can
handle those correctly with the same LinSpace type. This plus the
complexity of handling corner cases like lengths zero and one, makes
me wonder if it's really a good idea to make this a type at all.
This forces the LinSpace divisor – either n*e or just n – to be in
the range [1/2,0), which, in particular, prevents overflow even if
the old computation could overflow. Needs more testing of corner
cases, but this passes all of the existing tests at least.
E.g. linspace(realmin(), realmax(), maxintfloat())
The issue was that doing `one(p) << p` where p is Int produces zero
on a 32-bit system if p is larger than 32, which it can be for some
Float64 arguments.

Getting tests to pass for 32-bit also required avoiding taking the
length of LinSpace objects whose lengths are 2^31 or longer. This is
a bit of a more general problem since we can represent FloatRange
and LinSpace objects whose length cannot be represented as an Int on
both 64- and 32-bit platforms, although it is, of course, easier to
encounter the problem on 32-bit systems. In fact, on 64-bit systems
it is not actually possible to construct LinSpace{Float64} objects
whose length cannot be represented as an Int. The issue does exist
for FloatRange on 64-bit systems, however. One possible solution is
to introduce length(T,x) that gives lengths in a particular type and
just define length(x) as length(Int,x) as the default behavior.
StefanKarpinski added a commit that referenced this pull request Apr 13, 2015
linspace: try to "lift" linspace the float ranges are [close #9637]
@StefanKarpinski StefanKarpinski merged commit 22a73f9 into master Apr 13, 2015
@timholy
Copy link
Sponsor Member

timholy commented Apr 16, 2015

What should the policy be for these?

julia> [0.0,0.5,1.0] == linspace(0.0, 1.0, 3)
false

julia> convert(Vector, linspace(0.0, 1.0, 3))
ERROR: MethodError: `convert` has no method matching convert(::Type{Array{T,1}}, ::LinSpace{Float64})
This may have arisen from a call to the constructor Array{T,1}(...),
since type constructors fall back to convert methods.
Closest candidates are:
  convert{T}(::Type{T}, ::T)
  convert{T<:Integer}(::Type{T<:Integer}, ::Rational{T<:Integer})
  convert{T<:Integer}(::Type{T<:Integer}, ::Float16)
  ...

@timholy
Copy link
Sponsor Member

timholy commented Apr 17, 2015

Bump @StefanKarpinski

@StefanKarpinski
Copy link
Sponsor Member Author

So we have the same issue with FloatRanges:

julia> [0.0,0.5,1.0] == 0:0.5:1
false

julia> convert(Vector, 0:0.5:1)
ERROR: MethodError: `convert` has no method matching convert(::Type{Array{T,1}}, ::FloatRange{Float64})
This may have arisen from a call to the constructor Array{T,1}(...),
since type constructors fall back to convert methods.
Closest candidates are:
  convert{T}(::Type{T}, ::T)
  convert{T<:Tuple{Any,Vararg{Any}}}(::Type{T<:Tuple{Any,Vararg{Any}}}, ::Tuple{Any,Vararg{Any}})
  convert{T<:Integer}(::Type{T<:Integer}, ::Rational{T<:Integer})
  ...

@StefanKarpinski
Copy link
Sponsor Member Author

Which doesn't really answer the question, but at least it's consistent. The same answer should apply.

@timholy
Copy link
Sponsor Member

timholy commented Apr 17, 2015

convert(Vector{Float64}, obj) works but not convert(Vector, obj)---is this deliberate? Can we add convert{T}(::Type{Vector}, r::LinSpace{T}) = convert(Vector{T}, r)?

@BobPortmann
Copy link
Contributor

Shouldn't this have become a new version of linrange (and LinRange) instead of linspace (and LinSpace) since it returns a range object. Then linspace would be removed?

@StefanKarpinski
Copy link
Sponsor Member Author

Linspace is the traditional name for this kind of functionality.

@BobPortmann
Copy link
Contributor

OK, I'll take your word for it. I've never run across the name in many years using Fortran and IDL and a bit of Mathematica and R. Julia was the first place I saw it.

@tkelman tkelman deleted the sk/linspace branch April 18, 2015 01:36
@StefanKarpinski
Copy link
Sponsor Member Author

It's from Matlab

@StefanKarpinski
Copy link
Sponsor Member Author

I'd be fine with renaming it to linrange, but I'd like some input from others before doing that. Any opinions?

@nalimilan
Copy link
Member

Since you made it return a range-like object, it makes sense to me to rename it to linrange. (Coming from R, I admit linspace was not familiar to me either -- but nor was linrange of course :-).

@sebastiang
Copy link

So now linspace returns a cleverer object, but if people are expecting it to be an array as once it was, what should they do now? I'm running into this trying to fix Gadfly issues with v0.4

@StefanKarpinski
Copy link
Sponsor Member Author

An obvious solution is to do collect(linspace(start, stop, len)). What's the specific issue you're encountering?

@sebastiang
Copy link

A call to distinguishable_colors in Gadfly assumes linspace will match the ::Vector constraint required, as do apparently the default values provided in the library. To be honest it looks like the constraints could be lifted from distinguishable_colors itself and the code would work correctly, but you hate to dig too deep.

function distinguishable_colors{T<:ColorValue}(n::Integer,
                        seed::Vector{T};
                        transform::Function = identity,
                        lchoices::Vector = linspace(0, 100, 15),
                        cchoices::Vector = linspace(0, 100, 15),
                        hchoices::Vector = linspace(0, 340, 20))

@mbauman
Copy link
Sponsor Member

mbauman commented Apr 22, 2015

I'd just change it to be an AbstractVector, but you're right that it probably doesn't need to be typed at all.

(This is an impressive tour d'packages that you've embarked upon! Very nice work. 👍 )

@StefanKarpinski
Copy link
Sponsor Member Author

Note that if we have automatic conversion of default arguments to the type of the argument then this wouldn't be a problem at all since the LinSpace object would be automatically converted to Vector. Between this and the Nullable argument, I really think the evidence that automatic conversion is the right thing to do is piling up.

@sebastiang
Copy link

You've seen a lot more than I have, but it seems odd to have conversion for some arguments and not others. The type magic is already pretty daunting as it is.

@StefanKarpinski
Copy link
Sponsor Member Author

There's no magic, it would just mean that f(x::T=d) implicitly calls convert(T,d) for you.

@sebastiang
Copy link

Ah, not for callers, just for the guy who writes the default value into his function definition.

@quinnj
Copy link
Member

quinnj commented Apr 22, 2015

I think it makes a lot of sense as well @StefanKarpinski; particularly as we've made the push to have convert(T,x) be as conservative as possible and throw errors if something doesn't convert cleanly.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

10 participants