15. polynomial multiplication and Fast Fourier Transformation

ArticleDirectory

• 15. polynomial multiplication and Fast Fourier Transformation
• Preface
• Section 1 polynomial Multiplication
• Section 3, FFT

15. polynomial multiplication and Fast Fourier Transformation

Preface

ClassicAlgorithmResearchThe series has been written into chapter 15th. In this chapter, we will introduce polynomial multiplication and Fast Fourier Transformation algorithms. This blog has also introduced the Discrete Fourier Transform (refer to: 10. full understanding of the Fourier transform algorithm, top, and top 10, full understanding of the Fourier transform algorithm, bottom ), this time we start with polynomial multiplication and then introduce the principles and implementation of the FFT algorithm. At the same time, although this article involves many mathematical formulas and theorems (of course, I will try to remove some theorems or proofs irrelevant to the central content we will introduce in this article, only to ensure that the reader is easy to accept or understand), but try to ensure that the reader is easy to understand.

Some friends suggested that the algorithm should be specific, so it would be OK and there is no need to learn everything. But I can't restrain my interest. I just want to write it. When I cannot force myself not to write it, this series of classical algorithms has been written down.

Section 1 polynomial Multiplication

We know that there are two methods to represent polynomials, namely, coefficient representation and point value representation. What is a coefficient representation? For example, a (x) = 6x ^ 3 + 7X ^ 2-10x + 9, B (X) =-2x ^ 3 + 4x-5, C (x) = a (x) * B (X) is a common polynomial multiplication algorithm. The coefficient is multiplied by the coefficient. This is called the coefficient representation.

So what is the dot-value representation ?. Point value of Polynomial A (x) with N times represents a set of N point values :{(X0,Y0 ),(X1,Y1 ),...,(XN-1,YN-1 )}. All xks are different, and when K is 0, 1 ,......, N-1, y (K) = a (XK ).

A Polynomial can be expressed by many different point values. This is because of any n different points x0, X1 ,..., A collection of xn-1, can be seen as the expression of this value method. For a polynomial expressed in the form of coefficients, it is easy to calculate its point value in principle. We only need to select n different points x0, X1 ,..., xn-1, then to k = ,..., n-1, obtain a (XK), and then obtain the time required for the value of the N points according to the Kona law as O (N ^ 2 ).

However, later, you will see that if XK is cleverly selected, the point value can be used to indicate that the polynomial Multiplication can be completed in linear time. Therefore, if we can utilize the mutual transformation between the coefficient representation and the point value representation, we can accelerate the polynomial multiplication (the process will be described in detail below ), the optimum time complexity of O (N * logn) is achieved.

As mentioned above, when using a coefficient representation, that is, using a general algorithm to represent a polynomial, the multiplication of two n polynomials takes O (N ^ 2) to complete. However, in dot-value notation, the running time of Polynomial multiplication is only O (n ). What we need to do next is to use the mutual conversion between the coefficient notation and the point-value notation to implementFast Multiplication of O (N * logn) in polynomial coefficient notation.

Section 2: Fast Multiplication of Polynomials  

Before that, we have to understand two concepts,Evaluate and Interpolation. In layman's terms, the so-called evaluate (coefficient-> point value) refers to the form of a polynomial in the form of a coefficient to a point value. Interpolation (point value-> coefficient) is the inverse process of value calculation, that is, it is a known point value representation, instead, you need to convert it into a polynomial coefficient representation (using N point-value pairs can also uniquely determine a polynomial that does not exceed n-1, which is called interpolation ). The transformation between this coefficient notation and the point-value notation, that is, whether it is from the coefficient notation to the point-value notation, or from the point-value notation to the coefficient notation called interpolation, the time complexity of the two processes is O (n * logn ).

As we have mentioned above, in polynomial multiplication, if a polynomial is expressed by a coefficient notation, the complexity of Polynomial multiplication will be O (N ^ 2 ), for polynomials expressed in dot-valued notation, the complexity of Polynomial multiplication is linear complexity O (N), that is: the point value can be used to represent the polynomial multiplication in linear time. At this time, the reader can make full use of your imagination,Let us assume that the following process is used to calculate the multiplication of polynomials.(But before that, you must first expand the two polynomials A (x) and B (X) to the number of times or the polynomial with a coefficient of 2n ),We will get the expected results:

