[BOJ/백준] 1067 이동

### FFT를 쓰는걸 보고 연습하는 문제 ####
import sys
from cmath import exp,pi
def fft(a):
    N=len(a)
    if N==1:
        return a
    a_even=fft(a[0::2])
    a_odd=fft(a[1::2])
    w_N=[exp(2j*pi*n/N) for n in range(N//2)]
    return [a_even[n] +w_N[n]*a_odd[n] for n in range(N//2)] + [a_even[n]-w_N[n]*a_odd[n] for n in range(N//2)]

def ifft(a):
    N=len(a)
    if N==1:
        return a
    a_even=ifft(a[0::2])
    a_odd=ifft(a[1::2])
    w_N=[exp(-2j*pi*n/N) for n in range(N//2)]
    return [a_even[n] +w_N[n]*a_odd[n] for n in range(N//2)] + [a_even[n]-w_N[n]*a_odd[n] for n in range(N//2)]

M=int(sys.stdin.readline())
N=2*M
even=0
for i in range(18):
    if M==2**i:
        even=-100
        break
    elif N<2**i:
        even=i
        break
A=list(map(int,sys.stdin.readline().split()))
B=list(map(int,sys.stdin.readline().split()))
if even==-100:
    A=A[:]+A[:]
    B=B[-1::-1]+[0]*M
    C=[0]*N
    A_fft=fft(A)
    B_fft=fft(B)
    for i in range(N):
        C[i]=A_fft[i]*B_fft[i]

    C_ifft=ifft(C)
    for k in range(N):
        C_ifft[k]=round(C_ifft[k].real/N)
    max_number=max(C_ifft)
else:
    N_prime=2**i
    N,N_prime=N_prime,N
    A=A[:]+[0]*(N-N_prime//2)
    B=B[-1::-1]+[0]*(N-N_prime)+B[-1::-1]

    C=[0]*N
    A_fft=fft(A)
    B_fft=fft(B)
    for i in range(N):
        C[i]=A_fft[i]*B_fft[i]
    C_ifft=ifft(C)
    for k in range(N):
        C_ifft[k]=round(C_ifft[k].real/N)
    max_number=max(C_ifft)
print(max_number)

casterian.net/archives/297

고속 푸리에 변환 (Fast Fourier Transform, FFT) – The Casterian

 

 

blog.naver.com/kks227/221633584963

고속 푸리에 변환(Fast Fourier Transform) (수정: 2019-09-05)

 

FFT에 대한 자세한 설명은 위 2개의 블로그를 보면 될것 같다.

 

기본적으로 주어진 수가 2의 제곱수일때와 아닐때를 구분해서 FFT를 구현했고, A와 B를 FFT를 통해 DFT를 구해놓은 상태에서, 

 

이 이동의 값은 A,B의 합성곱을 나타낼수 있는데, 이걸 FFT를 통해 구해낸 DFT A와 DFT B의 곱셉으로 A와B의 합성곱의 DFT화 된 C를 구할 수 있다.

 

우리가 원하는 것은 실수부에 있는 것이므로, DFT화 된 C를 IFFT로 고속 역 푸리에변환을 해주고, 이때의 값을 반올림을 한 상태에서 최대값이 되는 것을 구하면 된다.

 

이 문제 자체를 합성곱으로 해서 풀어야한다는 것을 알았지만, 이 FFT를 어떻게 구현하는지에 대해서는 지식이 부족하여, 위 블로그들을 참조하여 구현했다.

 

제 설명보다는 위의 2 블로그를 참조하는 것을 추천합니다.