title: Getallen vermenigvuldigen in O(n log n) stappen.. Nieuw Archief voor Wiskunde, 2020. &#x27E8;hal-03013445&#x27E9;
text: Getallen vermenigvuldigen in O(n log n) stappen.
  Joris van der Hoeven
  
  To cite this version:
  Joris van der Hoeven. Getallen vermenigvuldigen in O(n log n) stappen..
  Wiskunde, 2020. �hal-03013445�
  
  Nieuw Archief voor
  
  HAL Id: hal-03013445
  https://hal.science/hal-03013445
  Submitted on 23 Nov 2020
  
  HAL is a multi-disciplinary open access
  archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from
  teaching and research institutions in France or
  abroad, or from public or private research centers.
  L’archive ouverte pluridisciplinaire HAL, est
  destinée au dépôt et à la diffusion de documents
  scientifiques de niveau recherche, publiés ou non,
  émanant des établissements d’enseignement et de
  recherche français ou étrangers, des laboratoires
  publics ou privés.
  Joris van der Hoeven
  
  Getallen vermenigvuldigen in O(n log n) stappen
  
  NAW 5/21 nr. 1 maart 2020
  
  55
  
  Joris van der Hoeven
  CNRS, Department of Mathematics
  Simon Fraser University, Burnaby, B.C.
  vdhoeven@lix.polytechnique.fr
  
  De oplossing
  
  Getallen vermenigvuldigen
  in O(n log n) stappen
  Onlangs werd door David Harvey en Joris van der Hoeven bewezen dat twee getallen van
  n cijfers kunnen worden vermenigvuldigd in O(n log n) stappen [12]. In 1971 werd het
  bestaan van zo’n methode door Schönhage en Strassen geuit als vermoeden [21]. Ook
  verwachtten zij dat het sneller niet gaat, iets dat we nog steeds niet weten. Joris van der
  Hoeven vat in dit artikel eerst een lange geschiedenis kort samen en geeft dan inzicht in
  het nieuwe bewijs.
   chine waarop we onze berekening doen,
  zijn we enkel geïnteresseerd in rekentijd
  op een constante factor na. Een honderd
  maal snellere methode is dus geen we­
  zenlijke verbetering, maar een log log log n
  
  Foto: Nina Ramaswamy
  
  Wat is de beste manier om twee gehele
  getallen met elkaar te vermenigvuldigen?
  Deze ogenschijnlijk eenvoudige vraag is
  misschien wel het oudste onopgeloste wis­
  kundeprobleem. Zo is het zeker tien maal
  ouder dan het vermoeden van Fermat, dat
  inmiddels de stelling van Wiles heet.
  Wél moet hier natuurlijk duidelijk wor­
  den gemaakt wat we met ‘beste’ bedoe­
  len. Dit werd eigenlijk pas mogelijk na de
  uitvinding van de Turingmachine: uitgerust
  met een precies rekenmodel kunnen we
  spreken over het aantal stappen M(n) dat
  nodig is om twee getallen van n cijfers met
  elkaar te vermenigvuldigen. De beste me­
  thode is dan diegene waarvoor M(n) zo
  langzaam mogelijk stijgt als n groter wordt.
  Neem bijvoorbeeld de methode van de
  lagere school. Dan vermenigvuldigen we
  elk cijfer van het eerste getal met elk cijfer
  van het tweede getal terwijl we de resul­
  taten op de juiste posities bij elkaar op­
  tellen. Dat geeft M(n) = O(n2): er bestaat
  een constante C > 0 met M(n) # Cn2 voor
  iedere n.
  Om onze hoofdvraag onafhankelijk te
  maken van de daadwerkelijke Turingma­
  
  maal sneller algoritme is een enorme stap
  voorwaarts.
  Waarschijnlijk formuleerde Kolmogorov
  de vraag voor het eerst op deze manier, tij­
  dens een beroemd seminarium in Moskou.
  In zijn enthousiasme lanceerde hij ook
  meteen het vermoeden dat n2 = O ^ M(n)h.
  Want, zo redeneerde hij, als er een betere
  methode bestond dan die van school, dan
  zou die sinds zesduizend jaar toch wel zijn
  uitgevonden.
   De nauwkeurige formulering van een
  probleem maakt het zoeken naar oplos­
  singen makkelijker. Maar het ontbreken
  daarvan sluit het zoeken naar oplossin­
  gen nou ook weer niet uit. De Babyloniërs
  berekenden 2 al tot op zestien decima­
  len achter de komma. Wellicht vinden we
  ooit een kleitablet waarop hogesnelheids­
  vermenigvuldigingsmethodes staan uitge­
  legd. Liefhebbers van het getal r, zoals de
  bij ons bekende Ludolph van Ceulen, zou­
  den daar vast wel wat aan kunnen hebben
  gehad (zie Figuur 1).
  Figuur 1 Gedenksteen in de Pieterskerk in Leiden, met
  een reconstructie van de oorspronkelijke tekst van de grafsteen van Ludolph van Ceulen (1540–1610), die vijfentwintig jaar van zijn leven besteedde aan het uitrekenen
  van 35 decimalen van r.
  Karatsuba
  Kolmogorovs vermoeden hield niet lang
  stand. Na drie weken werd hij verrast door
  een jonge student met een vermenigvuldi­
  gingsmethode in O(nlog 3/ log 2) stappen. Me­
  teen schreef hij het nieuwe algoritme op en
  publiceerde het onder de naam van de be­
  denker, Karatsuba, samen met een ander ar­
  tikel dat er niets mee te maken had [15, 16].
  56
  
  Getallen vermenigvuldigen in O(n log n) stappen
  
  NAW 5/21 nr. 1 maart 2020
  
  Joris van der Hoeven
  
  Nu kunnen we links altijd extra nullen aan
  een getal toevoegen zodat het aantal cij­
  fers een tweemacht wordt. Uiteindelijk be­
  wijzen we zo dat
  M(n) = O(nlog 3/ log 2).
  Andrej Kolmogorov
  
  Anatoly Karatsuba
  
  Laten we Karatsuba’s idee uitleggen met
  een voorbeeld:
  
  In het algemeen kost het vermenigvul­
  digen van twee maal zo lange getallen on­
  geveer drie maal zoveel tijd. Dit kan goed
  worden samengevat in een recurrente on­
  gelijkheid:
  
  u = 220629012020,
  v = 314159265358.
  Eerst splitsen we beide getallen in tweeën:
  a
  
  x
  
  b
  
  M(2n) # 3M(n) + O(n).
   644474448 64444744448 644474448
  u = 220629 # 1000000 + 012020 ,
  v = 314159 # 1000000 + 265358.
  144424443
  144424443
  
  Voor een zekere constante C en n = 2 r hebben we dus
  
  Daarna doen we twee optellingen en drie
  vermenigvuldigingen:
  
  # 9M(n/4) + _1 + 23 i Cn
  
  c
  
  d
  
  M(n) # 3M(n/2) + Cn
  
  # 27M(n/8) + _1 + 23 + 94 i Cn
  
  a + b = 232649,
  c + d = 579517,
  ac = 069312586011,
  bd = 003189603160,
  (a + b) (c + d) = 134824050533.
  h
  
  # 3 r M(1) + O __ 23 i n i
  r
  
  = O(3 r).
  Van getallen naar veeltermen
  Een van de ingrediënten van Karatsuba’s
  methode is het opsplitsen van getallen in
  twee delen. Dit maakte het mogelijk om u
  en v te zien als veeltermen ax + b en cx + d
  van graad < 2. In feite kwam dat simpel­
  weg neer op het werken in een talstel­
  sel met grondtal x = 1000 000 in plaats
  van 10.
  Dit idee om het rekenen met gehele ge­
  tallen te vervangen door het rekenen met
  veeltermen werd in the negentiende eeuw
  reeds geopperd door Kronecker. Zo kunnen
  we een getal van n cijfers opsplitsen in l
  brokken van p|= ^n/l h cijfers en daar weer
  een veelterm van maken.
  Door het aantal brokken l groter dan
  twee te kiezen, werden er in de zestiger
  jaren een reeks snellere varianten van Ka­
  ratsuba’s methode ontwikkeld [17, 20, 22];
  zie Tabel 1 voor een historisch overzicht.
  Zonder al te ver op details in te gaan, is
  elk van deze varianten gebaseerd op een
  veralgemening van Karatsuba’s rekentruc.
   In dit geval wordt de vermenigvuldiging
  van twee veeltermen van graad < l gere­
  duceerd tot 2l - 1 vermenigvuldigingen van
  coëfficiënten, waarna
  M(n) = O(nlog(2l - 1)/ log(l - 1)).
  Nu merken we op dat het slechts twee keer
  aftrekken vergt om
  
  -4000
  
  Onbekend
  
  O(n2)
  
  1962
  
  Karatsuba
  
  O(nlog 3/ log 2)
  
  ad + bc = (a + b) (c + d) - ac - bd
  = 62321861362
  te bepalen. Tot slot geldt
  
  1963
  
  Toom
  
  O(n25 log 2 )
  
  1966
  
  Schönhage
  
  uv = (ax + b) (cx + d)
  = acx2 + (ad + bc) x + bd,
  
  1969
  
  Knuth
  
  1971
  
  Pollard
  
  O(n2 log 2 log n)
  O(n log n log log ng)
  
  zodat we klaar zijn na een laatste optel­
  som
  
  1971
  
  Schönhage-Strassen
  
  O(n log n log log n)
  
  2007
  
  Fürer
  
  O(n log n 2O(log n))
  
  acx2
  069312586011000000000000
  +
  62321861362000000 (ad bc) x
  bd
  003189603160
  uv
  069312648332864551603160
  
  2014
  
  Harvey-vdH-Lecerf
  
  O(n log n 8log n)
  
  2017
  
  Harvey
  
  O(n log n 6log n)
  
  2017
  
  Harvey-vdH
  
  Al met al hebben we een vermenigvuldi­
  ging van twaalf cijfers teruggevoerd tot
  drie vermenigvuldigingen van zes cijfers en
  een stel optel- en aftreksommen.
  2018
  
  Harvey-vdH
  
  O(n log n ^4 2 h
  
  2019
  
  Harvey-vdH
  
  logn
  
  O(n2
  
  2 logn
  log 2
  
  3
  
  (log n) 2 )
  
  2 logn
  
  )
  
  )
  )
  
  log) n
  
  )
  
  )
  
  O(n log n 4log n)
  O(n log n)
  
  Tabel 1 Betere en betere bovengrenzen voor M(n) in de loop van de jaren... In de tabel is de functie log# gedefinieerd
  n#
  door log) x = min{n ! N : (log % f
  % log) (x) G 1}.
   Getallen vermenigvuldigen in O(n log n) stappen
  
  Joris van der Hoeven
  
  Van veeltermen naar kringtermen
  Vanaf nu gaan we het vooral hebben over
  het vermenigvuldigen van veeltermen
  en is het handig om te werken met coëf­
  ficiënten in een algemene ring R. De pre­
  cieze keus van R laten we voorlopig in het
  midden.
  Voor wat komen gaat, is het ook han­
  dig om met zogenaamde ‘kringtermen’
  in plaats van veeltermen te werken. Een
  kringterm van lengte l is een element van
  R[x] /(x l - 1).
  De relatie x l = 1 treedt pas in actie zodra
  de graad van een veelterm over de l gaat.
  Het uitrekenen van een product van twee
  veeltermen in R[x] van graad < l/2 kunnen
  we dus net zo goed in R[x] /(x l - 1) doen.
  Het voordeel van kringtermen is dat we
  een aantal nieuwe rekentrucs rijker wor­
  den. Veronderstel bijvoorbeeld dat l = 2
  en dat 2 omkeerbaar is in R. Dan is het
  mogelijk om Karatsuba’s methode verder
  te verbeteren: modulo x2 - 1 geldt
  (a1 x + a0) (b1 x + b0)
  
  at bt - at bt
  at bt + at bt
  = 0 0 2 1 1 x + 0 0 2 1 1,
  
  waar
  at0 = a0 + a1, bt0 = b0 + b1,
  at1 = a0 - a1, bt1 = b0 - b1.
  orde l en leidt tot de ontbinding
  xl - 1 =
  
  % (x - ~ k).
  0Gk<l
  
  Is l ook nog eens omkeerbaar in R, dan
  kan bovendien de Chinese reststelling wor­
  den toegepast op deze ontbinding:
  R[x] /(x l - 1) ,
  
  % R[x] /(x - ~ k).
  0Gk<l
  
  Van links naar rechts heet dit isomorfisme
  de discrete fouriertransformatie of DFT.
   Voor alle A ! R[x] en 0 # k < l hebben we
  x = ~ k en P(x) = P(~ k) modulo x - ~ k, zo­
  dat
  DFT(Ar ) = (A(1), A(~), f, A(~ l - 1)).
  De algebra R[x] /(x - ~ k) is op zijn beurt
  weer isomorf met R voor elke k, zodat
  R[x] /(x l - 1) , R l.
  Nu is rekenen aan de rechterkant goed­
  koop: één vermenigvuldiging in R l komt
  bijvoorbeeld neer op l vermenigvuldigin­
  gen in R. Als we zowel de DFT en de in­
  verse DFT-1 snel uit kunnen rekenen, dan
  hebben we dus een goede methode om
  kringtermen A, B ! R[x] /(x l - 1) te verme­
  nigvuldigen:
  AB = DFT-1 (DFT(A) DFT(B)).
  Dit heet FFT-vermenigvuldigen.
  2
  
  In R[x] /(x - 1) vergt het uitrekenen van
  een product dus slechts twee vermenig­
  vuldigingen in R. Ook zijn er een aantal
  optellingen, aftrekkingen en delingen door
  twee nodig.
  Bovenstaande relaties zijn het gevolg
  van
  (x2 - 1) = (x - 1) (x + 1)
  
  (1)
  
  en toepassing van de Chinese reststelling
  op deze ontbinding:
  R[x] /(x2 - 1) , R[x] /(x - 1) # R[x] /(x + 1),
  a1 x + a0 7 ( a0 + a1, a0 - a1).
  Dit maakt het mogelijk willekeurige bereke­
  ningen in R[x] /(x2 - 1) te vervangen door
  berekeningen in R[x] /(x - 1) # R[x] /(x + 1)
  , R 2.
  FFT-vermenigvuldigen
  De ontbinding (1) geldt voor elke ring R. In
  sommige ringen kan ook x l - 1 op soortge­
  lijke wijze ontbonden worden voor l > 2.
   Dat is het geval zodra R een element ~
  l-1
  bezit met / k = 0 (~ j) k = 0. Zo’n element ~
  heet een principiële eenheidswortel van
  
  De snelle fouriertransformatie
  Maar hoe rekenen we zo’n DFT snel uit?
  Het geval l = 2 bestudeerden we al eerder.
  Meer in het algemeen leidt de ontbinding
  x2l - 1 = (x l - 1) (x l + 1)
  voor even lengtes 2l tot een isomorfisme
  R[x] /(x2l - 1) , R[x] /(x l - 1)# R[x] /(x l + 1).
  Als A, B ! R[x] veeltermen van graad < l
  zijn, dan zendt dit isomorfisme A + Bx l
  naar ^ A + B, A - B h, waar we voor het ge­
  mak de modulostrepen hebben weggela­
  ten. Het uitrekenen van A ! B komt neer
  op l maal optellen en aftrekken in R. Voor
  het isomorfisme in de andere richting ko­
  men daar 2l delingen door twee bij.
  Als ~ een principiële eenheidswortel
  van orde 2l is, dan hebben we bovendien
  het volgende isomorfisme:
  R[x] /(x l + 1) , R[x] /(y l - 1),
  / ak x k 7 / (ak ~ k) y k.
  0Gk<l
  
  0Gk<l
  
  Het uitrekenen van dit isomorfisme komt
  neer op l vermenigvuldigingen met mach­
  
  NAW 5/21 nr. 1 maart 2020
  
  57
  
  ten van ~. Hier is het belangrijk om op
  te merken dat een vermenigvuldiging met
  een macht van ~ soms goedkoper is dan
  een willekeurige vermenigvuldiging in R.
   Al met al laat dit zien hoe een DFT
  van lengte 2l kan worden gereduceerd tot
  twee DFTs van lengte l. Schrijven we F(l)
  voor de tijd die het kost om een DFT van
  lengte l uit te rekenen, T! voor de tijd van
  een optelling of aftrekking in R, en T~ voor
  de tijd van een vermenigvuldiging met een
  macht van ~, dan hebben we
  F(2l) # 2F(l) + 2lT! + lT~.
  Passen we deze formule recursief toe in
  het geval dat l = 2lg l een tweemacht is,
  dan vinden we
  F(l) # 2F(l/2) + l _T! + 12 T~ i
  
  # 4F(l/4) + 2l _T! + 12 T~ i
  h
  
  # (l/2) F(2) + (lg l - 1) l _T! + 12 T~ i
  # l lg l _T! + 12 T~ i.
  (2)
  
  Voor de inverse transformatie komen daar
  ook nog l delingen door l bij.
  Intermezzo
  Voordat we verder gaan, is het een mooi
  moment voor wat opmerkingen. De FFT
  brak door in 1965, na de publicatie van
  Cooley en Tukey’s artikel [3]. Maar een
  soortgelijke methode werd eerder beschre­
  ven in ongepubliceerd werk van Gauss
  [6, 14].
  Daarnaast hebben we gezien dat een
  vermenigvuldiging van kringtermen kan
  worden gereduceerd tot twee directe plus
  één inverse DFT. In 1970 merkte Bluestein
  op dat een DFT van lengte l ook kan wor­
  den gereduceerd tot een kringproduct van
  dezelfde lengte l [2].
  Laten we, om dit uit te leggen, voor het
  gemak veronderstellen dat l even is en h
  een principiële eenheidswortel van orde 2l
  met h2 = ~. Voor gehele j, k geldt dan
  2
  
  2
  
  h(j + l) = h j h2l(j + l/2) = h j
  
  2
  
  en
  2
  
  2
  
  2
  
  ~ jk = h j h k h-(j - k).
   De j-de coëfficiënt van de DFT van een
  kringterm u0 + g + u l - 1 x l - 1 is dus gelijk
  aan
  
  / uk ~ jk = h j2 / `uk h k2 jh-(j - k)2.
  0Gk<l
  
  0Gk<l
  
  De rechtersom herkennen we nu als de
  
  58
  
  NAW 5/21 nr. 1 maart 2020
  
  Getallen vermenigvuldigen in O(n log n) stappen
  
  Joris van der Hoeven
  
  van ~: we schuiven gewoon de bits van
  het getal op terwijl we de relatie 2 m = - 1
  gebruiken zodra we de 2m passeren. Voor
  l ! {m, 2m} laten Schönhage en Strassen
  vervolgens zien hoe een vermenigvuldiging
  in Z/(2 lm/2 + 1) Z kan worden gereduceerd
  tot l vermenigvuldigingen in R. Schrijven
  we M'(m) voor de kosten van één verme­
  nigvuldiging in R, dan geeft dit
  M'(lm/2) # lM'(m) + O(ml log l).
  (3)
  
  De term O(ml log l) komt van de DFTs, waar
  we gebruik maken van het feit dat T~ = O(l)
  in R. De term lM'(m) komt van de ‘interne’
  l-1
  vermenigvuldiging in % k = 0 R[x] /(x - ~ k)
  l
  ,R.
  Omgerekend ten opzichte van n leidt (3)
  grofweg tot de ongelijkheid
  
  Foto: MFO, Konrad Jakobs
  
  M'(n) # 2n1/2 M'(n1/2) + Cn log n, (4)
  
  Volker Strassen (links) en Arnold Schönhage (rechts)
  
  j-de coëfficiënt van het product van de
  volgende twee kringtermen:
  V =
  
  / `uk h
  
  0Gk<l
  
  W=
  
  k2 j x k,
  
  / h-k2 x k.
  0Gk<l
  
  Terug naar de gehele getallen
  In mijn geboortejaar 1971 verschenen er
  maar liefst drie methodes om de FFT toe
  te passen op het vermenigvuldigen van ge­
  hele getallen [19, 21]. Elk van deze metho­
  des was gebaseerd op een aparte keuze
  van R.
   Het meest voor de hand ligt het om
  R = C en ~ = e2ri/l te nemen. In ons ge­
  val kunnen we natuurlijk alleen maar met
  benaderingen van complexe getallen wer­
  ken. Voor het vermenigvuldigen van getal­
  len van n cijfers kan worden bewezen dat
  benaderingen tot op C log n cijfers achter
  de komma voldoende zijn voor een zekere
  constante C > 0. Dat betekent dat we kun­
  nen werken met l = O(n/ log n) brokken
  van O(log n) cijfers. In combinatie met (2)
  vinden we dan
  M(n) = O ^l log lM(log n)h = O ^nM(log n)h.
  Laten we deze formule recursief op zichzelf
  los, dan vinden we
  
  M(n) = O(n log n log log n log log log ng).
  Dit was het ‘eerste algoritme’ in het artikel
  van Schönhage–Strassen [21].
  Het nulde algoritme werd echter al iets
  eerder door Pollard ontdekt [19]. Híj stel­
  de voor om R = Fp te nemen, waar p een
  priemgetal is van de vorm p = s2 r + 1; hoe
  kleiner s, hoe beter. Voor zulke priemgetal­
  len p weten we dat er primitieve eenheids­
  wortels van orde 2r bestaan. Opnieuw is
  het dan mogelijk om log p = O(log n) en
  l = O(n/ log n) te kiezen, en we vinden
  een soortgelijke bovengrens voor M(n) als
  net. Eerlijkheidshalve moet er wel worden
  bijgezegd dat dit laatste níet in Pollards
  artikel stond. Hij was meer uit op een prak­
  tisch algoritme en beschreef verschillende
  optimalisaties. Voor zeer grote getallen is
  zijn algoritme op hedendaagse computers
  het beste [8].
   Voor het ‘tweede algoritme’ van
  Schönhage–Strassen schakelen we over
  naar het binaire stelsel en nemen we
  R = Z/^2 m + 1 h Z, waar m een tweemacht
  is. In deze ring hebben we per definitie
  2 m = - 1, zodat ~|= 2 een principiële een­
  heidswortel van orde 2m is. Een ander
  voordeel is dat ~ een snelle eenheidswortel is. Daarmee bedoelen we dat we snel
  kunnen vermenigvuldigen met machten
  
  voor een zekere constante C. Deze formule
  kunnen we vervolgens ‘uitrollen’:
  M'(n) # 2n1/2 M'(n1/2) + Cn log n
  # 4n3/4 M'(n1/4) + 2Cn log n
  # 8n7/8 M'(n1/8) + 3Cn log n
  h
  # n log n M'(O(1)) + Cn log n log log n.
  Dit bewijst dat M(n) = O(n log n log log n),
  een bovengrens die vijfenveertig jaar stand
  hield.
  Hoe verder?
  Van de drie op de FFT gebaseerde me­
  thodes uit 1971 heeft de laatste methode
  één groot nadeel: omdat l # 2m wordt een
  vermenigvuldiging van lengte n ‘slechts’
  gereduceerd tot vermenigvuldigingen van
  lengte O^ n h in plaats van O(log n). In te­
  genstelling tot de twee andere methodes
  kosten de DFTs daarentegen bijna niets.
  John Pollard
  
  Getallen vermenigvuldigen in O(n log n) stappen
  
  Joris van der Hoeven
  
  Dit levert twee sporen op voor verdere ver­
  beteringen.
  Een eerste mogelijkheid is om de kos­
  ten van DFTs over C en / of Fp te drukken.
  In 2007 lukte Fürer dit voor het eerst [5].
  Zijn methode leidt tot een ongelijkheid van
  de vorm
  M(n)
  M(n')
  n log n # K n' log n' + O(1),
  n' = O((log n) 2),
  en de verbeterde bovengrens
  )
  
  M(n) = O(n log nKlog n),
  k # % log h (x) # 1 ,.
   log) x = min " k ! N: ^log % f
  Fürer besteedde geen aandacht aan de
  precieze ‘uitdijingsfactor’ K > 1. Vanaf 2014
  lukte het David Harvey, Grégoire Lecerf en
  mijzelf om deze factor steeds verder naar
  beneden te jagen [9, 10, 11, 13]; zie Tabel 1.
  Onze tweede optie is een verscherping
  van de ongelijkheid (4). Hier merken we
  op dat de factor 2 in 2n1/2 M'(n1/2) pre­
  cies wordt gecompenseerd door het feit
  dat log n = 12 log n: bij het uitrollen kost
  elke iteratie precies Cn log n stappen extra.
  Kunnen we deze factor 2 maar een fractie
  verbeteren, dan verloopt het uitrollen zich
  als volgt:
  M(n) # 1,98n1/2 M(n1/2) + Cn log n
  # 1,982 n3/4 M(n1/4)
  + (1 + 0,99) Cn log n
  # 1,983 n7/8 M(n1/8)
  + (1 + 0,99 + 0,992) Cn log n
  h
  # o(n log n) + 100Cn log n.
  Ho! Lezen we dat goed?
  M(n) = O(n log n).
  Op deze wens kunnen we voortborduren.
  Zo is het ook genoeg om te bewijzen dat
  M(n) # an(d - 1)/d M(n1/d) + O(n log n), (5)
  waar d $ 2 en a < 1/d.
  Op naar hogere dimensies
  Nu heeft de methode van Schönhage en
  Strassen inderdaad een frustrerend aspect:
  de kunstmatig aangemaakte eenheidswor­
  tel ~ = 2 in R wordt eigenlijk reuze ‘on­
  dergebruikt’. Om die reden is onze lengte­
  reductie n U n zeer bescheiden.
  Nu zijn er ook DFTs die in meer dan
  één richting werken. Stel opnieuw dat R
  een willekeurige ring is met een principië­
  le eenheidswortel ~ van orde l. Een d-di­
  
  mensionale DFT bewerkstelligt dan het
  isomorfisme
  R[x1, f, xd] /(x1l - 1, f, x dl - 1)
  ,
  
  % R[x1, f, xd] /(x1 - ~1k1, f, xd - ~dkd).
   0 # k1, f, kd < l
  
  Als we de DFTs in x2, f, xd nu vervangen
  door DFTs over de ring R[x1] /(x1l - 1), dan
  kosten ze veel minder, want x1 is een snel­
  le eenheidswortel in deze ring. Het syste­
  matische gebruik van dit soort DFTs werd
  voor het eerst voorgesteld door Nussbau­
  mer en Quandalle [18].
  Met dat idee maken we meteen al
  een grote stap in de richting van (5).
  Als n = l d, dan vergt een vermenigvuldiging
  in R[x1, f, xd] /(x1l - 1, f, x dl - 1) namelijk
  O(n log n) optellingen en aftrekkingen in
  R plus n(d - 1)/d vermenigvuldigingen in
  R[x1] /(x1l - 1).
  Er is echter één probleem: snel vermenig­
  vuldigen in R[x1, f, xd] /(x1l - 1, f, x dl - 1) is
  niet hetzelfde als snel vermenigvuldigen in
  d
  R[x] /(x l - 1). We hebben dus een manier
  nodig om kringtermen in één variabele x
  om te smeden tot kringtermen in d varia­
  belen x1, f, xd.
  Met hulp van de oude Chinezen
  Nu is het in bepaalde gevallen inderdaad
  mogelijk om van dimensie te veranderen.
  Stel dat l1, f, ld onderling ondeelbaar zijn.
  Volgens de Chinese reststelling hebben we
  dan
  Z/(l1 gld) Z , Z/l1 Z + g + Z/ld Z.
  Formeel gesproken hebben we dus ook
  xZ/(l1 gld) Z , x1Z/l1 Z # g # x dZ/ld Z.
  Nemen we vervolgens lineaire combinaties,
  dan zien we uiteindelijk dat
  R[x] /(x l1 gld - 1)
  , R[x1, f, xd] /(x1l1 - 1, f, x dld - 1).
  In verband met DFTs werd dit voor het
  eerst opgemerkt door Good [7]. Agarwal en
  Cooley gebruikten dit isomorfisme daarna
  voor het uitrekenen van convoluties [1].
   Wel zitten we nu weer met een nieuw
  probleem: in de vorige paragraaf vroegen
  we om een isomorfisme waarvoor alle li
  gelijk zijn. Maar ons isomorfisme werkt
  juist alleen in het geval wanneer l1, f, ld
  onderling ondeelbaar zijn.
  Wat er tot slot nog ontbreekt is een ma­
  nier om de lengte van een DFT iets te
  veranderen. Dit zou het mogelijk maken
  
  59
  
  NAW 5/21 nr. 1 maart 2020
  
  om een d-dimensionale DFT van leng­
  te (l1, f, ld) om te smeden tot eentje van
  lengte (l, f, l). Met behulp van een d-di­
  mensionale versie van Bluesteins algo­
  ritme kan zo’n DFT vervolgens worden
  gereduceerd tot een vermenigvuldiging in
  R[x1, f, xd] /(x1l - 1, f, x dl - 1). En dat kan
  uiteindelijk weer efficiënt worden gedaan
  met behulp van snelle eenheidswortels.
  Gaussiaans overmonsteren
  Hoe kunnen we een DFT van lengte s ver­
  vangen door een DFT van een iets grotere
  lengte t? Om dit te bereiken veronderstel­
  len we vanaf hier dat R = C.
  Neem een DFT van lengte s. In de
  signaaltheorie wordt de invoer gezien als
  een reeks monsters van een signaal. De
  frequentie van de monstername is hierbij
  evenredig aan s. Is het mogelijk het oor­
  spronkelijke signaal uit onze reeks mon­
  sters te reconstrueren? Dan zouden we
  daarna een nieuwe reeks monsters kunnen
  nemen, met een andere frequentie.
  De meest voor de hand liggende manier
  om een digitaal signaal analoog te ma­
  ken is via convolutie met een Gaussiaan
  2
  Ga (x) = e- ax . Hoe kleiner a, hoe gladder
  (en minder scherp) het analoge signaal.
  Een andere prettige eigenschap is dat de
  fouriertransformatie van Ga opnieuw een
  Gaussiaan is.
   Laten we deze ideeën in formules verta­
  len. In plaats van kringtermen van lengte s
  beschouwen we nu de bijbehorende vecto­
  ren van coëfficiënten u ! C s. Wel spreken
  we af dat uk + js / uk voor alle j ! Z. Gege­
  ven u ! C s definiëren we Fs (u) ! C s door
  (Fs u) k |=
  
  jk
  
  / uk e-2ri s.
  0Gk<s
  
  Dit is een variant op onze eerdere definitie
  van een DFT (nu nemen we ~ = e-2ri/s).
  Vervolgens definiëren we twee lineaire af­
  beeldingen S, T : C s " Ct door
  (Su)k |= a-1 / e- ra
  
  -2 2 k - j 2
  )
  s (
  
  t
  
  s
  
  j!Z
  
  (Tu)k |= / e- ra
  
  j
  -2 2 k
  t ( - )2
  t
  
  j!Z
  
  s
  
  u j,
  
  uj.
  Tot slot definiëren we twee permutaties
  Ps : C s " C s en Pt : Ct " Ct door
  (Ps u) j|= utj,
  (Pt u) k|= u-sk.
   In [12, Theorem 4.2] bewijzen we dan dat
  het volgende diagram commuteert:
  
  60
  
  Getallen vermenigvuldigen in O(n log n) stappen
  
  NAW 5/21 nr. 1 maart 2020
  
  −→
  
  Kronecker
  
  −→
  
  projectie zonder exces
  
  −→
  
  Chinese reststelling
  
  −→
  
  FFT-vermenigvuldigen
  
  −→
  
  Gaussiaans overmonsteren
  
  −→
  
  Bluestein reductie
  
  −→
  
  Z
  
  ×
  
  Nussbaumer FFT
  
  C[x]
  
  ×
  ×
  
  C[x]/(x − 1)
  
  ×
  
  C[x1 , . . . , xd ]/(xl11 − 1, . . . , xldd − 1)
  
  DFT
  
  C[x1 , . . . , xd ]/(xl11 − 1, . . . , xldd − 1)
  
  DFT
  
  C[x1 , . . . , xd ]/(xl1 − 1, . . . , xld − 1)
  
  ×
  
  C[x1 , . . . , xd ]/(xl1 − 1, . . . , xld − 1)
  
  ×
  
  C[x]/(xl − 1)
  
  L
  
  Figuur 2 Schematische weergave van de verschillende reducties in het nieuwe algoritme. Hier en daar hebben we ietwat
  gesmokkeld om alles simpel te houden.
  Cs
  S.
  Fs
  
  Cs
  
  Ps
  
  Cs
  .T
  
  Ct
  
  Ft
  
  Ct
  
  Pt
  
  Ct
  
  Dit is wat we gebruiken om het berekenen
  van Fs terug te voeren op het berekenen
  van Ft . Omdat t > s, zijn de matrices voor
  S en T niet vierkant. Wél is het zo dat de
  elementen die niet op de diagonaal liggen
  razendsnel afnemen. Dat komt doordat
  Gaussianen uit het centrum ook razendsnel
  afnemen. Door t - s goed gekozen rijen uit
  T weg te nemen blijft er om die reden
  
  Referenties
  1
  
  2
  
  3
  
  4
  
  5
  
  6
  
  7
  
  R. Agarwal en J. Cooley, New algorithms for
  digital convolution, IEEE Transactions on
  Acoustics, Speech, and Signal Processing
  25(5) (1977), 392–410.
  Leo I. Bluestein, A linear filtering approach
  to the computation of discrete Fourier trans­
  form, IEEE Transactions on Audio and Electroacoustics 18(4) (1970), 451–455.
   J. W. Cooley en J. W. Tukey, An algorithm for
  the machine calculation of complex Fourier
  series, Math. Computat. 19 (1965), 297–301.
  A. Dutt en V. Rokhlin, Fast Fourier trans­
  forms for nonequispaced data, SIAM J. Sci.
  Comput. 14(6) (1993), 1368–1393.
  M. Fürer, Faster integer multiplication, in
  Proceedings of the Thirty-Ninth ACM Symposium on Theory of Computing (STOC
  2007), San Diego, CA, 2007, pp. 57–66.
  C. F. Gauss, Nachlass: theoria interpolationis
  methodo nova tractata, in Werke, Vol. 3,
  Königliche Gesellschaft der Wissenschaften,
  Göttingen, 1866, pp. 265–330.
  I. J. Good, The interaction algorithm and
  practical Fourier analysis, Journal of the
  Royal Statistical Society, Series B. 20(2)
  (1958), 361–372.
  een bijna diagonale matrix over. Dit kan
  worden benut om T-1 Pt Ft S snel uit te
  rekenen.
  Kiezen we a en de precisie met beleid,
  dan kunnen we bewijzen dat Fs op deze
  manier met bijna evenveel precisie als Ft
  kan worden berekend en dat de rekentijd
  van S and T-1 verwaarloosbaar is ten op­
  zichte van die van Ft.
  Daarmee is de kous af en hebben we
  bewezen dat M(n) = O(n log n). In Figuur 2
  worden alle benodigde reducties nog eens
  samengevat.
  8
  
  9
  10
  
  11
  
  12
  
  13
  
  14
  
  D. Harvey, Faster arithmetic for numberthe­
  oretic transforms, J. Symbolic Comput. 60
  (2014), 113–119.
  D. Harvey, Faster truncated integer multipli­
  cation, 2017, arXiv:1703.00640.
  D. Harvey en J. van der Hoeven, Faster inte­
  ger and polynomial multiplication using cy­
  clotomic coefficient rings, Technical Report,
  2017, arXiv:1712.03693.
  D. Harvey en J. van der Hoeven, Faster inte­
  ger multiplication using short lattice vectors,
  in R. Scheidler en J. Sorenson, eds., Proc. of
  the 13-th Algorithmic Number Theory Symposium, Open Book Series 2, Mathematical
  Sciences Publishes, 2019, pp. 293–310.
   D. Harvey en J. van der Hoeven, Integer multi­
  plication in time O(n log n), Technical Report,
  HAL, 2019, http://hal.archivesouvertes.fr/hal02070778.
  D. Harvey, J. van der Hoeven en G. Lecerf,
  Even faster integer multiplication, J. Complexity 36 (2016), 1–30.
  M. T. Heideman, D. H. Johnson en C. S. Bur­
  rus, Gauss and the history of the FFT, IEEE
  Acoustics, Speech and Signal Processing
  Magazine 1 (1984), 14–21.
  Joris van der Hoeven
  
  Een variant van onze overmonsterings­
  methode werd voor het eerst gepubliceerd
  door Dutt en Rokhlin [4]. Deze variant is
  algemener en maakt het mogelijk DFTs te
  berekenen van onregelmatig gemonsterde
  signalen. Daarentegen werkt hun metho­
  de alleen voor log s = O(a) en is daarom
  net iets te traag voor ons uiteindelijke
  doel.
  Wat hebben we eraan?
  Praktisch gezien moeten we dat afwachten.
  Voor theoretische doeleinden is de func­
  tie M(n) van belang om de kosten van
  allerlei rekenkundige operaties nauwkeu­
  rig te beschrijven. Zo kost het delen van
  getallen van n cijfers O ^ M(n)h = O(n log n)
  en het uitrekenen van een grootste ge­
  mene deler O ^ M(n) log nh = O(n log2 n).
  Ook n decimalen van r kunnen nu in
  O ^ M(n) log nh = O(n log2 n) stappen wor­
  den berekend. Een DFT van lengte l over
  de complexe getallen kost O ^ M(lp)h als
  we rekenen met een precisie van p $ log l
  cijfers [13]. Dat bepaalt ook de minimale
  CO2-uitstoot die men kan verwachten bij
  grote berekeningen aan het weer. Al met
  al speelt M(n) als elementaire ‘rekensnel­
  heid’ dus een zelfde soort rol als de licht­
  snelheid c in de natuurkunde.
  Kan het nog sneller?
  Dat weten we niet!
  
  s
  
  15 A. A. Karatsuba, The complexity of compu­
  tations, Proc. of the Steklov Inst. of Math.
   211 (1995), 169–183. English translation;
  Russian original at pp. 186–202.
  16 A. Karatsuba en J. Ofman, Multiplication
  of multidigit numbers on automata, Soviet
  Physics Doklady 7 (1963), 595–596.
  17 D. E. Knuth, The Art of Computer Programming, Volume 2: Seminumerical Algorithms,
  Addison-Wesley, 1969.
  18 H. J. Nussbaumer en P. Quandalle, Compu­
  tation of convolutions and discrete Fourier
  transforms by polynomial transforms, IBM J.
  Res. Develop. 22(2) (1978), 134–144.
  19 J. M. Pollard, The fast Fourier transform in
  a finite field, Mathematics of Computation,
  25(114) (1971), 365–374.
  20 A. Schönhage, Multiplikation großer Zahlen,
  Computing 1(3) (1966), 182–196.
  21 A. Schönhage en V. Strassen, Schnelle
  Multiplikation großer Zahlen, Computing 7
  (1971), 281–292.
  22 A. L. Toom, The complexity of a scheme of
  functional elements realizing the multipli­
  cation of integers, Soviet Mathematics 4(2)
  (1963), 714–716.
  