Re: [julia-users] Surprising range behavior
Hmm, using Julia to test the accuracy of Matlab. Interesting. On Thu, Apr 24, 2014 at 10:42 PM, Peter Simon wrote: > +1 -- Thanks! > > > On Thursday, April 24, 2014 11:33:01 AM UTC-7, Simon Byrne wrote: >> >> On Thursday, 24 April 2014 19:10:56 UTC+1, Simon Byrne wrote: >>> >>> It probably would be useful to have a linspace-ish interface for >>> constructing FloatRanges >>> >> >> To follow up >> https://github.com/JuliaLang/julia/pull/6627 >> >> >> >>
Re: [julia-users] Surprising range behavior
+1 -- Thanks! On Thursday, April 24, 2014 11:33:01 AM UTC-7, Simon Byrne wrote: > > On Thursday, 24 April 2014 19:10:56 UTC+1, Simon Byrne wrote: >> >> It probably would be useful to have a linspace-ish interface for >> constructing FloatRanges >> > > To follow up > https://github.com/JuliaLang/julia/pull/6627 > > > >
Re: [julia-users] Surprising range behavior
On Thursday, 24 April 2014 19:10:56 UTC+1, Simon Byrne wrote: > > It probably would be useful to have a linspace-ish interface for > constructing FloatRanges > To follow up https://github.com/JuliaLang/julia/pull/6627
Re: [julia-users] Surprising range behavior
Julia's computation is more accurate because we use pairwise summation instead of naive summation. This gives the same results as MATLAB: x = 0.0; for i = 1:100; x += pi/100; end; x - pi We don't use summation in Range construction, though. The reason the range ends up different is that pi/(pi/100) is very slightly less than 100, and MATLAB fudges it because the difference is less than 3 ulp. Simon On Thursday, April 24, 2014 2:05:51 PM UTC-4, Hans W Borchers wrote: > > There is a strange difference in results between Matlab and Julia that > could be the reason why [0:pi/100:pi] ends up differently: > > matlab> x = ones(100, 1) * (pi/100); > matlab> sum(x) - pi > ans = > -4.4409e-15 > > > julia> x = ones(100)*(pi/100); > julia> sum(x) - pi > 1.3322676295501878e-15 > > > So it appears that summing pi/100 100-times is smaller than pi in Matlab > and larger than pi in Julia. > > > On Thursday, April 24, 2014 5:47:34 PM UTC+2, Stefan Karpinski wrote: >> >> On Thu, Apr 24, 2014 at 11:43 AM, andrew cooke wrote: >> >>> >>> sorry, rational is a stupid suggestion for irrational numbers... >>> >>> what i am trying to say is that there is no perfect solution to this >>> problem (that i know of). if matlab works better then it is either luck or >>> because they are actually fudging things. fudging things is appealing but >>> usually bites you somewhere down the line. >>> >> >> This is exactly how I feel about it. I still think it may be possible to >> improve our FloatRange behavior even further, but I'd like to avoid >> "fudging it". >> >
Re: [julia-users] Surprising range behavior
It probably would be useful to have a linspace-ish interface for constructing FloatRanges: I wrote a quick one for KernelDensity.jl: https://github.com/JuliaStats/KernelDensity.jl/blob/master/src/univariate.jl#L60 but there's probably a better way to do this. On Thursday, 24 April 2014 19:05:51 UTC+1, Hans W Borchers wrote: > > There is a strange difference in results between Matlab and Julia that > could be the reason why [0:pi/100:pi] ends up differently: > > matlab> x = ones(100, 1) * (pi/100); > matlab> sum(x) - pi > ans = > -4.4409e-15 > > > julia> x = ones(100)*(pi/100); > julia> sum(x) - pi > 1.3322676295501878e-15 > > > So it appears that summing pi/100 100-times is smaller than pi in Matlab > and larger than pi in Julia. > > > On Thursday, April 24, 2014 5:47:34 PM UTC+2, Stefan Karpinski wrote: >> >> On Thu, Apr 24, 2014 at 11:43 AM, andrew cooke wrote: >> >>> >>> sorry, rational is a stupid suggestion for irrational numbers... >>> >>> what i am trying to say is that there is no perfect solution to this >>> problem (that i know of). if matlab works better then it is either luck or >>> because they are actually fudging things. fudging things is appealing but >>> usually bites you somewhere down the line. >>> >> >> This is exactly how I feel about it. I still think it may be possible to >> improve our FloatRange behavior even further, but I'd like to avoid >> "fudging it". >> >
Re: [julia-users] Surprising range behavior
There is a strange difference in results between Matlab and Julia that could be the reason why [0:pi/100:pi] ends up differently: matlab> x = ones(100, 1) * (pi/100); matlab> sum(x) - pi ans = -4.4409e-15 julia> x = ones(100)*(pi/100); julia> sum(x) - pi 1.3322676295501878e-15 So it appears that summing pi/100 100-times is smaller than pi in Matlab and larger than pi in Julia. On Thursday, April 24, 2014 5:47:34 PM UTC+2, Stefan Karpinski wrote: > > On Thu, Apr 24, 2014 at 11:43 AM, andrew cooke > > wrote: > >> >> sorry, rational is a stupid suggestion for irrational numbers... >> >> what i am trying to say is that there is no perfect solution to this >> problem (that i know of). if matlab works better then it is either luck or >> because they are actually fudging things. fudging things is appealing but >> usually bites you somewhere down the line. >> > > This is exactly how I feel about it. I still think it may be possible to > improve our FloatRange behavior even further, but I'd like to avoid > "fudging it". >
Re: [julia-users] Surprising range behavior
On Thu, Apr 24, 2014 at 11:43 AM, andrew cooke wrote: > > sorry, rational is a stupid suggestion for irrational numbers... > > what i am trying to say is that there is no perfect solution to this > problem (that i know of). if matlab works better then it is either luck or > because they are actually fudging things. fudging things is appealing but > usually bites you somewhere down the line. > This is exactly how I feel about it. I still think it may be possible to improve our FloatRange behavior even further, but I'd like to avoid "fudging it".
Re: [julia-users] Surprising range behavior
sorry, rational is a stupid suggestion for irrational numbers...
what i am trying to say is that there is no perfect solution to this
problem (that i know of). if matlab works better then it is either luck or
because they are actually fudging things. fudging things is appealing but
usually bites you somewhere down the line.
andrew
On Thursday, 24 April 2014 12:41:28 UTC-3, andrew cooke wrote:
>
>
> interesting, so what is the type of pi/100 in matlab? is it a rational or
> something? or do they do some kind of hacky approximation in their
> (in)equality tests?
>
> andrew
>
>
>
> On Thursday, 24 April 2014 10:28:16 UTC-3, Peter Simon wrote:
>>
>> *Matlab behavior:*
>>
>> >> format long
>> >> x = 0:pi/100:pi;
>> >> length(x)
>>
>> ans =
>>
>>101
>>
>> >> x(1) - 0
>>
>> ans =
>>
>> 0
>>
>> >> x(end) - pi
>>
>> ans =
>>
>> 0
>>
>> >> diff([max(diff(x)), min(diff(x))])
>>
>> ans =
>>
>> -4.440892098500626e-16
>>
>>
>> *Julia behavior:*
>>
>> julia> x = 0:pi/100:pi;
>>
>> julia> length(x)
>> 100
>>
>> julia> x[1]-0
>> 0.0
>>
>> julia> x[end]-pi
>> -0.031415926535897754
>>
>> julia> diff([maximum(diff(x)), minimum(diff(x))])
>> 1-element Array{Float64,1}:
>> -4.44089e-16
>>
>>
>> I highlighted the important differences in red. IMO the Matlab behavior
>> is more intuitive. If you choose an increment that "very nearly"
>> (apparently within a few ulp, as stated by Stefan) evenly divides the
>> difference between the first and last entries in the colon notation, Matlab
>> assumes that you want the resulting vector to begin and end with exactly
>> these values. This makes sense in most situations, IMO. E.g. when using a
>> range to specify sample points in a numerical integration scheme, or when
>> using a range to specify x-coordinates in an interpolation scheme. In both
>> of these cases you don't want the last point to be adjusted away from what
>> you specified using the colon notation.
>>
>> --Pete
>>
>> --Peter
>>
>> On Thursday, April 24, 2014 5:28:19 AM UTC-7, andrew cooke wrote:
>>>
>>>
>>> how does matlab differ from julia here? (i though the original problem
>>> was related to "fundamental" issues of how you represent numbers on a
>>> computer, not language specific).
>>>
>>>
>>> On Thursday, 24 April 2014 02:24:59 UTC-3, Peter Simon wrote:
Thanks, Simon, that construct works nicely to solve the problem I
posed.
I have to say, though, that I find Matlab's colon range behavior more
intuitive and generally useful, even if it isn't as "exact" as Julia's.
--Peter
On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote:
>
> pi*(0:0.01:1) or similar should work.
>
> On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote:
>>
>> Thanks for the explanation--it makes sense now. This question arose
>> for me because of the example presented in
>> https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ:
>>
>> x = [0:0.01:pi]
>>
>> used as the set of x-coordinates for sampling a function to be
>> integrated (ideally over the interval (0,pi)). But the range defined in
>> x
>> has a last entry of 3.14, which will contribute to the inaccuracy of the
>> integral being sought in that example. As an exercise, I was trying to
>> implement the interpolation solution described later in that thread by
>> Cameron McBride: "BTW, another possibility is to use a spline
>> interpolation on the original data and integrate the spline evaluation
>> with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200)
>> for
>> the x values, because CoordInterpGrid will not accept an array as its
>> first
>> argument, so you have to use a range object. But the range object has a
>> built-in error for the last point because of the present issue. Any
>> suggestions?
>>
>> Thanks,
>>
>> --Peter
>>
>> On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote:
>>>
>>> The issue is that float(pi) < 100*(pi/100). The fact that pi is not
>>> rational – or rather, that float64(pi) cannot be expressed as the
>>> division
>>> of two 24-bit integers as a 64-bit float – prevents rational lifting of
>>> the
>>> range from kicking in. I worried about this kind of issue when I was
>>> working on FloatRanges, but I'm not sure what you can really do, aside
>>> from
>>> hacks where you just decide that things are "close enough" based on
>>> some ad
>>> hoc notion of close enough (Matlab uses 3 ulps). For example, you can't
>>> notice that pi/(pi/100) is an integer – because it isn't:
>>>
>>> julia> pi/(pi/100)
>>> 99.99
>>>
>>>
>>> One approach is to try to find a real value x such that
>>> float64(x/100) == float64(pi)/100 and float64(x) == float64(pi). I
Re: [julia-users] Surprising range behavior
interesting, so what is the type of pi/100 in matlab? is it a rational or
something? or do they do some kind of hacky approximation in their
(in)equality tests?
andrew
On Thursday, 24 April 2014 10:28:16 UTC-3, Peter Simon wrote:
>
> *Matlab behavior:*
>
> >> format long
> >> x = 0:pi/100:pi;
> >> length(x)
>
> ans =
>
>101
>
> >> x(1) - 0
>
> ans =
>
> 0
>
> >> x(end) - pi
>
> ans =
>
> 0
>
> >> diff([max(diff(x)), min(diff(x))])
>
> ans =
>
> -4.440892098500626e-16
>
>
> *Julia behavior:*
>
> julia> x = 0:pi/100:pi;
>
> julia> length(x)
> 100
>
> julia> x[1]-0
> 0.0
>
> julia> x[end]-pi
> -0.031415926535897754
>
> julia> diff([maximum(diff(x)), minimum(diff(x))])
> 1-element Array{Float64,1}:
> -4.44089e-16
>
>
> I highlighted the important differences in red. IMO the Matlab behavior
> is more intuitive. If you choose an increment that "very nearly"
> (apparently within a few ulp, as stated by Stefan) evenly divides the
> difference between the first and last entries in the colon notation, Matlab
> assumes that you want the resulting vector to begin and end with exactly
> these values. This makes sense in most situations, IMO. E.g. when using a
> range to specify sample points in a numerical integration scheme, or when
> using a range to specify x-coordinates in an interpolation scheme. In both
> of these cases you don't want the last point to be adjusted away from what
> you specified using the colon notation.
>
> --Pete
>
> --Peter
>
> On Thursday, April 24, 2014 5:28:19 AM UTC-7, andrew cooke wrote:
>>
>>
>> how does matlab differ from julia here? (i though the original problem
>> was related to "fundamental" issues of how you represent numbers on a
>> computer, not language specific).
>>
>>
>> On Thursday, 24 April 2014 02:24:59 UTC-3, Peter Simon wrote:
>>>
>>> Thanks, Simon, that construct works nicely to solve the problem I posed.
>>>
>>>
>>> I have to say, though, that I find Matlab's colon range behavior more
>>> intuitive and generally useful, even if it isn't as "exact" as Julia's.
>>>
>>> --Peter
>>>
>>> On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote:
pi*(0:0.01:1) or similar should work.
On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote:
>
> Thanks for the explanation--it makes sense now. This question arose
> for me because of the example presented in
> https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ :
>
> x = [0:0.01:pi]
>
> used as the set of x-coordinates for sampling a function to be
> integrated (ideally over the interval (0,pi)). But the range defined in
> x
> has a last entry of 3.14, which will contribute to the inaccuracy of the
> integral being sought in that example. As an exercise, I was trying to
> implement the interpolation solution described later in that thread by
> Cameron McBride: "BTW, another possibility is to use a spline
> interpolation on the original data and integrate the spline evaluation
> with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200)
> for
> the x values, because CoordInterpGrid will not accept an array as its
> first
> argument, so you have to use a range object. But the range object has a
> built-in error for the last point because of the present issue. Any
> suggestions?
>
> Thanks,
>
> --Peter
>
> On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote:
>>
>> The issue is that float(pi) < 100*(pi/100). The fact that pi is not
>> rational – or rather, that float64(pi) cannot be expressed as the
>> division
>> of two 24-bit integers as a 64-bit float – prevents rational lifting of
>> the
>> range from kicking in. I worried about this kind of issue when I was
>> working on FloatRanges, but I'm not sure what you can really do, aside
>> from
>> hacks where you just decide that things are "close enough" based on some
>> ad
>> hoc notion of close enough (Matlab uses 3 ulps). For example, you can't
>> notice that pi/(pi/100) is an integer – because it isn't:
>>
>> julia> pi/(pi/100)
>> 99.99
>>
>>
>> One approach is to try to find a real value x such that
>> float64(x/100) == float64(pi)/100 and float64(x) == float64(pi). If any
>> such value exists, it makes sense to do a lifted FloatRange instead of
>> the
>> default naive stepping seen here. In this case there obviously exists
>> such
>> a real number – π itself is one such value. However, that doesn't quite
>> solve the problem since many such values exist and they don't
>> necessarily
>> all produce the same range values – which one should be used? In this
>> case,
>> π is a good guess, but only because we know that's a special and
>> important
>> number.
Re: [julia-users] Surprising range behavior
Thanks, I hadn't been aware of the range() function (even though it is
right there in the docs). So this works:
julia> range(0.0, pi/100.0, 101)
0.0:0.031415926535897934:3.1415926535897936
The final value is sufficiently close to pi:
julia> range(0.0, pi/100.0, 101)[end] - pi
4.440892098500626e-16
This, and the solution suggested by Simon Kornblith both work fine for
setting up inputs to CoordInterpGrid. So the functionality is definitely
there. But I think that most new users (especially Matlab converts) are
going to try the colon notation first, and perhaps struggle a bit in
figuring out the right way to generate the desired range.
--Peter
On Thursday, April 24, 2014 7:19:59 AM UTC-7, Tobias Knopp wrote:
>
> Shouldn't it be possible to use Range{T}(start::T,step,len::Integer) for
> this?
>
> Am Donnerstag, 24. April 2014 16:09:51 UTC+2 schrieb Peter Simon:
>>
>> At present, the Grid package function CoordInterpGrid will not accept a
>> linspace-generated vector as the x-coordinates, making precise,
>> coordinate-based interpolation a little problematic.
>>
>> --Peter
>>
>> On Thursday, April 24, 2014 6:51:40 AM UTC-7, J Luis wrote:
>>>
>>> I also agree that the Matlab behavior is more intuitive, but even it can
>>> fail for the same reason discussed here. After being bitten by one of such
>>> case that took me a while to debug, when I need the exact ends I now always
>>> use linspace()
>>>
>>> Quinta-feira, 24 de Abril de 2014 14:28:16 UTC+1, Peter Simon escreveu:
*Matlab behavior:*
>> format long
>> x = 0:pi/100:pi;
>> length(x)
ans =
101
>> x(1) - 0
ans =
0
>> x(end) - pi
ans =
0
>> diff([max(diff(x)), min(diff(x))])
ans =
-4.440892098500626e-16
*Julia behavior:*
julia> x = 0:pi/100:pi;
julia> length(x)
100
julia> x[1]-0
0.0
julia> x[end]-pi
-0.031415926535897754
julia> diff([maximum(diff(x)), minimum(diff(x))])
1-element Array{Float64,1}:
-4.44089e-16
I highlighted the important differences in red. IMO the Matlab
behavior is more intuitive. If you choose an increment that "very nearly"
(apparently within a few ulp, as stated by Stefan) evenly divides the
difference between the first and last entries in the colon notation,
Matlab
assumes that you want the resulting vector to begin and end with exactly
these values. This makes sense in most situations, IMO. E.g. when using
a
range to specify sample points in a numerical integration scheme, or when
using a range to specify x-coordinates in an interpolation scheme. In
both
of these cases you don't want the last point to be adjusted away from what
you specified using the colon notation.
--Pete
--Peter
On Thursday, April 24, 2014 5:28:19 AM UTC-7, andrew cooke wrote:
>
>
> how does matlab differ from julia here? (i though the original
> problem was related to "fundamental" issues of how you represent numbers
> on
> a computer, not language specific).
>
>
> On Thursday, 24 April 2014 02:24:59 UTC-3, Peter Simon wrote:
>>
>> Thanks, Simon, that construct works nicely to solve the problem I
>> posed.
>>
>> I have to say, though, that I find Matlab's colon range behavior more
>> intuitive and generally useful, even if it isn't as "exact" as Julia's.
>>
>> --Peter
>>
>> On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote:
>>>
>>> pi*(0:0.01:1) or similar should work.
>>>
>>> On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote:
Thanks for the explanation--it makes sense now. This question
arose for me because of the example presented in
https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ:
x = [0:0.01:pi]
used as the set of x-coordinates for sampling a function to be
integrated (ideally over the interval (0,pi)). But the range defined
in x
has a last entry of 3.14, which will contribute to the inaccuracy of
the
integral being sought in that example. As an exercise, I was trying
to
implement the interpolation solution described later in that thread by
Cameron McBride: "BTW, another possibility is to use a spline
interpolation on the original data and integrate the spline evaluation
with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200)
for
the x values, because CoordInterpGrid will not accept an array as its
first
argument, so you have to use a range object. But the range object has
a
b
Re: [julia-users] Surprising range behavior
Shouldn't it be possible to use Range{T}(start::T,step,len::Integer) for
this?
Am Donnerstag, 24. April 2014 16:09:51 UTC+2 schrieb Peter Simon:
>
> At present, the Grid package function CoordInterpGrid will not accept a
> linspace-generated vector as the x-coordinates, making precise,
> coordinate-based interpolation a little problematic.
>
> --Peter
>
> On Thursday, April 24, 2014 6:51:40 AM UTC-7, J Luis wrote:
>>
>> I also agree that the Matlab behavior is more intuitive, but even it can
>> fail for the same reason discussed here. After being bitten by one of such
>> case that took me a while to debug, when I need the exact ends I now always
>> use linspace()
>>
>> Quinta-feira, 24 de Abril de 2014 14:28:16 UTC+1, Peter Simon escreveu:
>>>
>>> *Matlab behavior:*
>>>
>>> >> format long
>>> >> x = 0:pi/100:pi;
>>> >> length(x)
>>>
>>> ans =
>>>
>>>101
>>>
>>> >> x(1) - 0
>>>
>>> ans =
>>>
>>> 0
>>>
>>> >> x(end) - pi
>>>
>>> ans =
>>>
>>> 0
>>>
>>> >> diff([max(diff(x)), min(diff(x))])
>>>
>>> ans =
>>>
>>> -4.440892098500626e-16
>>>
>>>
>>> *Julia behavior:*
>>>
>>> julia> x = 0:pi/100:pi;
>>>
>>> julia> length(x)
>>> 100
>>>
>>> julia> x[1]-0
>>> 0.0
>>>
>>> julia> x[end]-pi
>>> -0.031415926535897754
>>>
>>> julia> diff([maximum(diff(x)), minimum(diff(x))])
>>> 1-element Array{Float64,1}:
>>> -4.44089e-16
>>>
>>>
>>> I highlighted the important differences in red. IMO the Matlab behavior
>>> is more intuitive. If you choose an increment that "very nearly"
>>> (apparently within a few ulp, as stated by Stefan) evenly divides the
>>> difference between the first and last entries in the colon notation, Matlab
>>> assumes that you want the resulting vector to begin and end with exactly
>>> these values. This makes sense in most situations, IMO. E.g. when using a
>>> range to specify sample points in a numerical integration scheme, or when
>>> using a range to specify x-coordinates in an interpolation scheme. In both
>>> of these cases you don't want the last point to be adjusted away from what
>>> you specified using the colon notation.
>>>
>>> --Pete
>>>
>>> --Peter
>>>
>>> On Thursday, April 24, 2014 5:28:19 AM UTC-7, andrew cooke wrote:
how does matlab differ from julia here? (i though the original problem
was related to "fundamental" issues of how you represent numbers on a
computer, not language specific).
On Thursday, 24 April 2014 02:24:59 UTC-3, Peter Simon wrote:
>
> Thanks, Simon, that construct works nicely to solve the problem I
> posed.
>
> I have to say, though, that I find Matlab's colon range behavior more
> intuitive and generally useful, even if it isn't as "exact" as Julia's.
>
> --Peter
>
> On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote:
>>
>> pi*(0:0.01:1) or similar should work.
>>
>> On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote:
>>>
>>> Thanks for the explanation--it makes sense now. This question arose
>>> for me because of the example presented in
>>> https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ:
>>>
>>> x = [0:0.01:pi]
>>>
>>> used as the set of x-coordinates for sampling a function to be
>>> integrated (ideally over the interval (0,pi)). But the range defined
>>> in x
>>> has a last entry of 3.14, which will contribute to the inaccuracy of
>>> the
>>> integral being sought in that example. As an exercise, I was trying to
>>> implement the interpolation solution described later in that thread by
>>> Cameron McBride: "BTW, another possibility is to use a spline
>>> interpolation on the original data and integrate the spline evaluation
>>> with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200)
>>> for
>>> the x values, because CoordInterpGrid will not accept an array as its
>>> first
>>> argument, so you have to use a range object. But the range object has
>>> a
>>> built-in error for the last point because of the present issue. Any
>>> suggestions?
>>>
>>> Thanks,
>>>
>>> --Peter
>>>
>>> On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski
>>> wrote:
The issue is that float(pi) < 100*(pi/100). The fact that pi is not
rational – or rather, that float64(pi) cannot be expressed as the
division
of two 24-bit integers as a 64-bit float – prevents rational lifting
of the
range from kicking in. I worried about this kind of issue when I was
working on FloatRanges, but I'm not sure what you can really do, aside
from
hacks where you just decide that things are "close enough" based on
some ad
hoc notion of close enough (Matlab uses 3 ulps). For example, you
can't
notice that pi/(
Re: [julia-users] Surprising range behavior
At present, the Grid package function CoordInterpGrid will not accept a
linspace-generated vector as the x-coordinates, making precise,
coordinate-based interpolation a little problematic.
--Peter
On Thursday, April 24, 2014 6:51:40 AM UTC-7, J Luis wrote:
>
> I also agree that the Matlab behavior is more intuitive, but even it can
> fail for the same reason discussed here. After being bitten by one of such
> case that took me a while to debug, when I need the exact ends I now always
> use linspace()
>
> Quinta-feira, 24 de Abril de 2014 14:28:16 UTC+1, Peter Simon escreveu:
>>
>> *Matlab behavior:*
>>
>> >> format long
>> >> x = 0:pi/100:pi;
>> >> length(x)
>>
>> ans =
>>
>>101
>>
>> >> x(1) - 0
>>
>> ans =
>>
>> 0
>>
>> >> x(end) - pi
>>
>> ans =
>>
>> 0
>>
>> >> diff([max(diff(x)), min(diff(x))])
>>
>> ans =
>>
>> -4.440892098500626e-16
>>
>>
>> *Julia behavior:*
>>
>> julia> x = 0:pi/100:pi;
>>
>> julia> length(x)
>> 100
>>
>> julia> x[1]-0
>> 0.0
>>
>> julia> x[end]-pi
>> -0.031415926535897754
>>
>> julia> diff([maximum(diff(x)), minimum(diff(x))])
>> 1-element Array{Float64,1}:
>> -4.44089e-16
>>
>>
>> I highlighted the important differences in red. IMO the Matlab behavior
>> is more intuitive. If you choose an increment that "very nearly"
>> (apparently within a few ulp, as stated by Stefan) evenly divides the
>> difference between the first and last entries in the colon notation, Matlab
>> assumes that you want the resulting vector to begin and end with exactly
>> these values. This makes sense in most situations, IMO. E.g. when using a
>> range to specify sample points in a numerical integration scheme, or when
>> using a range to specify x-coordinates in an interpolation scheme. In both
>> of these cases you don't want the last point to be adjusted away from what
>> you specified using the colon notation.
>>
>> --Pete
>>
>> --Peter
>>
>> On Thursday, April 24, 2014 5:28:19 AM UTC-7, andrew cooke wrote:
>>>
>>>
>>> how does matlab differ from julia here? (i though the original problem
>>> was related to "fundamental" issues of how you represent numbers on a
>>> computer, not language specific).
>>>
>>>
>>> On Thursday, 24 April 2014 02:24:59 UTC-3, Peter Simon wrote:
Thanks, Simon, that construct works nicely to solve the problem I
posed.
I have to say, though, that I find Matlab's colon range behavior more
intuitive and generally useful, even if it isn't as "exact" as Julia's.
--Peter
On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote:
>
> pi*(0:0.01:1) or similar should work.
>
> On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote:
>>
>> Thanks for the explanation--it makes sense now. This question arose
>> for me because of the example presented in
>> https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ:
>>
>> x = [0:0.01:pi]
>>
>> used as the set of x-coordinates for sampling a function to be
>> integrated (ideally over the interval (0,pi)). But the range defined in
>> x
>> has a last entry of 3.14, which will contribute to the inaccuracy of the
>> integral being sought in that example. As an exercise, I was trying to
>> implement the interpolation solution described later in that thread by
>> Cameron McBride: "BTW, another possibility is to use a spline
>> interpolation on the original data and integrate the spline evaluation
>> with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200)
>> for
>> the x values, because CoordInterpGrid will not accept an array as its
>> first
>> argument, so you have to use a range object. But the range object has a
>> built-in error for the last point because of the present issue. Any
>> suggestions?
>>
>> Thanks,
>>
>> --Peter
>>
>> On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote:
>>>
>>> The issue is that float(pi) < 100*(pi/100). The fact that pi is not
>>> rational – or rather, that float64(pi) cannot be expressed as the
>>> division
>>> of two 24-bit integers as a 64-bit float – prevents rational lifting of
>>> the
>>> range from kicking in. I worried about this kind of issue when I was
>>> working on FloatRanges, but I'm not sure what you can really do, aside
>>> from
>>> hacks where you just decide that things are "close enough" based on
>>> some ad
>>> hoc notion of close enough (Matlab uses 3 ulps). For example, you can't
>>> notice that pi/(pi/100) is an integer – because it isn't:
>>>
>>> julia> pi/(pi/100)
>>> 99.99
>>>
>>>
>>> One approach is to try to find a real value x such that
>>> float64(x/100) == float64(pi)/100 and float64(x) == float64(pi). If any
>>> such value exists, it makes sense to do a lifted FloatRange
Re: [julia-users] Surprising range behavior
I also agree that the Matlab behavior is more intuitive, but even it can
fail for the same reason discussed here. After being bitten by one of such
case that took me a while to debug, when I need the exact ends I now always
use linspace()
Quinta-feira, 24 de Abril de 2014 14:28:16 UTC+1, Peter Simon escreveu:
>
> *Matlab behavior:*
>
> >> format long
> >> x = 0:pi/100:pi;
> >> length(x)
>
> ans =
>
>101
>
> >> x(1) - 0
>
> ans =
>
> 0
>
> >> x(end) - pi
>
> ans =
>
> 0
>
> >> diff([max(diff(x)), min(diff(x))])
>
> ans =
>
> -4.440892098500626e-16
>
>
> *Julia behavior:*
>
> julia> x = 0:pi/100:pi;
>
> julia> length(x)
> 100
>
> julia> x[1]-0
> 0.0
>
> julia> x[end]-pi
> -0.031415926535897754
>
> julia> diff([maximum(diff(x)), minimum(diff(x))])
> 1-element Array{Float64,1}:
> -4.44089e-16
>
>
> I highlighted the important differences in red. IMO the Matlab behavior
> is more intuitive. If you choose an increment that "very nearly"
> (apparently within a few ulp, as stated by Stefan) evenly divides the
> difference between the first and last entries in the colon notation, Matlab
> assumes that you want the resulting vector to begin and end with exactly
> these values. This makes sense in most situations, IMO. E.g. when using a
> range to specify sample points in a numerical integration scheme, or when
> using a range to specify x-coordinates in an interpolation scheme. In both
> of these cases you don't want the last point to be adjusted away from what
> you specified using the colon notation.
>
> --Pete
>
> --Peter
>
> On Thursday, April 24, 2014 5:28:19 AM UTC-7, andrew cooke wrote:
>>
>>
>> how does matlab differ from julia here? (i though the original problem
>> was related to "fundamental" issues of how you represent numbers on a
>> computer, not language specific).
>>
>>
>> On Thursday, 24 April 2014 02:24:59 UTC-3, Peter Simon wrote:
>>>
>>> Thanks, Simon, that construct works nicely to solve the problem I posed.
>>>
>>>
>>> I have to say, though, that I find Matlab's colon range behavior more
>>> intuitive and generally useful, even if it isn't as "exact" as Julia's.
>>>
>>> --Peter
>>>
>>> On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote:
pi*(0:0.01:1) or similar should work.
On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote:
>
> Thanks for the explanation--it makes sense now. This question arose
> for me because of the example presented in
> https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ :
>
> x = [0:0.01:pi]
>
> used as the set of x-coordinates for sampling a function to be
> integrated (ideally over the interval (0,pi)). But the range defined in
> x
> has a last entry of 3.14, which will contribute to the inaccuracy of the
> integral being sought in that example. As an exercise, I was trying to
> implement the interpolation solution described later in that thread by
> Cameron McBride: "BTW, another possibility is to use a spline
> interpolation on the original data and integrate the spline evaluation
> with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200)
> for
> the x values, because CoordInterpGrid will not accept an array as its
> first
> argument, so you have to use a range object. But the range object has a
> built-in error for the last point because of the present issue. Any
> suggestions?
>
> Thanks,
>
> --Peter
>
> On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote:
>>
>> The issue is that float(pi) < 100*(pi/100). The fact that pi is not
>> rational – or rather, that float64(pi) cannot be expressed as the
>> division
>> of two 24-bit integers as a 64-bit float – prevents rational lifting of
>> the
>> range from kicking in. I worried about this kind of issue when I was
>> working on FloatRanges, but I'm not sure what you can really do, aside
>> from
>> hacks where you just decide that things are "close enough" based on some
>> ad
>> hoc notion of close enough (Matlab uses 3 ulps). For example, you can't
>> notice that pi/(pi/100) is an integer – because it isn't:
>>
>> julia> pi/(pi/100)
>> 99.99
>>
>>
>> One approach is to try to find a real value x such that
>> float64(x/100) == float64(pi)/100 and float64(x) == float64(pi). If any
>> such value exists, it makes sense to do a lifted FloatRange instead of
>> the
>> default naive stepping seen here. In this case there obviously exists
>> such
>> a real number – π itself is one such value. However, that doesn't quite
>> solve the problem since many such values exist and they don't
>> necessarily
>> all produce the same range values – which one should be used? In this
>> case,
>> π is a good guess, but
Re: [julia-users] Surprising range behavior
*Matlab behavior:*
>> format long
>> x = 0:pi/100:pi;
>> length(x)
ans =
101
>> x(1) - 0
ans =
0
>> x(end) - pi
ans =
0
>> diff([max(diff(x)), min(diff(x))])
ans =
-4.440892098500626e-16
*Julia behavior:*
julia> x = 0:pi/100:pi;
julia> length(x)
100
julia> x[1]-0
0.0
julia> x[end]-pi
-0.031415926535897754
julia> diff([maximum(diff(x)), minimum(diff(x))])
1-element Array{Float64,1}:
-4.44089e-16
I highlighted the important differences in red. IMO the Matlab behavior is
more intuitive. If you choose an increment that "very nearly" (apparently
within a few ulp, as stated by Stefan) evenly divides the difference
between the first and last entries in the colon notation, Matlab assumes
that you want the resulting vector to begin and end with exactly these
values. This makes sense in most situations, IMO. E.g. when using a range
to specify sample points in a numerical integration scheme, or when using a
range to specify x-coordinates in an interpolation scheme. In both of
these cases you don't want the last point to be adjusted away from what you
specified using the colon notation.
--Pete
--Peter
On Thursday, April 24, 2014 5:28:19 AM UTC-7, andrew cooke wrote:
>
>
> how does matlab differ from julia here? (i though the original problem
> was related to "fundamental" issues of how you represent numbers on a
> computer, not language specific).
>
>
> On Thursday, 24 April 2014 02:24:59 UTC-3, Peter Simon wrote:
>>
>> Thanks, Simon, that construct works nicely to solve the problem I posed.
>>
>> I have to say, though, that I find Matlab's colon range behavior more
>> intuitive and generally useful, even if it isn't as "exact" as Julia's.
>>
>> --Peter
>>
>> On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote:
>>>
>>> pi*(0:0.01:1) or similar should work.
>>>
>>> On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote:
Thanks for the explanation--it makes sense now. This question arose
for me because of the example presented in
https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ :
x = [0:0.01:pi]
used as the set of x-coordinates for sampling a function to be
integrated (ideally over the interval (0,pi)). But the range defined in x
has a last entry of 3.14, which will contribute to the inaccuracy of the
integral being sought in that example. As an exercise, I was trying to
implement the interpolation solution described later in that thread by
Cameron McBride: "BTW, another possibility is to use a spline
interpolation on the original data and integrate the spline evaluation
with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200) for
the x values, because CoordInterpGrid will not accept an array as its
first
argument, so you have to use a range object. But the range object has a
built-in error for the last point because of the present issue. Any
suggestions?
Thanks,
--Peter
On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote:
>
> The issue is that float(pi) < 100*(pi/100). The fact that pi is not
> rational – or rather, that float64(pi) cannot be expressed as the
> division
> of two 24-bit integers as a 64-bit float – prevents rational lifting of
> the
> range from kicking in. I worried about this kind of issue when I was
> working on FloatRanges, but I'm not sure what you can really do, aside
> from
> hacks where you just decide that things are "close enough" based on some
> ad
> hoc notion of close enough (Matlab uses 3 ulps). For example, you can't
> notice that pi/(pi/100) is an integer – because it isn't:
>
> julia> pi/(pi/100)
> 99.99
>
>
> One approach is to try to find a real value x such that float64(x/100)
> == float64(pi)/100 and float64(x) == float64(pi). If any such value
> exists,
> it makes sense to do a lifted FloatRange instead of the default naive
> stepping seen here. In this case there obviously exists such a real
> number
> – π itself is one such value. However, that doesn't quite solve the
> problem
> since many such values exist and they don't necessarily all produce the
> same range values – which one should be used? In this case, π is a good
> guess, but only because we know that's a special and important number.
> Adding in ad hoc special values isn't really satisfying or acceptable. It
> would be nice to give the right behavior in cases where there is only one
> possible range that could have been intended (despite there being many
> values of x), but I haven't figured out how determine if that is the case
> or not. The current code handles the relatively straightforward case
> where
> the start, step and stop values are all ratio
Re: [julia-users] Surprising range behavior
how does matlab differ from julia here? (i though the original problem was related to "fundamental" issues of how you represent numbers on a computer, not language specific). On Thursday, 24 April 2014 02:24:59 UTC-3, Peter Simon wrote: > > Thanks, Simon, that construct works nicely to solve the problem I posed. > > I have to say, though, that I find Matlab's colon range behavior more > intuitive and generally useful, even if it isn't as "exact" as Julia's. > > --Peter > > On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote: >> >> pi*(0:0.01:1) or similar should work. >> >> On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote: >>> >>> Thanks for the explanation--it makes sense now. This question arose for >>> me because of the example presented in >>> https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ : >>> >>> x = [0:0.01:pi] >>> >>> used as the set of x-coordinates for sampling a function to be >>> integrated (ideally over the interval (0,pi)). But the range defined in x >>> has a last entry of 3.14, which will contribute to the inaccuracy of the >>> integral being sought in that example. As an exercise, I was trying to >>> implement the interpolation solution described later in that thread by >>> Cameron McBride: "BTW, another possibility is to use a spline >>> interpolation on the original data and integrate the spline evaluation >>> with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200) for >>> the x values, because CoordInterpGrid will not accept an array as its first >>> argument, so you have to use a range object. But the range object has a >>> built-in error for the last point because of the present issue. Any >>> suggestions? >>> >>> Thanks, >>> >>> --Peter >>> >>> On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote: The issue is that float(pi) < 100*(pi/100). The fact that pi is not rational – or rather, that float64(pi) cannot be expressed as the division of two 24-bit integers as a 64-bit float – prevents rational lifting of the range from kicking in. I worried about this kind of issue when I was working on FloatRanges, but I'm not sure what you can really do, aside from hacks where you just decide that things are "close enough" based on some ad hoc notion of close enough (Matlab uses 3 ulps). For example, you can't notice that pi/(pi/100) is an integer – because it isn't: julia> pi/(pi/100) 99.99 One approach is to try to find a real value x such that float64(x/100) == float64(pi)/100 and float64(x) == float64(pi). If any such value exists, it makes sense to do a lifted FloatRange instead of the default naive stepping seen here. In this case there obviously exists such a real number – π itself is one such value. However, that doesn't quite solve the problem since many such values exist and they don't necessarily all produce the same range values – which one should be used? In this case, π is a good guess, but only because we know that's a special and important number. Adding in ad hoc special values isn't really satisfying or acceptable. It would be nice to give the right behavior in cases where there is only one possible range that could have been intended (despite there being many values of x), but I haven't figured out how determine if that is the case or not. The current code handles the relatively straightforward case where the start, step and stop values are all rational. On Wed, Apr 23, 2014 at 5:59 PM, Peter Simon wrote: > The first three results below are what I expected. The fourth result > surprised me: > > julia> (0:pi:pi)[end] > 3.141592653589793 > > julia> (0:pi/2:pi)[end] > 3.141592653589793 > > julia> (0:pi/3:pi)[end] > 3.141592653589793 > > julia> (0:pi/100:pi)[end] > 3.1101767270538954 > > Is this behavior correct? > > Version info: > julia> versioninfo() > Julia Version 0.3.0-prerelease+2703 > Commit 942ae42* (2014-04-22 18:57 UTC) > Platform Info: > System: Windows (x86_64-w64-mingw32) > CPU: Intel(R) Core(TM) i7 CPU 860 @ 2.80GHz > WORD_SIZE: 64 > BLAS: libopenblas (USE64BITINT DYNAMIC_ARCH NO_AFFINITY) > LAPACK: libopenblas > LIBM: libopenlibm > > > --Peter > >
Re: [julia-users] Surprising range behavior
Thanks, Simon, that construct works nicely to solve the problem I posed. I have to say, though, that I find Matlab's colon range behavior more intuitive and generally useful, even if it isn't as "exact" as Julia's. --Peter On Wednesday, April 23, 2014 7:17:23 PM UTC-7, Simon Kornblith wrote: > > pi*(0:0.01:1) or similar should work. > > On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote: >> >> Thanks for the explanation--it makes sense now. This question arose for >> me because of the example presented in >> https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ : >> >> x = [0:0.01:pi] >> >> used as the set of x-coordinates for sampling a function to be integrated >> (ideally over the interval (0,pi)). But the range defined in x has a last >> entry of 3.14, which will contribute to the inaccuracy of the integral >> being sought in that example. As an exercise, I was trying to implement >> the interpolation solution described later in that thread by Cameron >> McBride: "BTW, another possibility is to use a spline interpolation on the >> original data and integrate the spline evaluation with quadgk()". It >> seems that one cannot use e.g. linspace(0,pi,200) for the x values, because >> CoordInterpGrid will not accept an array as its first argument, so you have >> to use a range object. But the range object has a built-in error for the >> last point because of the present issue. Any suggestions? >> >> Thanks, >> >> --Peter >> >> On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote: >>> >>> The issue is that float(pi) < 100*(pi/100). The fact that pi is not >>> rational – or rather, that float64(pi) cannot be expressed as the division >>> of two 24-bit integers as a 64-bit float – prevents rational lifting of the >>> range from kicking in. I worried about this kind of issue when I was >>> working on FloatRanges, but I'm not sure what you can really do, aside from >>> hacks where you just decide that things are "close enough" based on some ad >>> hoc notion of close enough (Matlab uses 3 ulps). For example, you can't >>> notice that pi/(pi/100) is an integer – because it isn't: >>> >>> julia> pi/(pi/100) >>> 99.99 >>> >>> >>> One approach is to try to find a real value x such that float64(x/100) >>> == float64(pi)/100 and float64(x) == float64(pi). If any such value exists, >>> it makes sense to do a lifted FloatRange instead of the default naive >>> stepping seen here. In this case there obviously exists such a real number >>> – π itself is one such value. However, that doesn't quite solve the problem >>> since many such values exist and they don't necessarily all produce the >>> same range values – which one should be used? In this case, π is a good >>> guess, but only because we know that's a special and important number. >>> Adding in ad hoc special values isn't really satisfying or acceptable. It >>> would be nice to give the right behavior in cases where there is only one >>> possible range that could have been intended (despite there being many >>> values of x), but I haven't figured out how determine if that is the case >>> or not. The current code handles the relatively straightforward case where >>> the start, step and stop values are all rational. >>> >>> >>> On Wed, Apr 23, 2014 at 5:59 PM, Peter Simon wrote: >>> The first three results below are what I expected. The fourth result surprised me: julia> (0:pi:pi)[end] 3.141592653589793 julia> (0:pi/2:pi)[end] 3.141592653589793 julia> (0:pi/3:pi)[end] 3.141592653589793 julia> (0:pi/100:pi)[end] 3.1101767270538954 Is this behavior correct? Version info: julia> versioninfo() Julia Version 0.3.0-prerelease+2703 Commit 942ae42* (2014-04-22 18:57 UTC) Platform Info: System: Windows (x86_64-w64-mingw32) CPU: Intel(R) Core(TM) i7 CPU 860 @ 2.80GHz WORD_SIZE: 64 BLAS: libopenblas (USE64BITINT DYNAMIC_ARCH NO_AFFINITY) LAPACK: libopenblas LIBM: libopenlibm --Peter >>>
Re: [julia-users] Surprising range behavior
(Simon, you may also be amused to learn that Google thinks that your post is written in Latin and offers to translate it. Mirabile dictu!)
Re: [julia-users] Surprising range behavior
On Wednesday, April 23, 2014 10:50:57 PM UTC-4, Steven G. Johnson wrote: > On Wednesday, April 23, 2014 10:17:23 PM UTC-4, Simon Kornblith wrote: >> >> pi*(0:0.01:1) or similar should work. >> > > Actually, that may not work because of > https://github.com/JuliaLang/julia/issues/6364 > ...which of course you know about because you submitted that bug report, sorry. (People should just use their Github usernames everywhere, rather than polluting the namespace with all of these "real name" aliases.)
Re: [julia-users] Surprising range behavior
I believe that bug only applies to multiplication/division of integer ranges, but it will certainly be good to have it fixed. On Wednesday, April 23, 2014 10:50:57 PM UTC-4, Steven G. Johnson wrote: > > > > On Wednesday, April 23, 2014 10:17:23 PM UTC-4, Simon Kornblith wrote: >> >> pi*(0:0.01:1) or similar should work. >> > > Actually, that may not work because of > https://github.com/JuliaLang/julia/issues/6364 > > However, this should be fixable (and in fact I just submitted a PR for > it). >
Re: [julia-users] Surprising range behavior
On Wednesday, April 23, 2014 10:17:23 PM UTC-4, Simon Kornblith wrote: > > pi*(0:0.01:1) or similar should work. > Actually, that may not work because of https://github.com/JuliaLang/julia/issues/6364 However, this should be fixable (and in fact I just submitted a PR for it).
Re: [julia-users] Surprising range behavior
pi*(0:0.01:1) or similar should work. On Wednesday, April 23, 2014 7:12:58 PM UTC-4, Peter Simon wrote: > > Thanks for the explanation--it makes sense now. This question arose for > me because of the example presented in > https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ : > > x = [0:0.01:pi] > > used as the set of x-coordinates for sampling a function to be integrated > (ideally over the interval (0,pi)). But the range defined in x has a last > entry of 3.14, which will contribute to the inaccuracy of the integral > being sought in that example. As an exercise, I was trying to implement > the interpolation solution described later in that thread by Cameron > McBride: "BTW, another possibility is to use a spline interpolation on the > original data and integrate the spline evaluation with quadgk()". It > seems that one cannot use e.g. linspace(0,pi,200) for the x values, because > CoordInterpGrid will not accept an array as its first argument, so you have > to use a range object. But the range object has a built-in error for the > last point because of the present issue. Any suggestions? > > Thanks, > > --Peter > > On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote: >> >> The issue is that float(pi) < 100*(pi/100). The fact that pi is not >> rational – or rather, that float64(pi) cannot be expressed as the division >> of two 24-bit integers as a 64-bit float – prevents rational lifting of the >> range from kicking in. I worried about this kind of issue when I was >> working on FloatRanges, but I'm not sure what you can really do, aside from >> hacks where you just decide that things are "close enough" based on some ad >> hoc notion of close enough (Matlab uses 3 ulps). For example, you can't >> notice that pi/(pi/100) is an integer – because it isn't: >> >> julia> pi/(pi/100) >> 99.99 >> >> >> One approach is to try to find a real value x such that float64(x/100) == >> float64(pi)/100 and float64(x) == float64(pi). If any such value exists, it >> makes sense to do a lifted FloatRange instead of the default naive stepping >> seen here. In this case there obviously exists such a real number – π >> itself is one such value. However, that doesn't quite solve the problem >> since many such values exist and they don't necessarily all produce the >> same range values – which one should be used? In this case, π is a good >> guess, but only because we know that's a special and important number. >> Adding in ad hoc special values isn't really satisfying or acceptable. It >> would be nice to give the right behavior in cases where there is only one >> possible range that could have been intended (despite there being many >> values of x), but I haven't figured out how determine if that is the case >> or not. The current code handles the relatively straightforward case where >> the start, step and stop values are all rational. >> >> >> On Wed, Apr 23, 2014 at 5:59 PM, Peter Simon wrote: >> >>> The first three results below are what I expected. The fourth result >>> surprised me: >>> >>> julia> (0:pi:pi)[end] >>> 3.141592653589793 >>> >>> julia> (0:pi/2:pi)[end] >>> 3.141592653589793 >>> >>> julia> (0:pi/3:pi)[end] >>> 3.141592653589793 >>> >>> julia> (0:pi/100:pi)[end] >>> 3.1101767270538954 >>> >>> Is this behavior correct? >>> >>> Version info: >>> julia> versioninfo() >>> Julia Version 0.3.0-prerelease+2703 >>> Commit 942ae42* (2014-04-22 18:57 UTC) >>> Platform Info: >>> System: Windows (x86_64-w64-mingw32) >>> CPU: Intel(R) Core(TM) i7 CPU 860 @ 2.80GHz >>> WORD_SIZE: 64 >>> BLAS: libopenblas (USE64BITINT DYNAMIC_ARCH NO_AFFINITY) >>> LAPACK: libopenblas >>> LIBM: libopenlibm >>> >>> >>> --Peter >>> >>> >>
Re: [julia-users] Surprising range behavior
Thanks for the explanation--it makes sense now. This question arose for me because of the example presented in https://groups.google.com/d/msg/julia-users/CNYaDUYog8w/QH9L_Q9Su9YJ : x = [0:0.01:pi] used as the set of x-coordinates for sampling a function to be integrated (ideally over the interval (0,pi)). But the range defined in x has a last entry of 3.14, which will contribute to the inaccuracy of the integral being sought in that example. As an exercise, I was trying to implement the interpolation solution described later in that thread by Cameron McBride: "BTW, another possibility is to use a spline interpolation on the original data and integrate the spline evaluation with quadgk()". It seems that one cannot use e.g. linspace(0,pi,200) for the x values, because CoordInterpGrid will not accept an array as its first argument, so you have to use a range object. But the range object has a built-in error for the last point because of the present issue. Any suggestions? Thanks, --Peter On Wednesday, April 23, 2014 3:24:10 PM UTC-7, Stefan Karpinski wrote: > > The issue is that float(pi) < 100*(pi/100). The fact that pi is not > rational – or rather, that float64(pi) cannot be expressed as the division > of two 24-bit integers as a 64-bit float – prevents rational lifting of the > range from kicking in. I worried about this kind of issue when I was > working on FloatRanges, but I'm not sure what you can really do, aside from > hacks where you just decide that things are "close enough" based on some ad > hoc notion of close enough (Matlab uses 3 ulps). For example, you can't > notice that pi/(pi/100) is an integer – because it isn't: > > julia> pi/(pi/100) > 99.99 > > > One approach is to try to find a real value x such that float64(x/100) == > float64(pi)/100 and float64(x) == float64(pi). If any such value exists, it > makes sense to do a lifted FloatRange instead of the default naive stepping > seen here. In this case there obviously exists such a real number – π > itself is one such value. However, that doesn't quite solve the problem > since many such values exist and they don't necessarily all produce the > same range values – which one should be used? In this case, π is a good > guess, but only because we know that's a special and important number. > Adding in ad hoc special values isn't really satisfying or acceptable. It > would be nice to give the right behavior in cases where there is only one > possible range that could have been intended (despite there being many > values of x), but I haven't figured out how determine if that is the case > or not. The current code handles the relatively straightforward case where > the start, step and stop values are all rational. > > > On Wed, Apr 23, 2014 at 5:59 PM, Peter Simon > > wrote: > >> The first three results below are what I expected. The fourth result >> surprised me: >> >> julia> (0:pi:pi)[end] >> 3.141592653589793 >> >> julia> (0:pi/2:pi)[end] >> 3.141592653589793 >> >> julia> (0:pi/3:pi)[end] >> 3.141592653589793 >> >> julia> (0:pi/100:pi)[end] >> 3.1101767270538954 >> >> Is this behavior correct? >> >> Version info: >> julia> versioninfo() >> Julia Version 0.3.0-prerelease+2703 >> Commit 942ae42* (2014-04-22 18:57 UTC) >> Platform Info: >> System: Windows (x86_64-w64-mingw32) >> CPU: Intel(R) Core(TM) i7 CPU 860 @ 2.80GHz >> WORD_SIZE: 64 >> BLAS: libopenblas (USE64BITINT DYNAMIC_ARCH NO_AFFINITY) >> LAPACK: libopenblas >> LIBM: libopenlibm >> >> >> --Peter >> >> >
Re: [julia-users] Surprising range behavior
The issue is that float(pi) < 100*(pi/100). The fact that pi is not rational – or rather, that float64(pi) cannot be expressed as the division of two 24-bit integers as a 64-bit float – prevents rational lifting of the range from kicking in. I worried about this kind of issue when I was working on FloatRanges, but I'm not sure what you can really do, aside from hacks where you just decide that things are "close enough" based on some ad hoc notion of close enough (Matlab uses 3 ulps). For example, you can't notice that pi/(pi/100) is an integer – because it isn't: julia> pi/(pi/100) 99.99 One approach is to try to find a real value x such that float64(x/100) == float64(pi)/100 and float64(x) == float64(pi). If any such value exists, it makes sense to do a lifted FloatRange instead of the default naive stepping seen here. In this case there obviously exists such a real number – π itself is one such value. However, that doesn't quite solve the problem since many such values exist and they don't necessarily all produce the same range values – which one should be used? In this case, π is a good guess, but only because we know that's a special and important number. Adding in ad hoc special values isn't really satisfying or acceptable. It would be nice to give the right behavior in cases where there is only one possible range that could have been intended (despite there being many values of x), but I haven't figured out how determine if that is the case or not. The current code handles the relatively straightforward case where the start, step and stop values are all rational. On Wed, Apr 23, 2014 at 5:59 PM, Peter Simon wrote: > The first three results below are what I expected. The fourth result > surprised me: > > julia> (0:pi:pi)[end] > 3.141592653589793 > > julia> (0:pi/2:pi)[end] > 3.141592653589793 > > julia> (0:pi/3:pi)[end] > 3.141592653589793 > > julia> (0:pi/100:pi)[end] > 3.1101767270538954 > > Is this behavior correct? > > Version info: > julia> versioninfo() > Julia Version 0.3.0-prerelease+2703 > Commit 942ae42* (2014-04-22 18:57 UTC) > Platform Info: > System: Windows (x86_64-w64-mingw32) > CPU: Intel(R) Core(TM) i7 CPU 860 @ 2.80GHz > WORD_SIZE: 64 > BLAS: libopenblas (USE64BITINT DYNAMIC_ARCH NO_AFFINITY) > LAPACK: libopenblas > LIBM: libopenlibm > > > --Peter > >
[julia-users] Surprising range behavior
The first three results below are what I expected. The fourth result surprised me: julia> (0:pi:pi)[end] 3.141592653589793 julia> (0:pi/2:pi)[end] 3.141592653589793 julia> (0:pi/3:pi)[end] 3.141592653589793 julia> (0:pi/100:pi)[end] 3.1101767270538954 Is this behavior correct? Version info: julia> versioninfo() Julia Version 0.3.0-prerelease+2703 Commit 942ae42* (2014-04-22 18:57 UTC) Platform Info: System: Windows (x86_64-w64-mingw32) CPU: Intel(R) Core(TM) i7 CPU 860 @ 2.80GHz WORD_SIZE: 64 BLAS: libopenblas (USE64BITINT DYNAMIC_ARCH NO_AFFINITY) LAPACK: libopenblas LIBM: libopenlibm --Peter
