Upgrade to Pro — share decks privately, control downloads, hide ads and more …

競技プログラミングのための FFT - 多項式乗算の高速化ことはじめ

競技プログラミングのための FFT - 多項式乗算の高速化ことはじめ

スライド内のリンクを見たい方はこちら:https://hcpc-hokudai.github.io/archive/math_fft_001.pdf

Avatar for tsutaj

tsutaj

June 24, 2018
Tweet

More Decks by tsutaj

Other Decks in Programming

Transcript

  1. ڝٕϓϩάϥϛϯάͷͨΊͷ ..6 ଟ߲ࣜ৐ࢉͷߴ଎Խ͜ͱ͸͡Ί tsutaj (@tsutaj) Hokkaido University M1 June 24,

    2018 tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 1 / 14
  2. ຊεϥΠυͷ໨ඪ ଟ߲ࣜ৐ࢉ p ࣍ଟ߲ࣜ g(x) = ∑p i=0 cixi ͱɺq

    ࣍ଟ߲ࣜ h(x) = ∑q i=0 dixi ͷ৐ࢉΛ ҎԼͷΑ͏ʹఆٛ (g ∗ h)(x) = g(x)h(x) = p+q ∑ k=0 min(k,p) ∑ i=0 cidk−ixk (1) ྫ: g(x) = 2x + 3, h(x) = 3x2 + 4x + 1 (g ∗ h)(x) = (2x + 3)(3x2 + 4x + 1) = 6x3 + 17x2 + 14x + 3 g ΋ h ΋ n ࣍ࣜͰ͋Δͱͨ͠ͱ͖ɺ৐ࢉʹ͔͔Δܭࢉྔ͸ O(n2) ͩ ͕ɺFFT ʹΑͬͯ O(n log n) ʹߴ଎ԽͰ͖Δʂ ຊεϥΠυͷ໨ඪ: FFT Λར༻ͯ͠ଟ߲ࣜ৐ࢉΛߴ଎Խ͢Δ͜ͱ ΄͔ͷԠ༻ͳͲʹ͸৮Ε·ͤΜ tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 2 / 14
  3. جຊઓུ p ࣍ଟ߲ࣜ g(x) ͱ q ࣍ଟ߲ࣜ h(x) ͷੵ ˠ

    p + q ࣍ࣜ m ࣍ଟ߲ࣜ f(x) ʹରͯ͠ɺগͳ͘ͱ΋ m + 1 ݸͷ఺ x0, . . . , xm Ͱ ͷ஋ f(x0), . . . , f(xm) ͕෼͔͍ͬͯΕ͹ɺ͜ΕΒΛશͯ௨ΔΑ͏ͳଟ ߲͕ࣜҰҙʹఆ·Δ m = 1 ௚ઢ ͷ৔߹ɺ2 ఺ͷ஋͕෼͔͍ͬͯΕ͹Ұҙʹఆ·Δ m = 2 ์෺ઢ ͷ৔߹ɺ3 ఺ͷ஋͕෼͔͍ͬͯΕ͹Ұҙʹఆ·Δ جຊઓུ 1 ݩͷؔ਺ΛٻΊΔͷʹे෼ͳ਺ͷ఺ʹ͍ͭͯɺͦͷؔ਺஋ΛٻΊΔ 2 ઃఆͨ͠఺Λશͯ௨ΔΑ͏ͳؔ਺ΛԿΒ͔ͷܗͰٻΊΔ p + q + 1 ݸͷ఺ʹ͓͚Δؔ਺஋ ͕෼͔͍ͬͯΕ͹ɺ(g ∗ h)(x) ͕෼͔ Γͦ͏ʂ ʮ(g ∗ h)(xi ) ΛٻΊΔʯ=ʮg(xi )h(xi ) ΛٻΊΔʯ ͳͷͰܭࢉ͸؆୯ ܭࢉʹ͓͍ͯศརͳ఺ΛఆΊͯɺָʹݩͷؔ਺ (g ∗ h)(x) Λಘ͍ͨ ָ͕Ͱ͖Δ఺ͱ͸Կ͔ʁ tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 3 / 14
  4. ఺ͷબͼํ N := (p + q + 1) ≤ N

    Λຬͨ͢࠷খͷ 2 ͷ΂͖৐ ͱ͠ɺN ݸͷ఺Λ༻ҙ ͢Δ͜ͱΛߟ͑Δ ఺ͷબͼํ ఺ x0, . . . , xN−1 ͱͯ͠બͿͷ͸ɺ1 ͷ N ৐ࠜʂ ζN = exp(2π √ −1 N ) ͱ͓͘ͱ͖ɺxi = (ζN )i ͱ͢Ε͹Α͍ 1 ͷ N ৐ࠜΛ఺ͱͯ͠બͿͱɺԿ͕خ͍͠ͷʁ (ζN )i = (ζN )j ͳΒ͹ɺi = j mod N ͕੒ཱ (ζN )N = 1 ΑΓɺ(ζN )i = (ζN )kN+i (k ∈ Z+) ͳͷͰ ௚ަੑ͕͋Δ N−1 ∑ i=0 ( ζj N )i ( ζk N )i = N−1 ∑ i=0 ( ζj N )i ( ζ−k N )i = N−1 ∑ i=0 (ζN )i(j−k) (2) = { N (j = k mod N) 0 (otherwise) (3) tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 4 / 14
  5. ܎਺ͷ෮ݩ આ໌ͷͨΊ f(x) = (g ∗ h)(x) = ∑ N−1

    j=0 hjxj ͱ͓͘ ઌ΄ͲબΜͩ఺ʹରͯ͠ɺ࣮ࡍʹ f ( ζ0 N ) , . . . , f ( ζi N ) , . . . , f ( ζN−1 N ) ΛٻΊͨͱ͢Δ ͔͜͜Β f(x) ͷ܎਺ h0, . . . , hN−1 Λ෮ݩ͢Δʹ͸ʁ ཭ࢄϑʔϦΤม׵ (Discrete Fourier Transformation, DFT) ଟ߲ࣜ f(x) ʹରͯ͠ɺબΜ֤ͩ఺Ͱͷ஋Λ܎਺ʹ΋ͭଟ߲ࣜ ˆ f(t) Λ࡞Δ (Լઢ෦ͷ f ͷஔ͖׵͑ʹ஫ҙʂ ) ˆ f(t) = N−1 ∑ i=0 f ( ζi N ) ti = N−1 ∑ i=0 N−1 ∑ j=0 hj ( ζi N )j ti (4) = N−1 ∑ j=0 hj N−1 ∑ i=0 ( ζj N t )i (5) tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 5 / 14
  6. ܎਺ͷ෮ݩ DFT ʹΑͬͯಘΒΕͨ ˆ f ͔Β܎਺Λ෮ݩ͍ͨ͠ɾɾɾ ˆ f ( ζ−k

    N ) ΛٻΊͯΈΑ͏ ˆ f ( ζ−k N ) = N−1 ∑ j=0 hj N−1 ∑ i=0 ( ζj N ζ−k N )i = N−1 ∑ j=0 hj N−1 ∑ i=0 (ζN )i(j−k) (6) = N × hk (ζN ͷ௚ަੑΛར༻) (7) ˆ f ( ζ−k N ) ΛٻΊΔ͜ͱʹΑͬͯɺͳΜͱ ݩͷؔ਺ f(x) ͷ k ൪໨ͷ܎ ਺ (Λ N ഒͨ͠΋ͷ) ͕ٻΊΒΕΔʂ ͏Ε͍͠ʂ ˆ f ( ζ0 N ) , . . . , ˆ f ( ζ−(N−1) N ) ΛٻΊͯ΍Ε͹Αͦ͞͏ tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 6 / 14
  7. ཭ࢄϑʔϦΤٯม׵ ཭ࢄϑʔϦΤٯม׵ f(x) ͷ DFT ˆ f(t) Λར༻ͯ͠ɺݩͷؔ਺ f(x) ͸ҎԼͷΑ͏ʹٻΊΒΕΔ

    f(x) = 1 N N−1 ∑ i=0 ˆ f ( ζ−i N ) xi (∵ ˆ f ( ζ−i N ) = N × hi) (8) ཭ࢄϑʔϦΤม׵ (DFT) ˆ f(t) = N−1 ∑ i=0 f ( ζi N ) ti (9) ͱൺֱ͢Δͱɺม਺ x ͱ t ͕ೖΕସΘ͍ͬͯͯɺDFT Ͱ ζi N Ͱ͋ͬͨ ͱ͜Ζ͕ ζ−i N ΁ͱมΘ͍ͬͯΔ ࣜ (8) Λ཭ࢄϑʔϦΤٯม׵ (inverse DFT) ͱݺͿ tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 7 / 14
  8. ͜͜·Ͱͷ·ͱΊ ೖྗ: p ࣍ଟ߲ࣜ g(x) ͱ q ࣍ଟ߲ࣜ h(x) ग़ྗ:

    ৐ࢉͯ͠ಘΒΕΔଟ߲ࣜ f(x) = (g ∗ h)(x) = g(x)h(x) = ∑ N−1 j=0 hjxj ཭ࢄϑʔϦΤม׵ɾٯม׵Λ༻͍ͨղ๏ 1 N ← p + q + 1 ≤ N Λຬͨ͢࠷খͷ 2 ͷ΂͖৐ 2 1 ͷ N ৐ࠜ (ζN )0, . . . , (ζN )N−1 શͯʹର͠ɺf ( ζi N ) Λܭࢉ 3 ٻΊͨؔ਺஋Λ΋ͱʹ ˆ f(t) = ∑ N−1 i=0 f ( ζi N ) ti Λ࡞Δ (DFT) 4 ˆ f ( ζ−k N ) = N × hk Ͱ͋Δ͜ͱΛར༻ͯ͠ɺݩͷଟ߲ࣜΛ෮ݩ (inverse DFT) ࢒Δ՝୊͸ɺ DFT ͱ inverse DFT Λߴ଎ʹܭࢉ͢Δ͜ͱ ߴ଎ʹ DFT ΛٻΊΔΞϧΰϦζϜ ˠ ߴ଎ϑʔϦΤม׵ (Fast Fourier Transformation, FFT) tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 8 / 14
  9. ߴ଎ϑʔϦΤม׵ ଟ߲ࣜ f(x) = ∑ N−1 i=0 hixi (N ͸

    2 ͷ΂͖৐) ΛҎԼͷΑ͏ʹ 2 ͭʹ෼͚ ͯΈΑ͏ f0(x) = N 2 −1 ∑ i=0 h2ixi = h0x0 + h2x1 + h4x2 + . . . (10) f1(x) = N 2 −1 ∑ i=0 h2i+1xi = h1x0 + h3x1 + h5x2 + . . . (11) ݩͷଟ߲ࣜ͸ɺ f(x) = f0(x2) + xf1(x2) ͱදͤΔ f0, f1 ͸ͦΕͧΕ N 2 − 1 ࣍ҎԼͷଟ߲ࣜ tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 9 / 14
  10. ߴ଎ϑʔϦΤม׵ ઌ΄Ͳ ˆ f(t) ΛٻΊΔʹ͸ f ( ζ0 N )

    , . . . , f ( ζN−1 N ) ͷ஋͕ඞཁͩͬͨ ಉ༷ʹߟ͑Δͱɺ ˆ f0(t), ˆ f1(t) ΛٻΊΔʹ͸ɺ f(x) = f0(x2) + xf1(x2) ΑΓɺҎԼ͕ඞཁ f0 ( ζ0 N ) , f0 ( ζ2 N ) , . . . , f0 ( ζ2(N−1) N ) (12) f1 ( ζ0 N ) , f1 ( ζ2 N ) , . . . , f1 ( ζ2(N−1) N ) (13) tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 10 / 14
  11. ߴ଎ϑʔϦΤม׵ ζ2 N = exp ( 2 × 2π √

    −1 N ) = exp ( 2π √ −1 N 2 ) = ζN 2 ΑΓɺࣜ (12), (13) ͸ҎԼ ͷΑ͏ʹมܗՄೳ f0 ( ζ0 N 2 ) , f0 ( ζ1 N 2 ) , . . . , f0 ( ζN−1 N 2 ) (14) f1 ( ζ0 N 2 ) , f1 ( ζ1 N 2 ) , . . . , f1 ( ζN−1 N 2 ) (15) tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 11 / 14
  12. ߴ଎ϑʔϦΤม׵ ͜͜Ͱɺi = j mod N 2 ͳΒ͹ ζi N

    2 = ζj N 2 ࣜ (14), (15) ͸ͦΕͧΕલ൒ͱޙ൒͕ಉҰࢹͰ͖Δ Αͬͯલ൒෦෼ͷΈΛٻΊΕ͹Α͘ɺ۩ମతʹ͸ҎԼΛٻΊΕ͹ྑ͍ f0 ( ζ0 N 2 ) , f0 ( ζ1 N 2 ) , . . . , f0 ( ζ N 2 −1 N 2 ) (16) f1 ( ζ0 N 2 ) , f1 ( ζ1 N 2 ) , . . . , f1 ( ζ N 2 −1 N 2 ) (17) tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 12 / 14
  13. ߴ଎ϑʔϦΤม׵ ߴ଎ϑʔϦΤม׵ ໨ඪ: ˆ f(t) ΛٻΊΔͨΊʹඞཁͳ f ( ζ0 N

    ) , . . . , f ( ζN−1 N ) ͷ஋ΛٻΊΔ 1 f(x) = f0(x2) + xf1(x2) ͷΑ͏ʹɺf Λ 2 ͭͷؔ਺ f0, f1 ʹ෼͚Δ 2 2 ͭͷଟ߲ࣜ f0, f1 ʹରͯ͠ f0 ( ζ0 N 2 ) , f0 ( ζ1 N 2 ) , . . . , f0 ( ζ N 2 −1 N 2 ) (18) f1 ( ζ0 N 2 ) , f1 ( ζ1 N 2 ) , . . . , f1 ( ζ N 2 −1 N 2 ) (19) Λܭࢉ͢Δ͜ͱͰɺf ( ζ0 N ) , . . . , f ( ζN−1 N ) ΛٻΊΔ 3 ͜ΕΛ࠶ؼతʹ࣮ߦ (αΠζ͕൒෼ʹͳͬͨ΋ͷΛ 2 ͭղ͘͜ͱͷ܁ Γฦ͠) tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 13 / 14
  14. ߴ଎ϑʔϦΤม׵ લϖʔδͰઆ໌ͨ͠௨Γʹ࠶ؼతʹॲཧ͢Δͱɺ͜Ε͸ܭࢉྔ O(N log N) Ͱ͋Δʂ inverse DFT ʹؔͯ͠΋ɺDFT Ͱ

    ζi N Ͱ͋ͬͨͱ͜Ζ͕ ζ−i N ΁ͱม Θ͍ͬͯͯ 1 N ഒ͞Ε͍ͯΔ͚ͩͳͷͰɺߴ଎ϑʔϦΤม׵ͱಉ༷ʹߴ ଎ʹॲཧՄೳ ΑͬͯɺશମͰ O(N log N) Ͱղ͚ͨʂʂ ࣮૷ྫ Link - END - tsutaj (Hokkaido Univ.) ڝٕϓϩάϥϛϯάͷͨΊͷ FFT June 24, 2018 14 / 14