Tích chập, Biến đổi Fourier nhanh và đa thức (2022)

 (alvarorevuelta.com)
1 điểm bởi GN⁺ 2024-07-02 | 1 bình luận | Chia sẻ qua WhatsApp
  • Khi khai triển các đa thức bậc lớn theo cách học ở phổ thông, phải nhân mọi cặp hạng tử nên chi phí O(n²) nhanh chóng trở thành nút thắt cổ chai
  • Phép nhân các vector hệ số của đa thức tương đương với tích chập của tín hiệu rời rạc, và kết quả của [2, 3, 4] và [5, 6, 7] là [10, 27, 52, 45, 28]
  • DFT đưa tín hiệu rời rạc sang miền tần số, còn FFT tính cùng phép biến đổi đó trong O(n log n), tạo ra khác biệt rõ rệt với đầu vào lớn
  • Tích chập trong miền thời gian có thể chuyển thành phép nhân từng phần tử trong miền tần số, nên có thể xử lý phép nhân đa thức nhanh hơn bằng cách biến đổi bằng FFT, nhân, rồi đưa ngược lại bằng IFFT
  • Với bậc nhỏ, chi phí đi-về của FFT/IFFT có thể triệt tiêu lợi ích, nhưng khi bậc tăng lên thì cách dùng FFT hiệu quả hơn

Vì sao phép nhân đa thức trở nên chậm

  • Đa thức P(x) được biểu diễn bằng cách cộng các hạng tử là hệ số a_k nhân với lũy thừa của biến x
    • Ví dụ P(x)=5x²+2x+9 là một đa thức bậc 2
    • Vector hệ số có thể được viết là [5, 2, 9] hoặc [9, 2, 5] tùy theo quy ước ký hiệu
  • Phép cộng và phép trừ tương đối đơn giản vì chỉ cần cộng hoặc trừ các hạng tử cùng bậc
    • Trong Python, có thể duyệt từng hệ số bằng zip(p, q) rồi tính a + b hoặc a - b
    • Nếu bậc khác nhau, có thể dùng zip_longest
  • Phép nhân đòi hỏi phải nhân từng hạng tử với nhau rồi cộng lại các hạng tử cùng bậc, nên lượng tính toán tăng mạnh
    • Kết quả của (2x²+3x+4) × (5x²+6x+7) là 10x⁴+27x³+52x²+45x+28
    • Độ phức tạp của cách này là O(n²), và số phép nhân cần thực hiện tăng theo bậc đa thức

Vector hệ số và tích chập

  • Trong miền rời rạc, tích chập của hai tín hiệu p và q được định nghĩa là y[n]=Σ p[k]·q[n-k]
  • Cách tính là lật ngược q, rồi dịch nó trên p từ trái sang phải và cộng các tích của những phần tử đang chồng lên nhau
  • Các tín hiệu ví dụ như sau
    • p = [2, 3, 4]
    • q = [5, 6, 7]
  • Khi lật và dịch q, từng hệ số đầu ra được tạo theo thứ tự sau
    • 2×5 = 10
    • 2×6 + 3×5 = 27
    • 2×7 + 3×6 + 4×5 = 52
    • 3×7 + 4×6 = 45
    • 4×7 = 28
  • Kết quả tích chập là y = [10, 27, 52, 45, 28]
    • Nó trùng với các hệ số của 10x⁴+27x³+52x²+45x+28 thu được từ phép nhân đa thức
    • Vì vậy, phép nhân đa thức có thể được xem là tích chập của các vector hệ số

Biến đổi Fourier và FFT

  • Biến đổi Fourier chuyển tín hiệu từ miền thời gian sang miền tần số
    • Ở góc nhìn thời gian, tín hiệu được xem như giá trị tại từng thời điểm cụ thể
    • Ở góc nhìn tần số, tín hiệu được diễn giải như tổng của các dao động ở những tần số khác nhau
  • Các tần số dao động được biểu diễn bằng sin và cos, mỗi thành phần có hệ số và pha riêng
  • Khi áp dụng FFT cho sóng sin thuần 5Hz, trong miền tần số nó hiện ra như một delta tại vị trí 5Hz
    • Điều này cho thấy sóng sin trong miền thời gian có thể được biểu diễn chỉ bằng một thành phần sin 5Hz
  • Các thuật ngữ liên quan được phân biệt như sau
    • Fourier Transform(FT): biến đổi Fourier được định nghĩa trên miền liên tục
    • Discrete Fourier Transform(DFT): biến đổi Fourier được định nghĩa cho tín hiệu rời rạc
    • Fast Fourier Transform(FFT): thuật toán tính DFT trong O(n log n) thay vì O(n²)
  • DFT biến tín hiệu thời gian rời rạc x[n] thành X[k] trong miền tần số
    • Mỗi X[k] được tính bằng cách nhân các mẫu đầu vào với số phức biểu diễn một tần số cụ thể rồi cộng lại

