Editing Schönhage–Strassen algorithm (section)

== Description==
This section has a simplified version of the algorithm, showing how to compute the product <math>ab</math> of two natural numbers <math>a,b</math>, modulo a number of the form <math>2^n+1</math>, where <math>n=2^kM</math> is some fixed number.  The integers <math>a,b</math> are to be divided into <math>D=2^k</math> blocks of <math>M</math> bits, so in practical implementations, it is important to strike the right balance between the parameters <math>M,k</math>.  In any case, this algorithm will provide a way to multiply two positive integers, provided <math>n</math> is chosen so that <math>ab < 2^n+1</math>.

Let <math>n=DM</math> be the number of bits in the signals <math>a</math> and <math>b</math>, where <math>D=2^k</math> is a power of two.  Divide the signals <math>a</math> and <math>b</math> into <math>D</math> blocks of <math>M</math> bits each, storing the resulting blocks as arrays <math>A,B</math> (whose entries we shall consider for simplicity as arbitrary precision integers).

We now select a modulus for the Fourier transform, as follows.  Let <math>M'</math> be such that <math>DM'\ge 2M+k</math>.  Also put <math>n'=DM'</math>, and regard the elements of the arrays <math>A,B</math> as (arbitrary precision) integers modulo <math>2^{n'}+1</math>.  Observe that since <math>2^{n'} + 1 \ge 2^{2M+k} + 1 = D2^{2M}+1</math>, the modulus is large enough to accommodate any carries that can result from multiplying <math>a</math> and <math>b</math>.  Thus, the product <math>ab</math> (modulo <math>2^n+1</math>) can be calculated by evaluating the convolution of <math>A,B</math>.  Also, with <math>g=2^{2M'}</math>, we have <math>g^{D/2}\equiv -1\pmod{2^{n'}+1}</math>, and so <math>g</math> is a primitive <math>D</math>th root of unity modulo <math>2^{n'}+1</math>.

We now take the discrete Fourier transform of the arrays <math>A,B</math> in the ring <math>\mathbb Z/(2^{n'}+1)\mathbb Z</math>, using the root of unity <math>g</math> for the Fourier basis, giving the transformed arrays <math>\widehat A,\widehat B</math>.  Because <math>D=2^k</math> is a power of two, this can be achieved in logarithmic time using a [[fast Fourier transform]].

Let <math>\widehat C_i=\widehat A_i\widehat B_i</math> (pointwise product), and compute the inverse transform <math>C</math> of the array <math>\widehat C</math>, again using the root of unity <math>g</math>.  The array <math>C</math> is now the convolution of the arrays <math>A,B</math>.  Finally, the product <math>ab\pmod{2^n+1}</math> is given by evaluating
<math display="block">ab\equiv \sum_j C_j2^{Mj}\mod{2^n+1}.</math>

This basic algorithm can be improved in several ways.  Firstly, it is not necessary to store the digits of <math>a,b</math> to arbitrary precision, but rather only up to <math>n'+1</math> bits, which gives a more efficient machine representation of the arrays <math>A,B</math>.  Secondly, it is clear that the multiplications in the forward transforms are simple bit shifts.  With some care, it is also possible to compute the inverse transform using only shifts.  Taking care, it is thus possible to eliminate any true multiplications from the algorithm except for where the pointwise product <math>\widehat C_i=\widehat A_i\widehat B_i</math> is evaluated.  It is therefore advantageous to select the parameters <math>D,M</math> so that this pointwise product can be performed efficiently, either because it is a single machine word or using some optimized algorithm for multiplying integers of a (ideally small) number of words.  Selecting the parameters <math>D,M</math> is thus an important area for further optimization of the method.