[C++] BOJ #1067 - Moving

This problem is classified as Platinum II difficulty.

To solve this problem, you need to understand convolution codes. Convolution codes are broadly categorized into non-cyclic and cyclic codes. This problem is related to cyclic convolution codes. Though convolution codes are widely used in communication, they are also an important concept in understanding systems in general.

You are given two sequences X and Y, each containing N numbers. The goal is to create a cyclic convolution code of these sequences and find the maximum value from the result.

---

To do this, you need to understand what a convolution code is. Let’s consider two sequences {1, 3, 2, 4} and {1, 2, 3, 4}. If you write these two sequences using convolution code (non-cyclic convolution), it proceeds as follows:

Acyclic Convolution Code

You fix the first sequence and reverse the second sequence, sliding it one element at a time while computing dot products. For instance, multiplying corresponding elements and summing them gives:
• First position: \(1 \times 1 = 1\)
• Second: \(1 \times 2 + 3 \times 1 = 5\)
• and so on…

In fact, convolution is similar to the multiplication of multi-digit numbers. Think of multiplying 4231 by 4321 — the process is quite analogous:
• Units digit: \(1 \times 1 = 1\)
• Tens digit: \(1 \times 2 + 3 \times 1 = 5\)
• Eventually, this gives a result structured like 21000 + 1100 + 50 + 1

When you compute the convolution of two sequences of length 4, the result is a sequence of 7 elements. This process effectively mirrors the multiplication of two numbers.

 

convolution

Using FFT

To solve this problem efficiently, Fast Fourier Transform (FFT) is used. If the numbers in the sequences are large or the number of elements is substantial, rounding errors may arise, preventing exact results. However, since this problem involves numbers below 100 and around 60,000 elements, the precision should be manageable. Using double-precision floating points should ensure accurate results up to hundreds of billions. The maximum product can be up to \(6 \times 10^6 \times 6 \times 10^6 = 36 \times 10^{12}\), so rounding errors might occur in edge cases.

Using FFT reduces the time complexity to \(O(N \log N)\) for both FFT and inverse FFT (iFFT). Even if the number of data points is doubled for cyclic convolution, it is still much faster than the naive \(O(N^2)\) convolution approach, which would time out for data size around 60,000.

Steps for Cyclic Convolution Using FFT

To compute cyclic convolution using FFT:
1. Compute FFT of both sequences:
\(X(z) = FFT(x(k)),\quad Y(z) = FFT(y(k))\)
2. Multiply the FFT results element-wise:
\(Z(z) = X(z) \cdot Y(z)\)
3. Apply inverse FFT:
\(z(k) = iFFT(Z(z))\)

One user reports that their implementation using this method completed in 52ms, though some achieved 16ms. This suggests there may be room for optimization in the FFT implementation.

The FFT implementation they used was borrowed from an external source. Apart from that, the rest of the implementation follows the above formula and is relatively straightforward.

//----------------------------------------------------------
//    baekjoon #1067 - Moving
//        - by Aubrey Choi
//        - created at 2020-02-13
//----------------------------------------------------------
#include <stdio.h>
#include <complex>
#define    MAX    (1<<17)

std::complex<double> ca[MAX], cb[MAX];
void fft(std::complex<double> a[], int n, bool inv)
{
    int i, j, k, bit;
    for(i=1,j=0;i<n;i++) 
    {
        for(bit=n>>1;!((j^=bit)&bit);bit>>=1);
        if(i<j) swap(a[i], a[j]);
    }
    for(i=1;i<n;i<<=1) 
    {
        double x = (inv ? 1 : -1) * M_PI / i;
        std::complex<double> w={cos(x),sin(x)}, th, tmp;
        for(j=0;j<n;j+=i<<1) for(k=0,th=1;k<i;k++)
        {
            tmp = a[i+j+k]*th;
            a[i+j+k] = a[j+k]-tmp;
            a[j+k] += tmp;
            th *= w;
        }
    }
    if(inv) for(i=0;i<n;i++) a[i]/=n;
}
void cycleconv(int a[], const int b[], int n)
{
    int i, k;
    for(k=1;n>k;k<<=1); k<<=1;
    for(i=0;i<n;i++) ca[i]=a[i];
    for(i=1,cb[0]=b[0];i<n;i++) cb[k-i]=cb[n-i]=b[i];
    fft(ca, k, false); fft(cb, k, false);
    for(i=0;i<k;i++) ca[i]*=cb[i];
    fft(ca, k, true);
    for(i=0;i<n;i++) a[i]=(int)(ca[i].real()+0.5);
}
int main()
{
    int n, i, a[60000], b[60000], max;
    scanf("%d",&n);
    for(i=0;i<n;i++) scanf("%d",a+i);
    for(i=0;i<n;i++) scanf("%d",b+i);
    cycleconv(a, b, n);
    for(i=1,max=a[0];i<n;i++) max=(max<a[i])?a[i]:max;
    printf("%d\n", max);
}