Thanks Devon for your suggestions. mS is a more or less direct
translation of Twan v Laarhovens mulSquareAt. Though for an accurate
explanation of mS I think it's better to read the blog Tracy Harms
pointed at. :-)
My approach to this problem was mainly translating the Haskell into J.
Reading the whole story with mathematical 'eyes' will probably lead
expert J-ers to a better J-solution.
A 'condensed' mS with your suggestions can be:
(one step away from tacit?)
mS2=: 4 : 0
*/ x([: +/ [*]{y,0$~[:+/(#y)<:])"1 (IJ+<.-:#y)-/i:<.-:#x
)
But it doesn't bring a (significant) speed advantage.
Is there a better way for the filling zeros? They are needed for out of
range indexing of y (again a direct interpretation/translation of the
Haskell source).
NoKPs =: 4 : 0
IJ=: y
ep=. (0x&(,+,~)^: x, 1) * _1&*^:(<1+x) _1x^x
pa=. 1 1 0 1 1x & (+//.@(*/)) &.>^: (<x+1) <1
(2x^x)%~ +/ ep * > pa mS2 &.> |. (*1 _1$~#)&.>pa
)
ep is my interpretation of the sign alternating multinomials used in the
Haskell source:
mn =: (!@(+/)% */@:!)`0: @. (+./@(0&>))
((_1x ^...@-) * [: mn"1 (,. |.)@:i.) 5+1x
_1 5 _10 10 _5 1
((_1x ^...@-) * [: mn"1 (,. |.)@:i.) 6+1x
1 _6 15 _20 15 _6 1
This will give some speed advantage for 100 NoKPs 4 4:
(0x&(,+,~)^: 6, 1) * _1&*^:(<1+6) _1x^6
1 _6 15 _20 15 _6 1
powers_of_a: pa are boxed polynomials and
|. (*1 _1$~#)&.>pa
powers_of_b: their 'counterpolynomials' so to say.
BTW: Twan v Laarhovens array-based solution runs compiled in ~125 ms on
my computer compared to:
ts '100 NoKPs 4 4'
0.363078 7499328
=@@i
Devon McCormick schreef:
> Just a few notes based on only a cursory examination of this code: it does
> not appear to use arrays well. Working on higher-level arrays rather than
> scalars simplifies your code and potentially speeds it up.
>
> For instance, in the function "mS" below, nearly every line has a duplicate,
> doing the same thing on slightly different data. Also,it's hard to tell
> from a single example, but in "NoKPs" it looks like you're breaking apart
> your two-element vector argument into two global scalars - not a good
> practice for a number of reasons. This breaking the one thing - a vector -
> into two separate scalars - is at the root of the duplicated code in "mS".
>
> To clean this up,here's an example of how this first part of "mS":
>
> ny=. #y
> ax=. i:<.-: #x
> hy=. <.-: ny
> xi=: hy + I - ax
> xj=: hy + J - ax
>
> could be simplified to
>
> xij=. (IJ+<.-:#y)-/i:<.-:#x
>
> (where I re-combine the two global scalars back into a vector ( IJ=: I,J) -
> only an interim solution as there's no need for a global at all, especially
> if it's only ever used once.)
>
> This results in a single, two-row matrix "xij" - instead of the two vectors
> "xi" and "xj" - and eliminates five temporary variables that don't get used
> again.
> The two pairs of duplicate lines after this could also be similarly
> simplified.
>
> Basically, it looks like "mS" is really a one- or two-line function. Once
> you've saved all those lines of code, it wouldn't hurt to give back some of
> the typing saved with a descriptive comment and perhaps a more meaningful
> name for the function; the function name doesn't even need to be longer if
> you explain its meaning in a comment.
>
> I'd be interested in finding out if reducing the code like this gives any
> performance improvement; it almost certainly would if it allows you to come
> up with a tacit version though that may be too much work in this case.
>
> Good luck,
>
> Devon
>
> On Sun, Jan 4, 2009 at 3:12 AM, Aai <[email protected]> wrote:
>
>
>> Cleaning up a bit. And reusing resources gives a small improvement:
>>
>>
>> mS=: 4 : 0
>> ny=. #y
>> ax=. i:<.-: #x
>> hy=. <.-: ny
>> xi=: hy + I - ax
>> xj=: hy + J - ax
>> p1=. +/ x * xi{y,0$~+/ny<: xi
>> p2=. +/ x * xj{y,0$~+/ny<: xj
>> p1*p2
>> )
>>
>> NoKPs =: 4 : 0
>> 'I J'=: y
>> ep=. ((0x&(,+,~)^: (]`1:))*(_1x&*^:(<@(1&+)`(_1x&^)))) x
>> pa=. 1 1 0 1 1 & (+//.@(*/)) &.>^: (<x+1) <1
>> (2x^x)%~ +/ ep * > pa mS &.> |. (*1 _1$~#)&.> pa
>> )
>>
>>
>> Aai schreef:
>>
>>> Update:
>>>
>>> ah=: 1 1 0 1 1x
>>> bh=: 1 _1 0 _1 1x
>>>
>>> ppr=: +//.@(*/)
>>>
>>> mulSqAt=: 4 : 0
>>> nx=. #x
>>> ax=. i:<.-: nx
>>> ay=. <.-: #y
>>> aymi=. ay + I - ax
>>> aymj=. ay + J - ax
>>> nuli=. 0$~+/-.(#y)>aymi
>>> nulj=. 0$~+/-.(#y)>aymj
>>> p1=. +/ x * (aymi){y,nuli
>>> p2=. +/ x * (aymj){y,nulj
>>> p1*p2
>>> )
>>>
>>> pathstr =: 4 : 0
>>> 'I J'=: y
>>> ep=: ((0x&(,+,~)^: (]`1:))*(_1x&*^:(<@(1&+)`(_1x&^)))) x
>>> tp=: (ah&ppr&.>^: (<x+1) <1) mulSqAt &.> |. bh&ppr&.>^: (<x+1) <1
>>> (2x^x)%~ +/ ep * > tp
>>> )
>>>
>>> 100 pathstr 4 4
>>>
>>>
>> 2422219241769802380469882122062019059350760968380804461263234408581143863208781993964800
>>
>>> ts '100 pathstr 4 4'
>>> 0.535013 10054144
>>>
>>> Time to speed up things!
>>>
>>> =@@i
>>>
>>>
>>> Aai schreef:
>>>
>>>
>>>> Trying to translate some parts of pathstensor from part IV:
>>>>
>>>> The part with 'multinomial' * sign:
>>>> multinomial [na,nb] * (-1)^nb
>>>>
>>>> can be interpreted as:
>>>>
>>>> mnpd=: (0&(,+,~)^: (]`1:))*(_1&*^:(<@(1&+)`(_1&^)))
>>>>
>>>> (not really a first shot for) emulating mulArray and powers_of:
>>>>
>>>> a=: _2]\ _2 1 _1 1 0 0 1 1 2 1
>>>> b=: _2]\ _2 1 _1 _1 0 0 1 _1 2 1
>>>>
>>>> mularr=: [:({...@{.,+/@{:)@|:/.({.,{:)@(+,*)"1"1 _
>>>>
>>>> powers_of=: 4 : 0
>>>> if. 0=y do. 0 1 return. end.
>>>> if. 1=y do. x return. end.
>>>> x&mularr^: (y-1) x
>>>> )
>>>>
>>>> Some examples:
>>>>
>>>> a powers_of 0
>>>> 0 1
>>>> a powers_of 1
>>>> _2 1
>>>> _1 1
>>>> 0 0
>>>> 1 1
>>>> 2 1
>>>> a powers_of 2
>>>> _4 1
>>>> _3 2
>>>> _2 1
>>>> _1 2
>>>> 0 4
>>>> 1 2
>>>> 2 1
>>>> 3 2
>>>> 4 1
>>>>
>>>> However they will be way too slow for pathstensor 100 (4,4). So
>>>> something to chew on.
>>>>
>>>> ... still trying to wrap my brain around 'mulSquareAt arr arr (i,j)' and
>>>> all the other stuff :-)
>>>>
>>>>
>>>>
>>>> =@@i
>>>>
>>>>
>>>>
>>>> Tracy Harms schreef:
>>>>
>>>>
>>>>
>>>>> There is a lot of array-oriented thinking in a series of four blog
>>>>> posts by Twan van Laarhoven. The ideas include tensor product of rings
>>>>> and arrays as polynomials. Although his code is in Haskell, some here
>>>>> may share my interest in what he's written. I won't be surprised if
>>>>> somebody finds something fun to translate into J, or to do differently
>>>>> in J.
>>>>>
>>>>> Part one is here:
>>>>> http://twan.home.fmf.nl/blog/haskell/Knight1.details
>>>>>
>>>>> The overall blog post series (with all four parts) is here:
>>>>> http://twan.home.fmf.nl/blog/
>>>>>
>>>>> The problem examined is: A knight is placed at the origin of a
>>>>> chessboard that is infinite in all directions. How many ways are there
>>>>> for that knight to reach cell (i,j) in exactly n moves?
>>>>>
>>>>> --
>>>>> Tracy
>>>>> ----------------------------------------------------------------------
>>>>> For information about J forums see http://www.jsoftware.com/forums.htm
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>> ----------------------------------------------------------------------
>>>> For information about J forums see http://www.jsoftware.com/forums.htm
>>>>
>>>>
>>>>
>>>>
>>> ----------------------------------------------------------------------
>>> For information about J forums see http://www.jsoftware.com/forums.htm
>>>
>>>
>>>
>>>
>> ----------------------------------------------------------------------
>> For information about J forums see http://www.jsoftware.com/forums.htm
>>
>>
>
>
>
>
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