1. The coefficient notation is converted to the dot-value notation. First, a polynomial is expressed by the coefficient notation, and then the polynomial is evaluated, that is, the polynomial is changed from the coefficient notation to the point notation, and the time complexity is O (n * logn );
2. The point-value representation is used to calculate polynomial multiplication. After a polynomial is expressed in dot-valued notation, the linear complexity of the polynomial is O (n );
3. Point-value notation is converted to coefficient notation. Finally, the polynomial represented by the dot-value notation is converted into a polynomial represented by the coefficient notation, Which is O (n * logn ).

In summary, the above three operations are represented by the coefficient representation.The total time complexity of Polynomial multiplication has been reduced to O (N * logn + N * logn) = O (N * logn).

As shown in, this process is completed from left to right.Polynomial multiplication calculation process. However, there are two methods to complete this process. One is the common method mentioned in the first section, that is, the multiplication operation is directly performed on the polynomials represented by the coefficient representation, the complexity is O (n ^ 2), which is reflected in the ordinary multiplication process. There is alsoThree steps1. convert a polynomial from a coefficient representation to a point value representation (point value process); 2. Use a point value representation to complete polynomial multiplication; 3. Convert the dot-value notation to the coefficient notation (Interpolation Process ). After three steps, the total time complexity isO (N * logn). It is embodied in the three merging processes of evaluation + Pointwise multiplication + interpolation.

You have seen it. We skillfully use the linear time Multiplication Algorithm of polynomials in the form of point values to speed up the operation of Polynomial multiplication expressed in the form of coefficients.In this process, the most important thing in our work is to quickly convert a polynomial from the form of coefficients to the form of point values (evaluate, we call it FFT), and then convert it From the point value form to the coefficient form (interpolation, we call it FFT inverse) at O (N * logn) time ). The time complexity of the rapid multiplication of polynomials expressed in the form of final coefficients is O (n * logn ).. It is a tough task.

Yes, it must be noted that the item w2n is a 2n unit Compound Root. It cannot be ignored that the unit Compound Root is used in both the above two representations, namely, the coefficient representation and the point value representation, before the O (N * logn) complete the evaluation and interpolation operations within the time (for unit Compound Root, please refer to the relevant information. Because in order to ensure the plain understanding of this article, we have no intention to introduce a large number of formulas or theorems ).

Section 3, FFT

Note: a considerable amount of text in this section comes from the introduction to algorithm in the Chinese version 2 and Wikipedia.

Before introducing FFT, we need to be familiar with how the half theorem is: if n> 0 and n are even numbers, then, the square of unit root for n times is equal to N/2 n/2. We already know that by using a method called Fast Fourier Transform (FFT), we can calculate dftn (A) in the O (N * logn) Time ), if a direct method is used, the complexity is O (n ^ 2 ). FFT utilizes the special properties of the unit Compound Root.

As you can see, the FFT method uses a divide-based policy, which is called divide-based governance, that is, divide-and-conquer. In the process of multiplying two polynomials A (x) and B (X), FFT uses an even number of lower-mark coefficient in a (x) and an odd number of lower-mark coefficient, two new polynomials A [0] (X) and a [1] (x) with the number limit of n/2 are defined respectively ):

A [0] ( x ) = A 0 + A 2 x + A 4 x 2 + : : : + A n -2 x n /2-1,

A [1] ( x ) = A 1 + A 3 x + A 5 x 2 + : : : + A n -1 x n /2-1.

Note that a [0] contains the coefficients of all even numbers of the lower mark in a (the last digit of the corresponding binary number of the lower mark is 0 ), A [1] contains the coefficients of all the odd-numbered sub-numbers in a (the last digit of the corresponding binary number of the sub-mark is 1 ). The following formula is available:

In this way, the problem where a (x) is located is converted:

1) calculate the value of the polynomials A [0] And a [1] at the point where the number limit is n/2, and then

2) combine the above results.

Next, we useUnit Root n timesWN.

WNNature:

1. Periodicity,WNHas a period of N.
2. Symmetry:.

For the sake of simplicity, we will set the sequence length to be transformed belowN= 2R. Based on the symmetry of the above unit root, the summation interval can be divided into two parts when the series is obtained:

