Editing Schönhage–Strassen algorithm (section)

==Implementation details==

===Why ''N'' = 2{{sup|''M''}} + 1 mod ''N'' ===

In Schönhage–Strassen algorithm, <math> N = 2^M + 1 </math>. This should be thought of as a binary tree, where one have values in <math> 0 \le  \text{index} \le 2^{M}=2^{i+j} </math>. By letting <math> K \in [0,M] </math>, for each {{mvar|K}} one can find all <math> i+j = K </math>, and group all <math>(i,j)</math> pairs into M different groups. Using <math> i+j = k </math> to group <math>(i,j)</math> pairs through convolution is a classical problem in algorithms.<ref>{{cite book |last1=Kleinberg |first1=Jon |last2=Tardos |first2=Eva |title=Algorithm Design |date=2005 |publisher=Pearson |isbn=0-321-29535-8 |page=237 |edition=1|url=https://archive.org/details/AlgorithmDesign1stEditionByJonKleinbergAndEvaTardos2005PDF/page/n259/mode/2up}}</ref>

Having this in mind, <math> N = 2^M + 1 </math> help us to group <math> (i,j) </math> into <math> \frac{M}{2^k} </math> groups for each group of subtasks in depth {{mvar|k}} in a tree with <math> N = 2^{\frac{M}{2^k}} + 1  </math>

Notice that <math> N = 2^M + 1 = 2^{2^L} + 1 </math>, for some L. This makes N a [[Fermat number]]. When doing mod <math> N = 2^M + 1 = 2^{2^L} + 1 </math>, we have a Fermat ring.

Because some Fermat numbers are Fermat primes, one can in some cases avoid calculations.

There are other ''N'' that could have been used, of course, with same prime number advantages. By letting <math> N = 2^k - 1 </math>, one have the maximal number in a binary number with <math> k+1 </math> bits.
<math> N = 2^k - 1 </math> is a Mersenne number, that in some cases is a Mersenne prime. It is a natural candidate against Fermat number <math> N = 2^{2^L} + 1 </math>

==== In search of another ''N'' ====

Doing several mod calculations against different {{mvar|N}}, can be helpful when it comes to solving integer product. By using the [[Chinese remainder theorem]], after splitting {{mvar|M}} into smaller different types of {{mvar|N}}, one can find the answer of multiplication {{mvar|xy}} <ref>{{cite web |last1=Gaudry |first1=Pierrick |last2=Alexander |first2=Kruppa |last3=Paul |first3=Zimmermann |title=A GMP-based implementation of Schönhage-Strassen's large integer multiplication algorithm |page=6 | year=2007|url=https://inria.hal.science/inria-00126462/file/fft.final.pdf}}</ref>

Fermat numbers and Mersenne numbers are just two types of numbers, in something called generalized Fermat Mersenne number (GSM); with formula:<ref>{{cite web |last1=S. Dimitrov |first1=Vassil |last2=V. Cooklev |first2=Todor |last3=D. Donevsky |page=2|year=1994|first3=Borislav |title=Generalized Fermat-Mersenne Number Theoretic Transform |url=https://www.researchgate.net/publication/3324542}}</ref>

:<math> G_{q,p,n} = \sum_{i=1}^p q^{(p-i)n} = \frac{q^{pn}-1}{q^n-1} </math>

:<math> M_{p,n} =  G_{2,p,n} </math>

In this formula, <math> M_{2,2^k} </math> is a Fermat number, and <math>  M_{p,1} </math> is a Mersenne number.

This formula can be used to generate sets of equations, that can be used in CRT (Chinese remainder theorem):<ref>{{cite web |last1=S. Dimitrov |first1=Vassil |last2=V. Cooklev |first2=Todor |last3=D. Donevsky |page=3|year=1994|first3=Borislav |title=Generalized Fermat-Mersenne Number Theoretic Transform |url=https://www.researchgate.net/publication/3324542}}</ref>

:<math>g^{\frac{(M_{p,n}-1)}{2}} \equiv -1  \pmod {M_{p,n}}</math>, where {{mvar|g}} is a number such that there exists an {{mvar|x}} where  <math> x^2 \equiv g \pmod {M_{p,n}}</math>, assuming <math> N = 2^n </math>

Furthermore; <math>g^{2^{(p-1)n}-1} \equiv a^{2^n -1}  \pmod {M_{p,n}}</math>, where {{mvar|a}} is an element that generates elements in <math> \{1,2,4,...2^{n-1},2^n\} </math> in a cyclic manner.

If <math> N=2^t </math>, where <math> 1 \le t \le n </math>, then <math> g_t = a^{(2^n-1)2^{n-t}} </math>.

=== How to choose ''K'' for a specific ''N'' ===

The following formula is helpful, finding a proper {{mvar|K}} (number of groups to divide {{mvar|N}} bits into) given bit size {{mvar|N}} by calculating efficiency :<ref>{{cite web |last1=Gaudry |first1=Pierrick |last2=Kruppa |first2=Alexander |last3=Zimmermann |first3=Paul |title=A GMP-based Implementation of Schönhage-Strassen's Large Integer Multiplication Algorithm | page = 2 | date=2007|url=https://inria.hal.science/inria-00126462/file/fft.final.pdf}}</ref>

<math> E = \frac{\frac{2N}{K}+k}{n} </math> {{mvar|N}} is bit size (the one used in <math> 2^N + 1 </math>) at outermost level. {{mvar|K}} gives <math> \frac{N}{K} </math> groups of bits,  where <math> K = 2^k </math>.

{{mvar|n}} is found through {{mvar|N, K}} and {{mvar|k}} by finding the smallest {{mvar|x}}, such that  <math>  2N/K +k \le n = K2^x </math>

If one assume efficiency above 50%, <math> \frac{n}{2} \le \frac{2N}{K}, K \le n </math> and {{mvar|k}} is very small compared to rest of formula; one get

:<math> K \le 2\sqrt{N} </math>

This means: When something is very effective; {{mvar|K}} is bound above by <math> 2\sqrt{N} </math> or asymptotically bound above by <math> \sqrt{N} </math>

===Pseudocode===

Following algorithm, the standard Modular Schönhage-Strassen Multiplication algorithm (with some optimizations), is found in overview through <ref>{{cite web | year=2014 | page=28 |last1=Lüders |first1=Christoph |title=Fast Multiplication of Large Integers: Implementation and Analysis of the DKSS Algorithm |url=https://www.researchgate.net/publication/273701188}}</ref> 
{{olist
|1= Split both input numbers {{mvar|a}} and {{mvar|b}} into n coefficients of s bits each.
Use at least {{tmath|K + 1}} bits to store them,
to allow encoding of the value {{tmath|2^{K}.}}
|2= Weight both coefficient vectors according to (2.24) with powers of {{mvar|θ}} by performing cyclic
   shifts on them.
|3= Shuffle the coefficients {{tmath|a_i}} and {{tmath|b_j}} .
|4= Evaluate {{tmath|a_i}} and {{tmath|b_j}} . Multiplications by powers of ω are cyclic shifts.
|5= Do {{mvar|n}} pointwise multiplications {{tmath|1=c_k := a_kb_k}} in {{tmath|Z/(2^K + 1)Z}}. If SMUL is used recursively,
   provide {{mvar|K}} as parameter. Otherwise, use some other multiplication function like T3MUL
   and reduce modulo {{tmath|2^{K} + 1}} afterwards.
|6= Shuffle the product coefficients {{tmath|c_k}}.
|7= Evaluate the product coefficients {{tmath|c_k}}.
|8= Apply the counterweights to the {{tmath|c_k}} according to (2.25). Since {{tmath|\theta^{2n} \equiv 1}} it follows that
   {{tmath|\theta^{-k} \equiv \theta^{n-k} }}
|9=  Normalize the {{tmath|c_k}} with {{tmath|1/n \equiv 2^{-m} }} (again a cyclic shift).
|10= Add up the {{tmath|c_k}} and propagate the carries. Make sure to properly handle negative coefficients.
|11= Do a reduction modulo {{tmath|2^{N} + 1}}.
}}
* T3MUL = Toom–Cook multiplication 
* SMUL  = Schönhage–Strassen multiplication
* Evaluate = FFT/IFFT

=== Further study ===
For implemantion details, one can read the book ''Prime Numbers: A Computational Perspective''.<ref name="crandall">R. Crandall & C. Pomerance. ''Prime Numbers – A Computational Perspective''. Second Edition, Springer, 2005. Section 9.5.6: Schönhage method, p. 502. {{ISBN|0-387-94777-9}}</ref> This variant differs somewhat from Schönhage's original method in that it exploits the [[discrete weighted transform]] to perform [[negacyclic convolution]]s more efficiently. Another source for detailed information is [[Donald Knuth|Knuth]]'s ''[[The Art of Computer Programming]]''.<ref>{{cite book |last=Knuth |first=Donald E. |title=[[The Art of Computer Programming]] |year=1997 |volume=2: Seminumerical Algorithms |edition=3rd |publisher=Addison-Wesley |isbn=0-201-89684-2 |chapter=§ 4.3.3.C: Discrete Fourier transforms |pages=305–311 |chapter-url=https://archive.org/details/artofcomputerpro0002knut_u2o0/page/305/ }}</ref>