Chuyển thành phép nhân trong miền tần số

  • Ưu điểm cốt lõi của DFT và miền tần số là có thể biến tích chập thành phép nhân từng phần tử
    • Tích chập hai tín hiệu trong miền thời gian tương đương với việc nhân hai tín hiệu trong miền tần số
    • Phép nhân có thể được tính nhanh hơn tích chập
  • Quy trình thực hiện phép nhân đa thức nhanh gồm các bước sau
    • Biến đổi đa thức sang miền tần số bằng FFT: O(n log n)
    • Nhân từng phần tử trong miền tần số: O(n)
    • Biến đổi kết quả ngược lại miền thời gian bằng IFFT: O(n log n)
  • Tổng thể, dùng FFT cho phép thực hiện phép nhân đa thức với độ phức tạp O(n log n)
  • Với các đa thức lớn, cách này nhanh hơn phép nhân O(n²) kiểu phổ thông

Cài đặt Python và benchmark

  • multiply_naive dùng vòng lặp lồng nhau để nhân mọi cặp hệ số rồi cộng vào vị trí kết quả i + j
    • Độ dài kết quả là len(p) + len(q) - 1
    • Độ phức tạp là O(n²)
  • multiply_fft thực hiện phép nhân hệ số dựa trên FFT/IFFT
    • Nó tính độ dài là lũy thừa của 2 đủ để chứa kết quả có kích thước ít nhất len(p) + len(q) - 1
    • Dùng np.pad để đệm hai đầu vào
    • Nhân từng phần tử sau khi biến đổi bằng np.fft.fft
    • Dùng np.fft.ifft để biến đổi ngược, rồi làm tròn phần thực để chuyển thành hệ số nguyên
  • Với đầu vào ví dụ p = [2, 3, 4], q = [5, 6, 7], cả hai cách đều trả về [10, 27, 52, 45, 28]
  • Trong benchmark, thay vì multiply_naive, bài viết so sánh cách FFT với multiply_convolve dùng np.convolve
    • Lý do là multiply_naive dùng vòng lặp Python nên chậm, khó so sánh trực tiếp với cách FFT dùng np.fft.fft
    • np.convolve thực hiện cùng phép toán bằng mã C mức thấp
  • Bậc đa thức được tăng trong khoảng range(1, 30000, 1000), và ở mỗi bậc, hai đa thức được tạo với các hệ số ngẫu nhiên từ 1 đến 999999
    • Mỗi cách được đo thời gian trung bình với n_runs = 5
    • Ở bậc thấp, chi phí biến đổi qua lại FFT/IFFT có thể khiến cách FFT không có lợi thế
    • Khi bậc tăng lên, cách dùng FFT cho kết quả hiệu quả hơn nhiều

Bài viết liên quan