FOdd(K)AndFEven(K)There are two transformations about the n/2 points of the sequence's odd and even numbers respectively. This formula can only calculateYKFor the first n/2 points of the last n/2 points, noteFOdd(K)AndFEven(K)They are all functions with a period of n/2 and are based on the symmetry of the unit root. Therefore, there is a conversion formula as follows:

• .

In this way, an n point transformation is divided into two n/2 point transformations, so that we can continue to break down. This isKuli-Fig-Based Fast Fourier TransformationThe basic principle of the algorithm. At this time, it is not difficult to analyze the time complexity of the algorithm.O(NLogN).

Okay, I didn't expect that there are so many formulas in this article (No way, it's just a math job to do this FFT. The same is true for the basis of wavelet and Fourier analysis I bought earlier. It is a collection of mathematical formulas and theorems. However, if you can understand it better, you cannot understand it. You can only understand it a little ).

Now, let's take a simple look at the butterfly Algorithm in FFT. This article will end. As shown in:

 

Someone explained this.Butterfly Algorithm: For discrete signals of n (x power of 2) points, divide them into two sequences, namely, 0, 2, 4 ,...., 2 k (as a) and 1, 3, 5 ,...., 2k-1 (recorded as B), the first half of the Fourier transformation of N points can be derived from the Fourier transformation result is a + B * rotation factor, and the second half is the A-B * rotation factor. So the Fourier transformation of the N point is transformed into the Fourier transformation of the two n/2 point sequences respectively, for each n/2 point sequence, the first step of recursion, until there is only a two-point sequence, it will only become a simple addition and subtraction relationship. Connect the addition and subtraction relationships of these points with a line. It looks like a butterfly. Okay. For more information, see Chapter 30th about algorithms.

For example, we know that for a 4x4 matrix operation, if we use a conventional algorithm, we still need to perform 64 multiplication operations and 48 addition operations. This takes a lot of time, so in H. in 264, there is an improved algorithm (Butterfly algorithm) that can reduce the number of operations. This kind of matrix calculation algorithm is very clever. It can be completely constructed by using the integer nature and symmetry of the constructed matrix.Convert a multiplication operation to an addition operation. As shown in:

The followingCodeFrom a digital image processing bookSource code:

Void winapi FFT (complex <double> * TD, complex <double> * FD, int R) <br/> {<br/> // number of points to be converted by the payment option <br/> long count; </P> <p> // cyclic variable <br/> int I, j, k; </P> <p> // intermediate variable <br/> int bfsize, P; </P> <p> // angle <br/> double angle; </P> <p> complex <double> * w, * X1, * X2, * X; </P> <p> // calculate the number of points in the pay-as-you-go transformation <br/> COUNT = 1 <R; </P> <p> // allocate memory required for computation <br/> W = new complex <double> [count/2]; <br/> X1 = new complex <double> [count]; <br/> X2 = new complex <double> [count]; </P> <p> // calculate the weighting coefficient <br/> for (I = 0; I <count/2; I ++) <br/> {<br/> angle =-I * pI * 2/count; <br/> W [I] = complex <double> (COS (angle ), sin (angle); <br/>}</P> <p> // write time domain points to X1 <br/> memcpy (x1, TD, sizeof (complex <double>) * count); </P> <p> // use the butterfly algorithm to perform a quick token transform <br/> for (k = 0; k <r; k ++) <br/> {<br/> for (j = 0; j <1 <K; j ++) <br/> {<br/> bfsize = 1 <(r-k); <br/> for (I = 0; I <bfsize/2; I ++) <br/>{< br/> P = J * bfsize; <br/> X2 [I + P] = x1 [I + P] + X1 [I + P + bfsize/2]; <br/> X2 [I + P + bfsize/2] = (x1 [I + P]-X1 [I + P + bfsize/2]) * W [I * (1 <k)]; <br/>}< br/> X = x1; <br/> X1 = x2; <br/> X2 = x; <br/>}</P> <p> // sort again <br/> for (j = 0; j <count; j ++) <br/>{< br/> P = 0; <br/> for (I = 0; I <r; I ++) <br/>{< br/> If (J & (1 <I) <br/>{< br/> P + = 1 <(r-i-1 ); <br/>}< br/> FD [J] = x1 [p]; <br/>}</P> <p> // release the memory <br/> Delete W; <br/> Delete x1; <br/> Delete X2; <br/>}

I will continue to consider and supplement it later..