{'fraction_non_alphanumeric': 0.031597172125234455, 'fraction_numerical': 0.01271719191211328, 'mean_word_length': 4.640971980002325, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Winograd Convolution for Deep Neural Networks: Efficient Point Selection

Convolutional neural networks (CNNs) have dramatically improved the accuracy of tasks such as object recognition, image segmentation and interactive speech systems. CNNs require large amounts of computing resources because ofcomputationally intensive convolution layers. Fast convolution algorithms such as Winograd convolution can greatly reduce the computational cost of these layers at a cost of poor numeric properties, such that greater savings in computation exponentially increase floating point errors. A defining feature of each Winograd convolution algorithm is a set of real-value points where polynomials are sampled. The choice of points impacts the numeric accuracy of the algorithm, but the optimal set of points for small convolutions remains unknown. Existing work considers only small integers and simple fractions as candidate points. In this work, we propose a novel approach to point selection using points of the form {-1/c , -c, c, 1/c } using the full range of real-valued numbers for c. We show that groups of this form cause cancellations in the Winograd transform matrices that reduce numeric error. We find empirically that the error for different values of c forms a rough curve across the range of real-value numbers helping to localize the values of c that reduce error and that lower errors can be achieved with non-obvious real-valued evaluation points instead of integers or simple fractions. We study a range of sizes for small convolutions and achieve reduction in error ranging from 2% to around 59% for both 1D and 2D convolution. Furthermore, we identify patterns in cases when we select a subset of our proposed points which will always lead to a lower error. Finally we implement a complete Winograd convolution layer and use it to run deep convolution neural networks on real datasets and show that our proposed points reduce error, ranging from 22% to 63%.


Introduction
Resource constrained embedded devices play an important role as edge devices in Internet-of-things (IoT) network and autonomous systems such as driverless cars. To truly enable these technologies, there is a critical need to execute resource-intensive algorithms like convolutional neural networks (CNNs) on edge devices. CNNs are used extensively for applications like facial recognition, image segmentation, speech recognition and object detection [1,2,3,4]. For such resource intensive algorithms, computing at the edge is essential because sending large quantities of data over a network to cloud based infrastructure consumes much energy and may be too slow for real-time processing of data. It is also a challenge to use resource constrained devices for such resource intensive applications due to the limited processing, memory and/or battery capacity of many edge devices [5,6] Table 1: The number of scalar multiply operations used by the Hadamard product step of 2D Winograd convolution for different input block sizes. Note that an output block size of 1 × 1 corresponds to simple direct convolution.
The computational cost of a CNN is primarily derived from convolution layers. A large array of input data is convolved with a large number of much smaller convolution kernels. A simple implementation of convolution requires O(kn) operations to convolve an input of size n with a size k kernel. So-called fast convolution algorithms can greatly reduce the operation count. For example, the input and kernel can be transformed to the Fourier domain using the fast Fourier transform (FFT) with only O(nlog 2 (n) and O(klog 2 (k)) operations respectively. In the Fourier domain, convolution requires only n + k − 1 complex-number multiply operations. Thus, the dominant term in FFT convolution is typically O(nlog 2 (n)) operations, and FFT has been successfully used for CNN convolution [7].
An alternative to FFT is Winograd's fast convolution algorithm. Like FFT convolution, Winograd's approach applies transforms to the inputs and output. However, unlike FFT, Winograd's algorithm performs the core of the convolution in n + k − 1 real-number multiply operations. Complex-number multiplication requires at least four times as many arithmetic operations * as real-number multiplication. Thus, Winograd convolution can be significantly faster than FFT convolution, particularly for the small convolution kernels found in CNNs. Indeed, Winograd's approach is optimal with respect to the number of multiply operations in the inner core of the algorithm [8]. Thus, with an input of length n and a kernel of length k, Winograd requires only a theoretical minimum of n+k−1 general multiplications (Hadamard product operations) [9,10]. Note that Lavin and Gray's seminal paper on Winograd convolution for CNNs actually uses a closely-related approach, Toom-Cook convolution. All Toom-Cook convolution algorithms have an equivalent Winograd algorithm, and in the remainder of the paper we will follow the convention of the literature of using the terms interchangeably.
In Winograd convolution, the input and kernel are sampled at a given set of points using transform matrices. This is followed by a Hadamard product (element-wise multiplication) which is then converted back to the original space by an inverse transform. In CNNs the kernels are typically small, most commonly 3 × 3 or 5 × 5, but the inputs are large. The common practice is to break the input into sub-blocks and convolve the kernel with each sub-block. Winograd gains its efficiency from computing multiple output points at once. Table 1 show the number of pairwise multiplication operations needed for various input block sizes and for kernels of size 3 × 3 and 5 × 5. Large block sizes need more multiplies in total for the block, but result in less computation per output point.
Maximizing the input block size minimizes the computation needed per point in the output. However, Winograd convolution suffers from numeric accuracy issues: the floating point (FP) error increases exponentially with the size of the convolution [11]. Thus, we cannot simply use a large block size to minimize computation. Our block size must also be small enough to keep the numerical error low. However, if we can find ways to reduce the numerical error, it will allow us to use larger block sizes, and require less computation. Our work aims to reduce numeric errors in Winograd convolution so that larger tile sizes can be used for greater efficiency.
In this work, we propose a new approach to improving the numerical accuracy of Winograd convolution, with the goal of reducing the computation and energy needed for CNNs. We make the following contributions: searching sets of rational numbers with small numerators and denominators, the smooth error curve allows us to find rational or irrational real-number interpolation points that reduce the numeric error.
• For larger block sizes, extend our approach to two variables c, d with the same form, i.e., • We evaluate our proposed symmetric points strategy and find that for most important cases for CNNs it outperforms the best known existing point selection strategy.
• We also evaluate our proposed symmetric points on state-of-the-art deep convolution networks using real world data set and show that it outperforms for the majority of cases the existing point selection strategy.
This article is organized as follows: Section 2 provides a literature review and how our work differs from it. Section 3 gives an introduction to implementing convolution in deep neural networks using fast convolution methods while Section 4 describes Toom-Cook algorithm, a form of fast convolution in some detail. The proposed point selection and all the relevant details are given in Section 5 while the results comparing our proposal to the state-of-the-art is given in Section 6. Section 7 outlines some additional experiments performed in an attempt to improve the numeric accuracy of the algorithm and Section 8 outlines future research direction. Section 9 concludes the article.

Related Work
The Toom-Cook algorithm treats the input and kernel as polynomials and computes the product of the polynomials using Lagrange interpolation. The approach was originally designed for multiplication of large multi-word integers [12,13], but it can also be used for convolution. Later, Winograd [8] generalized Toom-Cook using the Chinese remainder theorem, and proved that it is optimal with respect the number of pairwise multiplications. In Toom-Cook the transforms correspond to sampling and interpolating polynomials at points on the real number line. The choice of points determines the values in the transform matrices. By choosing appropriate points, the cost of computing the transforms can be reduced. In particular, using integer and simple fraction points can simplify the transforms by making some matrix values zero, ±1, or constants that allow fixed-point multiplication using shifts and adds. A great deal of work on Toom-Cook and Winograd convolution focuses on finding transform matrices that minimize the number of operations needed for the transforms [14,15,16]. Bodrato [17] presents a search algorithm for finding transforms that require a minimum number of integer add and shift operations.
Lavin and Gray were the first to apply Winograd (Toom-Cook) convolution to CNNs and show that Winograd convolution is faster than direct convolution on all layers of the VGG network [18]. Whereas much signal processing has traditionally been performed in fixed-point, CNNs are commonly implemented using floating point. Lavin and Gray note that Lagrange interpolation has poor numerical properties that can result in large floating point errors.
There are two main existing approaches when selecting interpolation points to minimize the error in floating point Winograd convolution. First, the Chebyshev nodes define a set of n interpolation points defined by the equation x k = cos( 2k−1 2n π), k = 1, ..., n. The Chebyshev nodes are designed to reduce Runge's phenomenon, that is inaccuracies around the boundaries of the interpolated range.
Second, Barabasz et al. [10] present a theoretical error analysis of floating point (FP) accuracy for Winograd algorithms and derive a bound on the worst-case FP error. They show that the error bound for Toom-Cook algorithm grows exponentially with the size of the convolution and propose methods to reduce the error. They propose a canonical evaluation ordering of summation based on Huffman coding and find that it reduces FP error in Winograd based convolution by 12% to 14%. They also study mixed-precision and pairwise summation algorithms. They also consider the problem of selecting suitable interpolation points to reduce the numerical error. The existing best practice in the signal processing literature was that small numerators and denominators are the best points for fixed-point Winograd convolution. Barabasz et el. empirically evaluated many combinations of these simple fractional points and found good sets of points for the most common sizes used in CNNs.
Our work differs from Barabasz et al. [10] in three main ways. First, unlike prior which considers only simple fractions as suitable points, we consider the full range of rational and irrational numbers as potential points. Second, we propose combinations of points in a specific form that causes terms in the transform matrices to cancel one another: −c, − 1 c , 1 c , c. We show that points of this form have a mostly smooth error curve, that allows us to find close to optimal values of c within a manageable search time. The points proposed by Barabasz et al. are purely empirical in nature, whereas we propose points of a specific form with a particular error distribution that allows us localize points that reduce the convolution error. Finally, we evaluate our proposed points on synthetic and real CNN data and show that they result in significantly lower errors.
Among other contributions is the one by Marques et al. which presents a Winograd-aware quantized network using two different output sizes of 4×4 and 2×2 which incorporates the numerical inaccuracies of the Winograd transformations to the learning of model parameters in a CNN [19]. Their Winograd-aware ResNet-18 network results in a 2.66× speed-up for the CIFAR-10 dataset.
Zhao et al. present a combination of the Winograd convolution and Strassen matrix multiplication algorithms to reduce computational complexity [20]. For different matrix sizes, they show that a combination of Winograd and Strassen reduces the number of multiplications as compared to conventional convolution, at the cost of some increase in addition operations. However, the number of extra addition operations is much less than the reduction in the number of multiplications. This reduction saves 75% of the execution as compared to other alternatives To address the problem of an increased numerical error when Winograd convolution algorithms are applied to large input and kernel tiles, Huang et al. present a novel way to decompose large kernel tiles into several small kernels [21] at the cost of some increased computation. This allows the Winograd algorithm to be applied to a wide range of general convolutions and achieves a speed-up of approximately 2× without affecting numerical accuracy significantly. Podili et al. present an implementation of Winograd convolution in an FPGA-based accelerator using double-buffering to reduce memory latency [24]. They also propose a data layout to reduce memory bandwidth and claim to achieve 1.2× improvement in throughput by using 3× fewer multipliers and 2× less on-chip memory. Although they claim that their proposed implementation does not impact accuracy, they do not present data to support this claim. Xiao et al. present a fused layer architecture with dynamic programming to determine the structure of fusion and an automatic tool-flow from Caffe to FPGA using Vivado HLS [25]. They fuse multiple layers in CNNs by reusing the intermediate data to save memory transfer, which is important to performance on their memory bandwidth-limited FPGA.
Vincent et al. [26] propose to reduce the condition number of the Vandermonde matrices by scaling the convolution matrices. They find that this approach allows them to reduce the error in just one specific case: when convolving a 5 × 5 kernel to create a 9 × 9 output block for AlexNet and Inception v3.

Deep Convolutional Neural Networks based on Fast Convolution
A deep convolutional neural network consists of many different layers, of which the convolutional layer is the most compute intensive. The conventional convolution layer computes an individual element of the output feature map by multiplying and accumulating the corresponding input feature map with filters. However, the convolution involved in a CNN are short with popular filter sizes of 3 × 3, 5 × 5 and 11 × 11 [1,3,4]. It was demonstrated by Winograd in [9] that the existing method proposed by Toom-Cook [12,13] generates optimal convolution algorithms in terms of minimum number of general multiplications for fixed filter kernel and input tile sizes. These convolution algorithms were also shown by Lavin and Gray [18] to be around 2× as fast as direct convolution.
The generic 2D convolution can be represented as follows: which reads a 2D input feature map of size H × W with C channels and convolves it with a bank of K filters with C channels. Each filter kernel has a dimension of R × S. The family of Winograd convolution algorithms perform this convolution by (1) transforming the input and kernel tile into the Winograd domain, (2) performing a point-wise multiplication (the actual convolution), and (3) transforming the resultant product back to the original domain using an inverse transform. Here, convolution in the spatial domain is transformed to element-wise multiplication in the Winograd domain.
The Winograd family of convolution algorithms consists of a wide variety of algorithms with different trade-offs. To generate these algorithms, one requires a polynomial which may be linear or super-linear. Linear polynomials generate algorithms that are equivalent to those generated by the Toom-Cook method; they guarantee a theoretical minimum number of element-wise multiplication operations [11]. In this paper we focus on these optimal Winograd convolution algorithms that minimize the element-wise multiplication, and are equivalent to the set of Toom-Cook convolution algorithms. A wider range of non-optimal Winograd algorithms are investigated by Barabasz and Gregg in [11].

Toom-Cook Algorithm
The Toom-Cook algorithm or method is a linear convolution algorithm based on representing convolution as a polynomial product. The linear convolution of g of size M and d of size N can be represented as:  (2), direct computation of Equation (2) is equivalent to the polynomial product [15]: To compute this polynomial product, assuming both polynomials are of degree N , Toom-Cook first evaluates each polynomial at 2N − 1 distinct points α i . This is equivalent to transforming the two polynomials to the Toom-Cook (or Winograd) domain. Point-wise multiplications of g(α i ) and d(α i ) to produce y(α i ) = g(α i )d(α i ) follows (the actual convolution). Finally, using the Lagrange interpolation, y(x) is recovered as [27,28]: The Lagrange interpolation states that given a set of n distinct points (α n ) and corresponding values of y(α n ), the polynomial y(x) can be uniquely determined. The degree of the polynomial y(x) must be less than or equal to n.
Polynomial evaluation can be determined by applying the following Vandermonde Matrix to a given polynomial [28]: where V is the linear transform matrix of the Chinese remainder theorem (CRT) which maps In this way, Toom-Cook method is a special case of the Winograd family of algorithms which are based on the Chinese remainder theorem (CRT).
Applying the inverse Vandermonde matrix (V −1 ) returns the original coefficients and thus V −1 corresponds to interpolation and can either be computed using the Lagrange's formulation of Equation (4) [28] or CRT applied on polynomials [15]. Mathematically, the whole transform can be represented as [10]: where ⊙ represents point-wise (Hadamard) multiplication.
Applying the matrix exchange theorem, we obtain the following formulation: Representing V T d as A T , V g as G and V −T as B T gives us the following formulations for 1D transform: The minimal 2D algorithm (F (m × m, r × r)) can be obtained by nesting the 1D algorithm (F (m, r)) with itself to obtain:

Toom-Cook: Matrix Construction
The transform matrices, A/A T , B/B T and G/G T , identified in Equations (8) and (9) are based on n different real points α 0 , . . . , α n−1 . These are the real-valued points where the two input polynomials are evaluated. Typically, in a conventional Winograd convolution, an input of length n is convolved with a kernel of length k to produce an output of length m = n + k − 1 using n + k − 1 general multiplications. However, since the formulations of Equations (8) and (9) involve matrix interchange, Lavin and Gray expressed their convolution in terms of the output size m which is computed using a kernel of length k and input size n = m + k − 1 using m + k − 1 general multiplications.
Matrices A T and G are Vandermonde matrices of sizes m × n and n × k, respectively. The inverse Vandermonde matrix B T is of size n × n. The Lagrange formulation of Equation (4) shows that B T matrix will have a scaling factor of N i = j =i (α i − α j ) associated with each row. This scaling factor can either be applied to the output of the pairwise multiplication of Toom-Cook convolution or embedded into the matrix G (the preferred choice) [10].
The generalized form of the matrices are thus: and is also referred to as a reducing polynomial and the choice of a good reducing polynomial affects the efficiency of the algorithm [15]. For example, if the Lagrange interpolation points are 0, ±1, then M (x) = x(x − 1)(x + 1) with m i (x) being x, x − 1 and x + 1. Based on these polynomials, the convolution operation can also be mathematically expressed as [15]: The FP error in the Toom-Cook convolution is largely dependent on the point selection and the size of input/output tiles. With regards to the interpolation points, different points give different FP error, both in rounding and representation [10]. Furthermore, larger output and filter sizes reduce the number of element-wise multiplications per computed output, but they also are a source of numerical instability [26]. The growth in the number of additions and constant multiplications in the transforms of input, kernel, and output is quadratic with the tile size. However, each input tile is convolved with many different kernels, and each kernel is convolved with many different input tiles. Thus, although the cost of the transforms is quadratic, this cost is amortized over many uses of the same transformed input tile or kernel.
Note that the point selection may have a small impact on the execution time of Winograd convolution for CNNs. Where the transform matrix contains a zero value, no arithmetic is needed to implement that part of the transform. Similarly, where an element of the transform matrix has the value ±1, only addition and no multiplication is needed to implement that part of the transform. Thus, the point selection affects both the floating point error and the number of operations needed to implement the transform matrices. Indeed zero values in the transform matrices tend to reduce both the computational cost and the arithmetic error. However, the overall computation cost of the transforms is relatively low for Winograd convolution in CNNs. As noted above, each input tile is convolved with many kernels and vice versa. Thus, the computational cost of Winograd convolution for CNNs is typically dominated by the elementwise multiplication of many pairs of input tile and kernel, not by the transforms which are performed just once for each input tile and kernel.

Modified Toom-Cook Algorithm
The Toom-Cook algorithm reduces the number of element-wise multiplications at the cost of an increase in the number of additions needed to implement the transforms. To reduce the number of additions, while keeping the number of multiplications same, a modified form of the Toom-Cook algorithm is used. It has been shown by Barabasz et al. [10], that this modified form reduces floating point (FP) error in both linear transforms and the convolution operation.
The main idea behind the modified form is to have an input of size n − 1 instead of n with the same kernel size. Essentially, we are solving for one size smaller problem with the interpolation points ranging from α 0 . . . α n−2 and the matrices A T , G and B T are constructed based on these points.
In the modified algorithm, the reducing polynomial is taken as and we compute: is a polynomial of degree one less than of y(x) of Equation (12), which is recovered by: The modified Toom-Cook algorithm is represented as: The three transform matrices are initially generated for n − 1 points and extra rows/columns are added based on the following formulations [10]: where G (n−2) , A (n−2) and B (n−2) are the matrices for original method with n − 1 points and G m(n−1) , A m(n−1) and B m(n−1) are the matrices for the modified method.
Essentially, this entails that a (n − 1)th row is added to both the G and A matrix with all zeros and a 1 at the last position. For matrix B, a final row and column are added where the last row contains all zeros bar a 1 at the last position and the last column contains consecutive coefficients of polynomial M ′ (x). Thus [10]: With the modified Toom-Cook technique resulting in less additions and also significantly reducing the FP error, we will use this technique to perform all our analysis.

Proposed Point Selection for Toom-Cook Algorithm
The overall performance of Toom-Cook convolution algorithms, in terms of FP error, is dependent on the points on which the polynomials are evaluated. For example, in order to limit additional general multiplications with m = 2 and r = 3 (F (2, 3)), the evaluation points of [0, 1, −1] generate the matrices that only involve additions, subtractions and shifts by 1 [16]. Two of these points, 0 and ∞ guarantees zeros in the three matrices of Equation (6) and thus remove any source of FP error [10].
There is no known systematic method for selecting the best points to minimize the error [26]. Barabasz  First, small integers and simple fractions such as ±1, ± 1 2 , ±2, ± 1 4 , ±4 work well because they have few significant binary digits and can be represent exactly. Second, each point c is raised to successive powers, such as c 1 , c 2 , c 3 , c 4 . . . . Point values that are close to 1.0 do not become too large or small when raised to powers. Third, and on the other hand, parts of the transform matrices depend on the difference (subtraction) between different interpolation points. If different points are too close together, the difference is a very small number which tends to result in increased numeric error. Fourth, selecting pairs of points that differ in sign, that is c and −c, and those that are reciprocal, that is c and 1 c often leads to lower numeric errors. Such pairs of points help to create zero or ±1 coefficients, while also helping with the conditioning of transform matrices [29,30]. Note that these four features are often in conflict with one another, and there is no known method for finding sets of points that are optimal with respect to the floating point error. In our proposed method, we select the following set of points in addition to 0 and ∞: for F (2, 3) (that is input block size 4, output block size 2, kernel size 3), for F (4, 3), and for F (8,3). For input points: n < 6, we empirically evaluate a subset of points for F (4, 3). For e.g., for F (3, 3) the following points are evaluated: For n > 6, we select all of c points, i.e., − 1 c , −c, c, 1 c and select a subset of d points, depending on the value of n. So, for e.g., for F (7, 3), we evaluate the following points, in addition to all c points, 0 and ∞: In addition to these points, we also construct points using a subset of those proposed by Barabasz et al. in [10]. For e.g., for F (4, 3), Barabasz et al. propose: [0, −1, 1, 1 2 , −3, ∞] and we construct the following points using these (in addition to 0 and ∞): In the following sections, we describe how selecting the proposed interpolation points results in simpler transform matrices as compared to selecting points with no such structure. Beginning with F (2, 3), we selected [0, c, −c, ∞] as the input evaluation points. As stated earlier, picking up points with opposite signs helps in getting a zero or one coefficient in the matrices. Two out of the three matrices generated by these input evaluation points are given in Equation ( 27). In order to compare, the two matrices generated if four disparate points were selected is given in Equation (26).
When we replace a and b with c and −c respectively, a number of terms in the two matrices cancel which leads to simpler computations and helps reduce the error. For e.g., −a − b is reduced to 0 in B T and 1/(a − b) is reduced to 1/ − 2c in G. Similar simplifications occur in the A T matrix.
Matrices generated for four disparate points [a, b, c, d] are given in Equations (29), (31) and (31). where When a, b, c and d are replaced with − 1 c , −c, 1 c , c, the simplification in matrices is evident. A number of values are either reduced to zero/one, or significantly simplified. Similar reductions are achieved for larger input/output sizes as well.
Note that the simplifications of matrix terms caused by our point selection may slightly reduce the computational cost of Winograd convolution for CNNs, as well as reducing the numerical error. Where the resulting transform matrices contain a zero value, no computation is needed to implement that part of the transform. Similarly, where the transform matrix contains the value ±1, only addition, not multiplication, is needed to implement part of the transform. Our proposed approach to point selection tends to result in simple values in the transform matrices, and thus it may reduce the computation cost of transforms compared to other point selections. But it is not specifically designed to maximize the number of zero and ±1 values in the transform matrices. In particular, it is sometimes possible to increase the number of ±1 values in the matrix with a different point selection.
However, the impact on execution time of slightly increasing the number of ±1 values in the transform matrices is minor. As noted in Section 4, each input tile is convolved with many kernels and vice versa, so the cost of the transforms is amortized across many uses. Thus, the largest computational component of Winograd convolution for CNNs is the element-wise Hadamard product, not the transforms. The Hadamard product of many input tiles and many kernels, and the summation of convolution results across channels is normally implemented using matrix multiplication. This part of Winograd convolution dominates the execution time, and is entirely independent of the interpolation points that are used in the transforms.

Experimental Evaluation
In this section we apply our method of point selection for Winograd convolution, and search for suitable interpolation points. The search space of floating point values is extremely large, but we find that the error curve for points in our form is smooth at a macro level, with a clear region of values that are lower than other parts of the curve. This allows us to focus our search on that region of the error curve where the error is low, and find real-valued points that reduce the reduce the numerical error. This smoothness allows us to find good interpolation points in an otherwise impossibly large search space. We evaluate our selected sets of points on two different types of dataset. First, we evaluate with uniform random values to demonstrate the value of the points in general. Second, we implement our method within the Caffe framework for deep neural networks, and evaluate our point selections on several real neural networks. Barabasz et al. in [10] propose using Huffman tree evaluation order that reduces the average error at no additional cost in computations. We also use the same evaluation order and present results based on it. Although simulations were carried out for all of simple Winograd, Winograd with Huffman and Winograd with Huffman using double precision floating point for the transforms and single precision for Hadamard product, we only mention these in cases where the conclusions reached are different than Huffman based Winograd convolution. For error calculation, different types of convolution algorithm were evaluated on 5000 randomly generated input data for each evaluation. This is done to cancel any effect of individual inputs and make the results statistically constant. L1 norm is then computed between the output of each type of convolution and a normal convolution using double precision floating point values. The selected points were then used to evaluate Winograd convolution on selected deep convolution networks with real world data.

The Error Curve
Typically, as suggested in literature, the points that reduce the error in Winograd/Toom-Cook convolution shall either be integers or simple fractions of the form a b , where both a and b are integers. However, in this work, instead of only focusing on these values, we are looking at the error curve to find optimal points that reduce the error.
In order to illustrate this, we take the example of F (4, 3) convolution and evaluate it on the points, {− 1 c , −c, 0, c, 1 c }. We vary the value of c from 1.1 to 2.5 with a step size of 10 −3 and calculate the error in both 1D and 2D convolution.    The error is plotted as a function of c and is shown in Fig. 1 for 1D convolution and in Fig. 2 for 2D convolution. All three types of convolution, simple, Huffman based and Huffman with mixed precision is shown.
Both the curve and all three different types of convolution show a similar behaviour, with minimum error values corresponding to values of c in the range of 1.6 -1.8 with error increasing either side of this range. The points considered are not in the form favoured in the literature but behave well as compared to points, for e.g., − 1 2 , −2, 0, 2, 1 2 . Before analysing points based on this methodology and those based on simple fractions, we analyse whether the points proposed by Barabasz et al. in [10] combine well with the our proposed point formulation, as shown in Equation (25), in the next section.

Searching for Optimality
Here we combine some of the points proposed in [10] with our proposal evaluation points and analyse the resulting convolution error. In other words, in addition to 0 and ∞, we fix more points while selecting a subset of our proposed points to constitute the remaining. We refer to these as base points. We present results for three cases, F (2, 3), F (3, 3) and (4,3). For the case of F (2, 3), results are shown in Table 2 Table 3 shows the same analysis for F (3,3). Choosing only from our proposed points does not result in minimum convolution error which is given when 0.5 or ±1 is chosen as the base points for 1D and 2D convolution, respectively  and gives at least a 16.67% and 20% improvement for 1D and 2D convolution, respectively. However, looking at the results of F (4, 3), our proposed points give the lowest error among all options, improving by at least 2.75% and 5.95% for 1D and 2D convolution, respectively. This difference between F (3, 3) and F (4, 3) is due to the fact that for F (3, 3) not all points from {− 1 c , −c, c, 1 c } can be chosen at one given time, reducing the opportunity for maximum cancellations and reductions to 0 and ±1.
As shown in Fig. 3, certain combination of points are better than others. Upon further investigation, we found that when combining base points with a subset of our proposed points, some combinations were always better than others. For example, as shown in Table 5 F (3, 3).

Points 1D
Value of c    Base points of [10] with Huffman summation tree, the pattern that gives minimum FP error is The patterns identified for 2D convolution are given in Table 6 for the same convolution size as those given in Table 6. The patterns are fairly similar for both 1D and 2D convolution, which suggests that the problem of finding points for 2D convolution is not a fundamentally different problem to finding good 1D points.

State-of-the-Art
Having identified the best point selections using our proposed method, we compare the floating point errors in 1D and 2D Winograd/Toom-Cook convolution against the state-of-the-art numbers given by Barabasz et al. [10]. The floating point error is defined as the L1 norm between the output of Winograd convolution and a normal convolution using double precision floating point values. Barabasz et al. analytically and experimentally evaluated a number of points and suggest optimal points that reduced the error. For this work, we experimentally evaluated hundreds of thousands of different combinations for various input polynomial evaluation points and present the result for 1D and 2D convolution in Tables 7 and 8, respectively. The first row of the table, n = 0, indicates direct convolution in single precision floating point and the error is computed using the same way as for the Winograd outputs.
Apart from the case of F (2, 3), our proposed method of picking points reduces the error compared to the best existing points [10]. We find improvements for 1D and 2D convolution ranging from 2.14% to 59.54% and from 17.5% to 35.75%, respectively. However, for the two cases of F (5, 3) and F (6, 3), we were not able to improve upon the evaluation points proposed by Barabasz et al. in [10] for either 1D and 2D convolution. It seems possible that the existing point selection for these convolution sizes may already be optimal.
This shows that selecting integer or simple fraction as polynomial evaluation points is not necessary and it is important to analyze the error curve to identify points to reduce the floating point error in 1D and 2D convolution. Table 6: Optimal pattern with minimum error for 2D convolution.

Deep Convolution Networks using Winograd Convolution
As mentioned earlier, short length Winograd/Toom-Cook convolution lends itself favourably towards convolution in deep neural networks. In order to realize the applicability of the proposed points on real deep convolution networks, Winograd convolution was implemented in Caffe v1.0. Caffe [31] is a deep learning framework which allows describing a neural network in a modular way and is developed by Berkeley AI Research (BAIR) and by community contributors. Specifically for this work, a customized distribution of Caffe developed by the Software Tools Group at Trinity College Dublin is used.
A complete convolution layer is implemented using Winograd for those layers where the kernel size is three. This is because all the previous evaluation is limited to this kernel size. However, the Winograd library of functions developed for the previous evaluation were not optimized to be used in an actual network. This limited the evaluation on real networks using only normal Winograd convolution and not the one where order of summation is based on Huffman coding. Furthermore, running a full network on all the images in the dataset proved to be prohibitively long.
Therefore, in order to present a proper evaluation, the networks are run on 1024 input images. Performance of Winograd convolution, for different evaluation points, against normal convolution is evaluated by computing the L1 norm of distances between the floating point values at the output of individual layers and normalize the norm by the number of images and layer in a particular network. The choice of 1024 images was motivated after running for a number of different input images, ranging from 2 to 1024, and seeing little difference in the L1 norm.
Two datasets are evaluated for a number of networks. For datasets, we chose the CIFAR10 and CIFAR100 dataset [32] while for the networks, we use a pre-trained GoogleNet [33], ResNet-20 [34] and SqueezeNet [35] with the CIFAR10 dataset and ResNet-20 [34] with the CIFAR100 dataset. Among the three networks, SqueezeNet is most suitable for embedded devices because of its smaller size and fewer parameters.  In future, all layers can be replaced by Winograd layers that cater to different kernel dimensions, stride factor and group convolution. It was not possible to run SqueezeNet using Winograd convolution for n = 10 since the dimension of the feature maps was too small for the given number of interpolation points.
The results are presented in Table 9 where the points are the same as those shown in Table 8. The best result, indicated by lowest L1 norm, for each network and number of Winograd input points, are highlighted in bold and underlined, for readability. A number of points proposed as part of point selection scheme in this work produce much better results than state-of-the-art, with a maximum of around 60% improvement. Furthermore, the results presented in Table 9 are fairly consistent with those shown in Table 8 for a single 2D convolution. The proposed points produce less error, not only for larger Winograd input size but also for most of the smaller input size.
The case of n = 6 is particularly interesting, because it corresponds to the set of points {0, ∞, − 1 c , −c, 1 c , c} which contains all four symmetric variants of c. This is the size of convolution where we might expect our approach to perform best, and indeed we see good results. The case of n = 10 is another sweet spot for our approach, where all four variants of c and d are included in the set of points alongside 0 and ∞.
Some of the largest reductions in error for our method arise where n = 9. It is not obvious why our method is so well suited to this case. However, when we compare the errors for n = 8 and n = 9 we see that the method that chooses simple fractions performs exceptionally poorly for n = 9. The error for n = 9 is almost ten times as large as the error for n = 8 when using simple fractions. In contrast, our approach manages to find better real-number values for this awkward-sized convolution case. Table 9 also present results for the Chebyshev nodes. Chebyshev nodes improve the conditioning of polynomial interpolation [36], which is an essential step of Winograd/Toom-Cook convolution. Looking at the numbers in Table 9, Chebyshev nodes produce good results. In fact, for GoogleNet and SqueezeNet, with n = 5, these nodes are able to achieve the lowest L1 norm and for ResNet-20 (on CIFAR100 dataset) they outperform the simple fractional points for n = 4 and n = 5. However, the L1 norm of these nodes progressively get more worse as compared to the other two schemes, as the number of input interpolation points is increased.

Other Experiments
Apart from experimenting with the proposed evaluation points, details and results about which are given in Sections 5 and 6, we experimented with the following methods in an attempt to further reduce the FP error: • Implementing Winograd transforms using Kronecker product for 2D convolution • Using pairwise summation while performing addition during the Winograd transforms • Using different base points for the two separate transforms used to implement 2D convolution The Kronecker product of two matrices is a generalization of the outer product of vectors. The Kronecker product of A and B, also referred to as direct or tensor product [37], is an (m 1 m 2 ) × (n 1 n 2 ) matrix, where A and B are m 1 × n 1 and m 2 × n 2 matrices and the product is written as C = A ⊗ B. The Kronecker product has one important property which is relevant to the 2D convolution and is given here: where Vec(X) indicates vectorized form of the matrix X.
The property shown in Equation (32) matches exactly with the terms GgG T and B T dB of Equation (9). So, instead of implementing 2D Winograd transforms using two multiple matrix multiplications, one can simply calculate the Kronecker product between G and G (similarly, B and B) and then do a vector multiplication. Another idea that led to this experiment was that multiplying the transform matrix together before doing a vector multiplication with the input/coefficient vector will result in a reduced FP error. However, the number of additions and multiplications involved in the operation defined by Equation (32) is much greater than the standard matrix multiplication of Equation (9). We found that FP errors obtained using the Kronecker product were greater than that of the original case. Similarly, we experimented with pairwise summation when multiplying the Kronecker product (G ⊗ G) with the vectorized form of g by first doing a point wise multiplication and then adding in pairs. However, the fundamental problem of the Kronecker product resulting in more multiplications and additions caused the error to remain large despite an improved summation method.
Finally, we experimented with using different sets of points for the first and second dimensions when doing 2D convolution. Typically, 2D convolution is implemented by nesting two 1D convolutions using identical points for both dimensions. We experimentally evaluated having two different sets of interpolation points for the two 1D convolutions that we use to compute 2D convolution. However, the initial results were poor so we did not experiment further with this approach.

Discussion and Future Work
This paper addresses the problem of selecting interpolation points for Winograd convolution to minimize the numerical error. Reducing the numerical error allows larger input block sizes to be used, which reduces the computational cost of Winograd convolution, as shown in Figure 1. Our proposed method allows us to abandon the existing practice of using just small integer and fractional points, and instead search the full space of real-valued points. This allows us to reduce the numerical error significantly for several important convolution sizes. As described in Section 7, we also studied other strategies to reduce the numerical error that were not successful.
The ideal solution to the problem of selecting interpolation points for Winograd convolution would be an analytic mathematical solution. The closest existing such approach is the Chebyshev nodes, which are designed to minimize the error. The Chebyshev nodes may be close to optimal in reducing the asymptotic error. However, for the smallsized convolutions found in CNNs, the complex interaction between the terms in the transform matrices result in small errors that outweigh the asymptotic advantages of the Chebyshev nodes. Given that many small interactions between point selections for constant-sized small convolutions, it seems possible that the structure of an analytic solution might be different for each small size of convolution.
One possible area of significant improvement is finding interpolation points that work well for a particular CNN. CNNs are typically trained ahead of time in large data centres, and the trained network is then deployed on edge devices. Thus, at the time that the CNN is deployed, its weights are already known. In the current paper, we search for interpolation points that reduce the numerical error for any set of weights and inputs. Inevitably, these points that are suitable across a wide range of weights and inputs are compromises between multiple conflicting goals. When the set of weights is known, fewer such compromises are necessary, and it may be possible to find interpolation points that work better for a specific set of weights.

Conclusion
Winograd and Toom-Cook convolution are efficient algorithms to compute short length convolutions that occur very regularly in deep convolutional neural networks. However, these convolution algorithms suffer from reduced floating point accuracy due to transforms to and from the Winograd/Toom-Cook domain. These convolution algorithms are based on polynomial interpolation using distinct points and research contributions have suggested using small integers or simple integers of the form a b , where both a and b are small integers to reduce the floating point computation errors. In this work, we propose a particular form of input points for the modified Toom-Cook algorithm, i.e., {− 1 c , −c, 0, c, 1 c }. We evaluate the error curve for these points and find that it is mostly smooth with a clear region of low errors. This allows us to find good values of c that reduce the numeric error without having to consider all possible floating point values. We find real-valued points whose form are not simple integers or fractions that reduce the error, especially for the case of F (4, 3) convolution with input and output tile size of n = 6 and m = 4. The reduction in error is 18% and 28% for 1D and 2D convolution, respectively.
We extended our method to the larger convolution size of F (8, 3) by repeating the pattern of F (4, 3), i.e., This allows us to achieve an error reduction of 59% and 28% for 1D and 2D convolutions, respectively compared to the best existing points. In addition to our proposed form of points, where the point 0 (and ∞) are fixed as base points, we have experimentally evaluated more base points in addition to 0 and identified that the F (3, 3) convolution benefits with a different base, achieving 20% and 36% improvement for 1D and 2D convolution.
Finally, we implemented a complete Winograd convolution layer in Caffe v1.0 and used it to run GoogleNet and SqueezeNet on the CIFAR10 dataset and ResNet-20 on both CIRAR10 and CIFAR100 dataset. We evaluated our proposed points along with points using existing strategy of only picking simple rationals and the Chebyshev nodes, using 1024 images. For evaluation, we computed the distance between output of each convolution layer of Winograd and normal convolution, i.e., L1 norm. The L1 norm was then normalized to the number of convolution layers and images. The points proposed in this work achieve reduction in L1 norm in majority of cases with improvement ranging from 22% to almost 63%.