Sorry - a few days have passed....
Interesting!
As I had put up this idea (using fft multiplication and a special case
for fft "division"), I feel I can
wonder if this helps the Rosetta Code exposition of the virtues of J! I
like my idea, of course,
and it's of potential interest to the J community, but are Rosetta
coders interested in tweaks
to improve performance?
REBoss's contribution is quite different, and very compact, so worth
including.
There's no chance the contribution will reach the length of the Java
one, say, but casual readers
not familiar with J might get the impression that all the code presented
is required.
Just a thought!
Anyway, I had a look at saving some of the continued redefinition of
roots of unit and reuse of "cube", given
the repeated use of fft in the cyclotomic function,by defining a
function fft1 which would be used instead of
fft"1 for multiples of lists of polynomials all of the same extended
length, but it doesn't help much here,
and is of course more complicated, so I won't share it.
However, if the fft multiplication version is to be included in the
Rosetta Code contribution, I think the
following offers some improvement, mainly from the point of view of
exposition. A bit of thinking reveals
that the maximum length of all the polynomials will be that of x^n - 1 .
I believe this block of code is self-contained, assuming addons include
the fftw library. unpad only needs
redefinition here for self-containedness. ctfft can be renamed of
course! Use of M. is a matter of taste,
I suppose; it does save having to define your own cache, as previously
discussed.
NB. Rosetta Coders don't need to worry about fft v fftw - that's for J-ers.
load '~addons/math/fftw/fftw.ijs'
fft =: fftw
ifft =: ifftw
roundreal =: [: (<.) 0.5 + 9&o.
unpad =: {.~ (1 + 0 i:~ 0 = ])
ctfft =: {{ assert.0<y
if. y = 1 do. _1 1 return. end.
'q p'=. __ q: y
if. 1=#q do.
,(y%*/q) {."0 q#1
elseif.2={.q do.
,(y%*/q) {."0 (* 1 _1 $~ #) ctfft */}.q
elseif. 1 e. 1 < p do.
,(y%*/q) {."0 ctfft */q
else.
ctlist =. ctfft "0 }:*/@>,{1,each q NB. ctlist is 2-d
table of polynomial divisors
dividend=. _1,(-y){.1 NB. coefficients of
(x^n) - 1
lg =. >.&.(2&^.) y NB. required length
of each polynomial for fft transforms here
fftdivis=. lg */ @:(fft @:{."1) ctlist NB. fft of product
of divisor polys = product of ffts of divisor polys
fftquot =. fftdivis %~ fft lg{.dividend NB. (fft of
(dividend%divisor) = (fft of dividend) % (fft of divisor)
unpad roundreal ifft"1 fftquot NB. clean up inverse
fft of (dividend&divisor)
end.
}} M.
Also, I noticed a few minor typos in the remarks between "If you take
all the divisors... " and "Task examples:"
Cheers,
Mike
On 02/03/2022 19:57, Raul Miller wrote:
I've updated the rosettacode page to include this implementation. And,
I plan on updating the j wiki polynomial division essay, also, once I
think through a bit how to describe this approach. (Descriptive
phrases and sentences would be welcome. Documentation needs good
perspectives, above all else.)
For the purpose of the rosettacode page, I used a different phrasing of rou:
rou =: {{(,j.) r,(%:0j1),j.+|.}.r=.^o.j.(i.@%&8 % -:) y}} NB.
roots of unity
This seems (to me) to be more concise, and straightforward, than the
purely tacit expression.
(That said, I'm not sure if this change in notation would be a good
fit for the FFT essay.)
Thanks,
--
This email has been checked for viruses by Avast antivirus software.
https://www.avast.com/antivirus
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
