I made that change during the whole g77-to-gfortran switch several years ago,
when it was hard to get slatec to compile on many platforms. I have no
objection to switching it back. MatrixOps was an attempt to put at least
*something* in place that would compile whenever the core of PDL would — but in
retrospect it was probably a mistake. Where it overlaps with Slatec or GSL,
Slatec or GSL invariable do better.
> On Mar 4, 2015, at 2:04 PM, Karl Glazebrook <karlglazebr...@mac.com> wrote:
>
> Hi Derek
>
> Yes it seems the switch from PDL::Slatec::matinv to PDL::MatrixOps::inv was
> at fault! Clearly the latter does not thread correctly.
>
> Further investigation suggests that for fitpoly1d($x,$y) with 2D $y then it
> fails when
> - $x is absent
> - $x is 2D
> But it WORKS when $x is 1D, with $y 2D which is probably the most common
> case and why it was not caught to now! (That and me not writing a test)
>
> Who made the change? Would anyone mind if we switched it back? The
> alternative is to make inv() thread correctly. I’d kind of like to correct
> this in 2.008. It’s my module :-)
>
> Karl
>
>
>> On 4 Mar 2015, at 8:05 am, Derek Lamb <de...@boulder.swri.edu
>> <mailto:de...@boulder.swri.edu>> wrote:
>>
>> Hi Karl,
>>
>> Maybe if we had a test for threaded fitpoly1d we would have caught this 4.5
>> years ago!
>>
>> Looks like the problem is with the change from PDL::Slatec::matinv to
>> PDL::MatrixOps::inv [8adb0b]
>> <http://sourceforge.net/p/pdl/code/ci/8adb0b9e6e88ded9153563509794ed253d4f50e1/>.
>> Running your example through the debugger I get:
>>
>> DB<34> p inv($C)
>>
>> [
>> [ 0.168 -0.064]
>> [ -0.32 0.16]
>> ]
>>
>> DB<35> p matinv($C)
>>
>> [
>> [
>> [ 0.6 -0.4]
>> [-0.4 0.4]
>> ]
>> [
>> [ 0.6 -0.4]
>> [-0.4 0.4]
>> ]
>> ]
>>
>> Wolfram Alpha suggests that matinv is correct. When running your
>> single-dimension example through, inv and matinv return identical results.
>>
>> My guess is PDL::MatrixOps::inv is not threadable. It certainly isn't
>> advertised as such in the docs, though other subs in MatrixOps are
>> advertised as such.
>>
>> cheers,
>> Derek
>>
>>
>>> On Mar 4, 2015, at 12:59 AM, Karl Glazebrook <karlglazebr...@mac.com
>>> <mailto:karlglazebr...@mac.com>> wrote:
>>>
>>> fitpoly1d() no longer seems to thread over extra dimensions:
>>>
>>> pdl> $x = pdl( [101,103,104,102,109] )
>>>
>>> pdl> p $x
>>> [101 103 104 102 109]
>>> pdl> ($xf, $c) = fitpoly1d($x, 2)
>>>
>>> pdl> p $xf
>>> [100.8 102.3 103.8 105.3 106.8]
>>> pdl>
>>> pdl> $x = pdl( [101,103,104,102,109],[101,103,104,102,109] )
>>>
>>> pdl> p $x
>>>
>>> [
>>> [101 103 104 102 109]
>>> [101 103 104 102 109]
>>> ]
>>>
>>> pdl> ($xf, $c) = fitpoly1d($x, 2)
>>>
>>> pdl> p $xf # WTF?!
>>>
>>> [
>>> [ 53.496 12.576 -28.344 -69.264 -110.184]
>>> [ 53.496 12.576 -28.344 -69.264 -110.184]
>>> ]
>>>
>>>
>>> I took a look at the Polynomial.pm and could no longer figure it out! I
>>> swear it used to work - anyone have any ideas?
>>>
>>> Karl
>>>
>>> ------------------------------------------------------------------------------
>>> Dive into the World of Parallel Programming The Go Parallel Website,
>>> sponsored
>>> by Intel and developed in partnership with Slashdot Media, is your hub for
>>> all
>>> things parallel software development, from weekly thought leadership blogs
>>> to
>>> news, videos, case studies, tutorials and more. Take a look and join the
>>> conversation now.
>>> http://goparallel.sourceforge.net/_______________________________________________
>>>
>>> <http://goparallel.sourceforge.net/_______________________________________________>
>>> pdl-devel mailing list
>>> pdl-devel@lists.sourceforge.net <mailto:pdl-devel@lists.sourceforge.net>
>>> https://lists.sourceforge.net/lists/listinfo/pdl-devel
>>> <https://lists.sourceforge.net/lists/listinfo/pdl-devel>
>>
>
> ------------------------------------------------------------------------------
> Dive into the World of Parallel Programming The Go Parallel Website, sponsored
> by Intel and developed in partnership with Slashdot Media, is your hub for all
> things parallel software development, from weekly thought leadership blogs to
> news, videos, case studies, tutorials and more. Take a look and join the
> conversation now.
> http://goparallel.sourceforge.net/_______________________________________________
> pdl-devel mailing list
> pdl-devel@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/pdl-devel
------------------------------------------------------------------------------
Dive into the World of Parallel Programming The Go Parallel Website, sponsored
by Intel and developed in partnership with Slashdot Media, is your hub for all
things parallel software development, from weekly thought leadership blogs to
news, videos, case studies, tutorials and more. Take a look and join the
conversation now. http://goparallel.sourceforge.net/
_______________________________________________
pdl-devel mailing list
pdl-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/pdl-devel
