I agree that the this is a surprising result---even
if correct by our current implementation. Maybe
the default for NV perl scalars should be type
PDL_Double. Maybe a perl scalar promotion
setting that could set a default value, double
or other specific type, smallest size, or ...
--Chris
On Wed, Nov 18, 2015 at 4:58 PM, Craig DeForest <defor...@boulder.swri.edu>
wrote:
> Yes, I believe this is a Bug in the conversion code logic. Currently we
> use the smallest possible type that can contain a value. Values that fit
> in a float get promoted to float and those that require a double get
> promoted to double. Here’s a four-atan trial, showing the autopromotion
> using a Perl scalar constructed to be similar to 2479/1024, that cannot fit
> in a float.
>
>
> sub try {
> my $val = shift;
> printf(“\nval: %.15f\nfloat: %.15f\ndefault: %.15f\ndouble:
> %.15f\n”,
> $val, atan(float($val)), atan($val), atan(double($val))
> );
> }
>
> foreach( float(2379.0)/1024.0, double(2379.0)/1024.0, 2379.0/1024.0,
> 2479.00000001/1024.0 ) { try($_) }
>
>
> yields
>
> val: 2.323242187500000
> float: 1.164332985877990
> default: 1.164332985877990
> double: 1.164332932171300
>
> val: 2.323242187500000
> float: 1.164332985877990
> default: 1.164332932171300
> double: 1.164332932171300
>
> val: 2.323242187500000
> float: 1.164332985877990
> default: 1.164332985877990
> double: 1.164332932171300
>
> val: 2.420898437509766
> float: 1.179073929786680
> default: 1.179073913765370
> double: 1.179073913765370
>
> where the top is the floating-point value, second is the double value (and
> both work as desired), third is the “bug” value (perl scalar is promoted to
> floating point), and fourth is the “antibug” value (perl scalar is promoted
> to double).
> Because a float can represent 2379/1024 exactly (after all, it’s just the
> integer 2379, shifted 10 bits to the right), you get a float. But 1024 is
> a pretty special number, and highly unusual. Changing the ratio even
> slightly yields promotion to double by default.
>
> pdl> try(2479/1025)
>
> val: 2.418536585365854
> float: 1.178729414939880
> default: 1.178729370918130
> double: 1.178729370918130
>
> Since the vast majority of numeric Perl scalars cannot be represented
> exactly in a float, it’s easy to get into the mindset that Perl scalars are
> autopromoted to double PDLs. By the Principle of Least Surprise, we should
> therefore always autopromote to the most common case (double), and require
> explicit cast to get to floating-point values.
>
> Sorry for the wall of text here, this is a pretty cool bug.
>
>
> On Nov 18, 2015, at 1:42 PM, Chris Marshall <devel.chm...@gmail.com>
> wrote:
>
> I think we may need something like
> an option to enable/disable forced
> conversion of perl scalars to double
> piddles.
>
> --Chris
>
> On Wed, Nov 18, 2015 at 3:03 PM, Chris Marshall <devel.chm...@gmail.com>
> wrote:
>
>> Hi Doug-
>>
>> I think the problem is only when
>> a perl scalar is the input to atan().
>> The top two cases are the expected
>> result. Maybe the logic in the type
>> conversion is off. You don't say
>> what version of PDL but I see the
>> same output from PDL-2.014_02.
>>
>> --Chris
>>
>>
>>
>> On Wed, Nov 18, 2015 at 1:29 PM, Douglas Hunt <dh...@ucar.edu> wrote:
>>
>>> Hi all: I just noticed a thing which caused me some trouble.
>>>
>>> The PDL 'atan' function defaults to 'float' precision instead of the
>>> more natural 'double':
>>>
>>> use PDL;
>>>
>>> my $a1 = atan(double(2379.0)/double(1024.0));
>>> my $a2 = atan(float(2379.0)/float(1024.0));
>>> my $a3 = atan(2379.0/1024.0);
>>>
>>> print "$a1, $a2, $a3\n";
>>>
>>> This prints:
>>> 1.1643329321713, 1.1643328666687, 1.1643328666687
>>>
>>> This feature causes a problem when porting to/from C as the C atan
>>> function operates on doubles.
>>>
>>> I think this traces down to Basic::Math in the math.pd file:
>>>
>>> my (@ufuncs1) = qw(acos asin atan cosh sinh tan tanh); # F,D only
>>> ...
>>> my $func;
>>> foreach $func (@ufuncs1) {
>>> pp_def($func,
>>> HandleBad => 1,
>>> NoBadifNaN => 1,
>>> GenericTypes => ['F','D'],
>>> Pars => 'a(); [o]b();',
>>> Inplace => 1,
>>> Doc => inplace_doc( $func ),
>>> Code => code_ufunc($func),
>>> BadCode => badcode_ufunc($func),
>>> );
>>> }
>>>
>>> Note that GenericTypes isF, D instead of D, F.
>>>
>>> I would go so far to say that this is a bug, as the default pdl type is
>>> double as in:
>>>
>>> my $pdl = pdl(4,5,6);
>>> p $pdl->info
>>> PDL: Double D [3]
>>>
>>> Regards,
>>>
>>> Doug Hunt
>>>
>>>
>>>
>>>
>>> ------------------------------------------------------------------------------
>>> _______________________________________________
>>> pdl-devel mailing list
>>> pdl-devel@lists.sourceforge.net
>>> https://lists.sourceforge.net/lists/listinfo/pdl-devel
>>>
>>
>>
>
> ------------------------------------------------------------------------------
> _______________________________________________
> pdl-devel mailing list
> pdl-devel@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/pdl-devel
>
>
>
------------------------------------------------------------------------------
_______________________________________________
pdl-devel mailing list
pdl-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/pdl-devel