1 bình luận

 
GN⁺ 2024-07-02
Ý kiến trên Hacker News

  • Điều luôn khiến tôi khó chịu trong các phần giải thích kiểu này là người ta thường quên mất sai số số học
    Không thể chỉ trừu tượng hóa phép nhân hệ số thành “thời gian hằng số”. Nếu làm vậy thì ngay từ đầu trừu tượng hóa toàn bộ phép nhân cũng chẳng khác gì. Nếu xét đến độ chính xác số học thì nó gần với O(n (log n)^3) hơn [1]
    [1]: http://numbers.computation.free.fr/Constants/Algorithms/fft....

    • Cận sai số được trích trong bài đó quá bi quan. Phiên bản mới nhất của Knuth đã có cận đúng, vì tôi đã báo cho ông ấy
    • Sẽ thật tốt nếu có thể dùng phép toán dựa trên quaternion trong bài của OP để giảm hoặc thậm chí loại bỏ sai số phép nhân [1],[2],[3]
      [1] One-Dimensional Quaternion Discrete Fourier Transform and an Approach to Its Fast Computation:
      https://www.mdpi.com/2079-9292/12/24/4974
      [2] Convolution Theorems for Quaternion Fourier Transform: Properties and Applications:
      https://onlinelibrary.wiley.com/doi/10.1155/2013/162769
      [3] On the Matrix Form of the Quaternion Fourier Transform and Quaternion Convolution:
      https://arxiv.org/abs/2307.01836
    • Nếu hệ số là số nguyên thì có thể dùng NTT với một modulus đủ lớn để thu được kết quả chính xác, và đặc biệt trên phần cứng thì thời gian nhân thậm chí còn có thể nhanh hơn
    • Đó cũng là lý do người ta phân biệt khoa học máy tính với kỹ thuật phần mềm :)

  • Có thể dùng cách này để nhân các số rất dài. Điểm cốt lõi là phép nhân đa thức tương đương với phép nhân số dài thông thường khi chưa xử lý phần nhớ (carry)
    Ví dụ, nếu có một số 1000 chữ số thì lấy từng chữ số làm hệ số của một đa thức có 1000 phần tử. Sau đó có thể nhân các đa thức này bằng cách FFT được mô tả trong bài. Để đổi kết quả trở lại thành số, cần xử lý carry. Nếu một phần tử lớn hơn 10 thì chuyển phần vượt quá sang chữ số kế tiếp, rồi biến các hệ số thành chữ số
    Ý tưởng cơ bản là vậy, nhưng có những điểm tinh tế liên quan đến độ chính xác cần thiết cho carry và việc bảo đảm rằng làm tròn kết quả FFT về số nguyên gần nhất vẫn cho ra đáp án đúng. Đây là cách GMP, thư viện tiêu biểu trong lĩnh vực này, thực hiện phép nhân số lớn

    • Đúng như bạn nói, số thập phân có thể biểu diễn thành đa thức tại x=10 nên khá dễ hiểu. Ví dụ 983 = 9x^2 + 8x + 3, tức là [9, 8, 3]
      Tôi tò mò là con số phải lớn đến mức nào thì cách này mới thực sự có ý nghĩa, và nó được dùng ở đâu

  • Nếu bạn chưa xem thì nên xem video này
    https://youtu.be/h7apO7q16V0?si=bmgUEMTQSqU3flIv
    Video này suy ra thuật toán FFT cho phép nhân đa thức rất xuất sắc. Tôi cứ khoảng 6 tháng lại xem lại một lần

  • Tính chất “tích chập trở thành phép nhân theo từng điểm” của FFT cũng đúng trong mọi nhóm nhân cyclic. Một cách suy diễn mang tính đại số hơn có ở https://www.sciencedirect.com/science/article/pii/S002200007...
    Đôi khi người ta gọi đây là “harmonic FFT”, và cũng có cả các FFT không harmonic: [LCH14] “additive NTT” trên GF(2^n), [HLP24] circle FFT trên đường tròn đơn vị X^2+Y^2=1 của trường hữu hạn, và [BCKL21] ecfft trên chuỗi isogeny của đường cong elliptic
    [LCH14]: https://arxiv.org/abs/1404.3458
    [HLP24]: https://eprint.iacr.org/2024/278
    [BCKL21]: https://arxiv.org/pdf/2107.08473

  • Ai là người đầu tiên đề xuất dùng FFT cho phép nhân đa thức nhanh hơn?
    Gần đây tôi tò mò nên thử tìm, dù không giỏi lần theo trích dẫn lắm nhưng đã lần ngược được tới bài báo năm 1995 của David Eppstein [0]. Ở đó, nó được dùng để giải hiệu quả bài toán tổng con sau các cập nhật tăng dần. Chắc hẳn trong TAOCP của Knuth còn có đề cập sớm hơn
    Việc phép nhân đa thức bằng FFT có thể giải bài toán subset sum chính xác có lặp trong thời gian cận mũ cũng khá gây sốc [1]. Điều quan trọng là thuật toán này là O(N log N), trong đó N không phải kích thước tập mà là phần tử lớn nhất, nên đây không phải kiểu phản ví dụ cho P ≠ NP
    [0] https://escholarship.org/content/qt6sd695gn/qt6sd695gn.pdf
    [1] https://x.com/festivitymn/status/1788362552998580473?s=46&t=...

  • Tôi nghĩ mọi học máy đều là việc giải các phương trình tích chập
    Bài báo này bàn về điều đó trong ngữ cảnh học tăng cường https://arxiv.org/abs/1712.06115, nhưng phần lớn các cách tiếp cận đều khớp với mô thức này

    • Về cơ bản chẳng phải là phương pháp kernel sao?

  • Tôi vừa triển khai một thuật toán dùng FFT để tính tích vô hướng trên một tập lớn các chuỗi con của chuỗi thời gian (matrix profile). Độ dài chuỗi thời gian n có thể lên tới hàng trăm triệu
    Nhờ tính tích chập nhanh bằng FFT, thời gian tính toán giảm từ O(n) xuống O(log n), và ở quy mô này mức cải thiện tốc độ là rất lớn. Nếu còn dùng GPU thì còn nhanh hơn nữa, kiểu xử lý 10 triệu điểm dữ liệu trong 0,1 giây trên laptop

  • Có vẻ “mẹo” cốt lõi của phép toán này là nhận ra điều sau:

    Nói cách khác, thực hiện tích chập của hai tín hiệu trong miền thời gian tương đương với nhân hai tín hiệu trong miền tần số.
    Đây là một bài viết hay vì đã chia một ý tưởng phức tạp thành những bước nhỏ hơn nhiều, đủ để cả người dở toán như tôi cũng có thể phần nào hiểu được. Nhưng có phải tôi đã bỏ lỡ một bước trung gian nào đó không? Hay tác giả cố tình để đó như một bài tập cho người đọc tự tìm hiểu? Đến đoạn đó tôi đã phải huy động hết sức mình về mặt toán học rồi, nên hơi có cảm giác kiểu “rồi giờ hãy vẽ nốt phần còn lại của con cú chết tiệt”. Có ai khác cũng thấy vậy không? Dù sao thì bài viết thực sự rất hay

    • Không biết điều này có giúp ích không: phép nhân hai đa thức mà ta học ở trường thực ra chính là tích chập
      Có tính chất rằng “thực hiện tích chập của hai tín hiệu trong miền thời gian tương đương với nhân hai tín hiệu trong miền tần số”, và FFT cho phép chuyển từ miền thời gian sang miền tần số. Vì vậy ta đưa các đa thức sang miền tần số bằng FFT, rồi trong miền đó chỉ cần nhân thôi. Việc này nhanh hơn tích chập. Không rõ như vậy đã làm sáng tỏ bước bị thiếu chưa; nếu vẫn còn phần nào thiếu, tôi có thể cập nhật bài viết

  • Vậy thì phân tích thừa số nguyên có phải là deconvolution rời rạc không? Nếu đặt biểu diễn FFT, tức phép nhân theo từng điểm, cạnh tableax, tức phép nhân dài thông thường/cộng dồn carry, thì sự phá vỡ đối xứng đó có cung cấp đủ thông tin cho một thuật toán nhanh hay không, tôi thấy khá tò mò

  • Tất nhiên, phép nhân đa thức kiểu ngây thơ sẽ chậm khi bậc đa thức tăng cao. Nhưng trên thực tế, khi nào ta thật sự cần xử lý hai đa thức bậc 100?
    Vì lý do đó, tôi từng có cảm giác rằng các hệ thống đại số máy tính không dùng cách này

    • Các hệ thống đại số máy tính, ví dụ như chebfun của Matlab, chuyển các hàm tùy ý thành đa thức bậc trên 100 để tìm nghiệm, cực trị, v.v. dễ hơn
    • Trong sửa lỗi và xử lý tín hiệu thì điều này rất phổ biến
      https://www.youtube.com/watch?v=CcZf_7Fb4Us
      https://en.wikipedia.org/wiki/Reed%E2%80%93Solomon_error_cor... là một ví dụ
    • Tôi từng tạo một chương trình[1] vì muốn dịch ngược tham số checksum CRC của các tệp lớn, bằng cách biến tệp thành đa thức GF(2) bậc hàng triệu và tính ước chung lớn nhất. Nếu không có phép nhân dựa trên FFT thì việc này không thể làm được trong thời gian hợp lý
      [1]: https://github.com/8051enthusiast/delsum
    • Góc nhìn tích chập này cùng với các kernel GPU nhanh cho FFT đã được dùng trong một số mô hình không gian trạng thái trước thời Mamba để mô hình hóa chuỗi dài; ở đây đa thức chính là chuỗi đầu vào
      Các bài blog của Hazy Research giai đoạn 2020~2023 có khá nhiều thông tin về hướng tiếp cận này
    • Xem https://news.ycombinator.com/item?id=40306339
      "(...) trong nghiên cứu vật lý, tôi từng xử lý các biểu thức dài gần 1 terabyte với hơn 100 triệu hạng tử"