1

I have an FFT and IFFT functions. And i know that

A*B = IFFT(FFT(A)*FFT(B))

Where

FFT(A)FFT(B)=[qw for q,w in zip(A,B)]

But when i input: 10 10 => output: [(0.5+0j), (0.5+0j)] What i`m doing wrong? Here is my code:

from cmath import exp,pi

def FFT(X):
    n = len(X)
    w = exp(-2*pi*1j/n)
    if n > 1:
        X = FFT(X[::2]) + FFT(X[1::2])
        for k in range(n//2):
            xk = X[k]
            X[k] = xk + w**k*X[k+n//2]
            X[k+n//2] = xk - w**k*X[k+n//2]
    return X

def IFFT(X):
    n = len(X)
    w = exp(2*pi*1j/n)
    if n > 1:
        X = IFFT(X[::2]) + IFFT(X[1::2])
        for k in range(n//2):
            xk = X[k]
            X[k] = (xk + w**k*X[k+n//2])/n
            X[k+n//2] = ((xk - w**k*X[k+n//2]))/n
    return X

s = input().split()
a1 = s[0]
b1 = s[1]
a = [int(x) for x in a1]
b = [int(x) for x in b1]
if (len(a)>len(b)):
   n = len(a)
   b.reverse()
   while (len(b)<n):
       b.extend([0])
   b.reverse()
else:
    n = len(b)
    a.reverse()
    while(len(a)<n):
        a.extend([0])
    a.reverse()
c = [q*w for q,w in zip(a,b)]

print (IFFT(c))
fferri
  • 18,285
  • 5
  • 46
  • 95
Reodont
  • 329
  • 1
  • 14
  • 4
    I am not sure I got your point. But A*B = IFFT(FFT(A)*FFT(B)) does not seem right to me, it rather goes as A \ast B = IFFT(FFT(A)*FFT(B)), where \ast is a convolution product. Check out [this wikipedia page](http://en.wikipedia.org/wiki/Convolution_theorem). – Tony Apr 26 '15 at 14:05
  • I thought that FFT(A*B) = FFT(A)*FFT(B). Thats mean that IFFT (FFT(A*B))=IFFT(FFT(A)*FFT(B)) and this is A*B = IFFT(FFT(A)*FFT(B)) – Reodont Apr 26 '15 at 14:17
  • 1
    I don't think FFT(A*B) = FFT(A)*FFT(B) is correct. You may want to check out the previous link. – Tony Apr 26 '15 at 14:20
  • It is first formula on this wikipedia page. Isn't it? – Reodont Apr 26 '15 at 14:23
  • 1
    No, the "star" represents a convolution product, while the dot represents the usual product. See [here](http://en.wikipedia.org/wiki/Convolution) for the definition of the convolution product. – Tony Apr 26 '15 at 14:26
  • Ok.. So if want to compute A*B i need to compute IFFT of convolution product of FFT(A) and FFT(B). That`s right? – Reodont Apr 26 '15 at 14:38
  • Yes, it is usually used the other way around: to compute the convolution product, you go through FFT to do it faster. – Tony Apr 26 '15 at 14:40
  • You are welcome. But I'm curious, what are you trying to do here ? Test your (I)FFT ? – Tony Apr 26 '15 at 14:52
  • Trying to multiplicate big integers as fast as coomputer can using FFT. – Reodont Apr 26 '15 at 15:02
  • Ok. Well it does not seem to be a very good plan to me if you are looking for efficiency. I don't think you can beat the plain multiplication with FFT (or anything else), it will always involve way more operations ... Nevertheless, you may consider [numpy arrays](http://docs.scipy.org/doc/numpy/reference/arrays.html) to do to the job for you in a very efficient way. Although it will use more or less the same strategy, numpy being written in C, it will be much faster here. – Tony Apr 26 '15 at 15:10
  • @Tony, Sorry for my ignorance but i have no idea how to write convoluation function. Do you know where i can find some code, without using of non-standart libraries? – Reodont Apr 26 '15 at 16:15
  • Well, as I said before, usually people use FFT to compute convolution products. I am not aware of any library providing that, perhaps Numpy. What is important here is that doing the plain multiplication (i.e. doing it term by term) is the best you can do ... Any other way would use a lot more operations ... If it is too slow for you purpose, you can make it much faster using library written in lower level languages (such as Numpy, written in C), which are much faster when it comes to loops. – Tony Apr 26 '15 at 16:31
  • @Tony FFT runs in `O(N*log(N))` time. Whereas the naive (term-by-term) method is `O(N^2)`. – Mysticial Apr 26 '15 at 16:59
  • Multiplying two arrays as he is doing is an O(N) operation (one multiplication per term). – Tony Apr 26 '15 at 17:00
  • @Wespocan Can you clarify whether you are trying to perform a convolution vs. an element-wise product of two arrays? – Mysticial Apr 26 '15 at 17:04
  • @Tony Rereading the comments, the OP is trying to multiply large integers. That means a convolution. So yes FFT makes sense. I assume this is an academic exercise since it will be difficult to do better than the built-in multiply operator. Python uses GMP for bignums which already use some form of FFT which is optimized for the purpose. – Mysticial Apr 26 '15 at 17:09
  • Hum, that is not what I understand from the comments. But indeed, assuming this is an academic exercise, it may be a convolution product, although I am not sure the OP is aware of that. If it is, I agree, the FFT makes sense. – Tony Apr 26 '15 at 17:15
  • @Mysticial: Nitpick: Python (at least CPython) doesn't use GMP for its int (long in Python 2) implementation, though that has been proposed a couple of times in the past. It implements Karatsuba multiplication, but nothing more exciting than that. – Mark Dickinson Apr 26 '15 at 18:05
  • 1
    Multiply big numbers is convolution, as you multiply every digit in one input by every digit in the other input when multiplying long hand, then do a linear summation of all n*m sub-products. To do linear convolution by FFT, don't forget to zero-pad to at least n+m-1 in length, or you will end up with a circular convolution instead. – hotpaw2 Apr 26 '15 at 18:39

1 Answers1

1

You should actually apply the FFT to your input sequences a and b. There seems to be an oversight after the extensive padding procedure.

c = IFFT( [ q*w for q,w in zip( FFT(a), FFT(b) ) ] )

What you are doing can be interpreted, in one of the many ways, as performing polynomial multiplication of

a[0]+a[1]*z+...+a[n-1]*z^(n-1)  mod (z^n-1) 

with

b[0]+b[1]*z+...+b[n-1]*z^(n-1)  mod (z^n-1) 

by first evaluating at n equidistant points on the unit circle, multiplying the point values and then interpolating the result polynomial from the multiplied values. The implementation and especially the distribution of the 1/n factor to the IFFT are conform to this procedure.

Lutz Lehmann
  • 25,219
  • 2
  • 22
  • 51