# [Rd] fft fails for lengths 392, 588, 968, 980 .... (PR#1429)

**brahm@alum.mit.edu
**
brahm@alum.mit.edu

*Thu, 4 Apr 2002 17:55:04 +0200 (MET DST)*

Gabriel Fricout <Gabriel.FRICOUT@cmm.ensmp.fr> reported a bug in fft(), and
Brian D. Ripley <ripley@stats.ox.ac.uk> narrowed it down to specific values of
N (the data vector length). For N <= 1000, BDR showed the failures occur at:
N = 392 = 7^2 * 2^2 * 2
N = 588 = 7^2 * 2^2 * 3
N = 968 = 11^2 * 2^2 * 2
N = 980 = 7^2 * 2^2 * 5
BDR's test is as follows:
R> for (i in 1:1000) {
R> X <- rnorm(i)
R> XX <- fft(fft(X), inverse=T)/i
R> if(max(Mod(XX-X)) > 1e-10) print(i)
R> }
In src/appl/fft.c, lines 820-822 set "maxf" in subroutine fft_factor. "maxf"
is supposed to be the "maximum factor size", and the algorithm basically grabs
the last non-squared factor and the last squared factor (there are "kt" squared
factors):
maxf = nfac[m_fac-kt-1];
if (kt > 0)
maxf = imax2(nfac[kt-1], maxf);
But the last squared factor (nfac[kt-1]) is not necessarily the largest, since
the squared factors are added to nfac in this order: 4, 3,5,7,..., 2 (see lines
769-791). Thus I believe we need to add these two lines after line 822:
if (kt > 1) maxf = imax2(nfac[kt-2], maxf);
if (kt > 2) maxf = imax2(nfac[kt-3], maxf);
I don't know for sure that the second line is necessary (I saw no error for
N = 1152 = 4^2 * 3^2 * 2^2 * 2). But I did rebuild R-devel (2002-04-03) with
these two lines added, and BDR's test ran without errors. I submit this as the
bug fix for PR#1429.
--
-- David Brahm (brahm@alum.mit.edu)
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._