Some Algorithms for Computing Short-Length Linear Convolution

Abstract

In this article, we propose a set of efficient algorithmic solutions for computing short linear convolutions focused on hardware implementation in VLSI. We consider convolutions for sequences of lengthN=2, 3, 4, 5, 6, 7, and 8. Hardwired units that implement these algorithms can be used as building blocks when designing VLSI -based accelerators for more complex data processing systems. The proposed algorithms are focused on fully parallel hardware implementation, but compared to the naive approach to fully parallel hardware implementation, they require from 25% to about 60% less, depending on the length N and hardware multipliers. Since the multiplier takes up a much larger area on the chip than the adder and consumes more power, the proposed algorithms are resource-efficient and energy-efficient in terms of their hardware implementation.

1. Introduction

Discrete convolution is found in many applications in science and engineering. Above all, it plays a key role in modern digital signal and image processing. In digital signal processing, it is the basis of filtering, multiresolution decomposition, and optimization of the calculation of orthogonal transform [1,2,3,4,5,6,7,8,9,10]. In digital image processing, convolution is a basic operation of denoising, smoothing, edge detection, blurring, focusing, etc. [11,12,13]. There are two types of discrete convolutions: the cyclic convolution and the linear convolution. General principles for the synthesis of convolution algorithms were described in [1,2,3]. The main emphasis in these works was made primarily on the calculation of cyclic convolution, while in many digital signal and image processing applications, the calculation of linear convolutions is required.

In recent years, convolution has found unusually wide application in neural networks and deep learning. Among the various kinds of deep neural networks, convolutional neural networks (CNNs) are most widely used [14,15,16,17]. In CNNs, linear convolutions are the most computationally intensive operations, since in a typical implementation, their multiple computations occupy more than 90% of the CNN execution time [14]. Only one convolutional level in a typical CNN requires more than two thousand multiplications and additions. Usually, there are several such levels in the CNN. That is why developers of such type of networks seek and design efficient ways of implementing linear convolution using the smallest possible number of arithmetic operations.

To speed up linear convolution computation, various algorithmic methods have been proposed. The most common approach to effective calculating linear convolution is dipping it in the space of a double-size cyclic convolution with the subsequent application of a fast Fourier transform (FFT) algorithm [15,16,17]. The FFT-based linear convolution method is traditionally used for large length finite impulse response (FIR) filters; however, modern CNNs use predominantly small length FIR filters. In this situation, the most effective algorithms used in the computation of a small-length linear convolution are the Winograd-like minimal filtering algorithms [18,19,20], which are the most used currently. The algorithms compute linear convolution over small tiles with minimal complexity, which makes them more effective with small filters and small batch sizes; however, these algorithms do not calculate the whole convolution. They calculate only two inner products of neighboring vectors formed from the current data stream by a moving time window of length N; therefore, these algorithms do not compute the true linear convolution.

At the same time, there are a number of CNNs in which it is necessary to calculate full-size small-length linear convolutions. In addition, in many applications of digital signal processing, there is the problem of calculating a one-dimensional convolution using its conversion into a multidimensional convolution. The algorithm thus obtained has a modular structure, and each module calculates a short-length one-dimensional convolution [21].

The most popular sizes of sequences being convoluted are sequences of length 2, 3, 4, 5, 6, 7, and 8. However, in the papers known to the authors, there is no description of resource-efficient algorithms for calculation of linear convolutions for lengths greater than four [1,4,6,21,22]. In turn, the solutions given in the literature for N = 2, N = 3, and N = 4 do not give a complete imagination about the organization of the linear convolution calculation process, since their corresponding signal flow graphs are not presented anywhere. In this paper, we describe a complete set of solutions for linear convolution of small length N sequences from 2 to 8.

2. Preliminaries

Let{_hm}and,{_xn},m=0,1,…,M−1,n=0,1,…,N−1be two finite sequences of length M and N, respectively. Their linear convolution is the sequence{_y1},i=0,1,…,M+N−2 defined by [1]:

_yi=∑n=0N−1_{hi−n xn},

where we take_hi−n=0, ifi−n<0.

As a rule, the elements of one of the sequences to be convolved are constant numbers. For definiteness, we assume that it will be a sequence{_hm}.

Because sequences{_xn}and{_hm} are finite length, then their linear convolution (1) can also be implemented as matrix-vector multiplication:

_{Y(N+M−1)×1}=_{H(N+M−1)×N XN×1}

where

_{H(N+M−1)×N}=_h0h1h0⋮_h1⋱_hM−1⋮⋯_h0hM−1⋯_h1⋱⋮_hM−1,

_XN×1=^{[_x0,_x1,…_xN−1]T},

_{Y(N+M−1)×1}=^{[_y0,_y1,…,_yN+M−2]T},

_HM×1=^{[_h0,_h1,…,_hM−1]T}.

In the future, we assume that_XN×1is a vector of input data,_{Y(N+M−1)×1}is a vector of output data, and_HM×1is a vector containing constants.

Direct computation of (3) takes MN multiplications and (M − 1)(N − 1) addition. This means that the fully parallel hardware implementation of the linear convolution operation requires MN multipliers and N + M − 3 multi-input adders with a different number of inputs, which depends on the length of the sequences being convolved. Traditionally, the convolution for which M = N is assumed as a basic linear convolution operation. Resource-effective cyclic convolution algorithms for benchmark lengths (N = 2, 3,...,16) have long been published [1,2,3,4,5,6,7,8,9]. For linear convolution, optimized algorithms are described only for the cases N = 2, 3, 4 [4,6,21,22]. Below we show how to reduce the implementation complexity of some benchmark-lengths linear convolutions for the case of completely parallel hardware their implementation. For completeness, we also consider algorithms for the sequences of lengths M = N = 2, 3, and 4.

So, considering the above, the goal of this article is to develop and describe fully parallel resource-efficient algorithms for N = 2, 3, 4, 5, 6, 7, 8. 3. Algorithms for Short-Length Linear Convolution

The main idea of presented algorithm is to transform the linear convolution matrix into circular matrix and two Toeplitz matrices. Then we can rewrite (3) in following form:

_Y(2N−1)×1=_{H(2N−1)×N XN×1}=_{HK×NH˘NHL×N}−_{0K×NHL×N01×NHK×N0L×NXN×1}

where_HK×N=_{TK(l)0K×(N−K)}and_HL×N=_{0L×(N−L)TL(r)}are matrices that are horizontal concatenations of null matrices and left-triangular or right-triangular Toeplitz matrices, respectively:

_TK×N(l)=_h00⋯0_h1h0⋯0⋮⋮⋱⋮_hK−1hK−2⋯_h0, _TL(r)=_hN−1hN−2⋯_hN−L−20_hN−1⋯_hN−L−1⋮⋮⋱⋮00⋯_hN−1,

which gives

_HK×N=_h00⋯000⋯0_h1h0⋯000⋮0⋮⋮⋱⋮⋮⋮⋱⋮_hK−1hK−2⋯_h000⋯0,

_HL×N=00⋯000⋯0⋮⋮⋱⋮00⋯0_hN−1hN−2⋯_hN−L−20_hN−1⋯_hN−L−1⋮⋮⋱⋮00⋯_hN−1.

A circulant matrix_H˘Nis a matrix of cyclic convolution H_N with rows cyclically shifted by n positions down:

_H˘N=_{IN←n HN}=_hKhK−1⋯_hK+1hK+1hK⋯_hK+2⋮⋮⋱⋮_hN−1hN−2⋯_h0h0hN−1⋯_h1⋮⋮⋱⋮_hK−1hK−2⋯_hK,

where_IN(←n)- is permutation matrix obtained from the identity matrix by cyclic shift of its columns by n positions to the left and:

_HN=_h0hN−1⋯_h1h1h0⋯_h2⋮⋮⋮⋮_hN−1hN−2⋯_h0.

The coefficients K and L are natural numbers arbitrary taken and fulfilling the dependence K + L = N − 1. These values are selected heuristically for each N separately.

The product_{H˘N XN×1}is calculated using the well-known fast convolution algorithm. The products of_{HK×N XN×1}and_{HL×N XN×1}are also calculated using fast algorithms for matrix-vector multiplication with Toeplitz matrices. We use all of the above techniques to synthesize the final short-length linear convolution algorithms with reduced multiplicative complexity.

3.1. Algorithm for N = 2

Let_X2×1=^[_x0,_x1]Tand_H2×1=^[_h0,_h1]Tbe 2-dimensional data vectors being convolved and_Y3×1=^{[_y0,_y1,_y2]T}be an input vector representing a linear convolution. The problem is to calculate the product

_Y3×1=_{H3×2 X2×1},

where

_H3×2=_h0h1h0h1.

Direct computation of (5) takes four multiplications and one addition. It is easy to see that the matrix_H3×2 possesses an uncommon structure. By the Toom–Cook algorithmic trick, the number of multiplications in the calculation of the 2-point linear convolution can be reduced [1,21].

With this in mind, the rationalized computational procedure for computing 2-point linear convolution has the following form:

_Y3×1=_{A3(2) D3 A3×2(2) X2×1}

where

_A3×2(2)=11−11, _A3(2)=11−111, _D3=diag(_h0,_h0−_h1,_h1),

_s0(2)=_h0, _s1(2)=_h0−_h1, _s2(2)=_h1.

Figure 1 shows a signal flow graph for the proposed algorithm, which also provides a simplified schematic view of a fully parallel processing unit for resource-effective implementing of 2-point linear convolution. In this paper, the all data flow graphs are oriented from left to right. Straight lines in the figures denote the data transfer (data path) operations. The circles in these figures show the operation of multiplication (multipliers in the case of hardware implementation) by a number inscribed inside a circle. Points where lines converge denote summation (adders in the case of hardware implementation) and dotted lines indicate the sign-change data paths (data paths with multiplication by −1). We use the usual lines without arrows on purpose, so as not to clutter the picture [23].

As it can be seen, the calculation of 2-point linear convolution requires only three multiplications and three additions. In terms of arithmetic units, a fully parallel hardware implementation of the processor unit for calculating a 2-point convolution will require three multipliers, one two-input adder, and one three-input adder. In terms of arithmetic units, a fully parallel hardware implementation of the processor unit for calculating a 2-point convolution will require three multipliers, one two-input adder, and one three-input adder instead of four multipliers and one two-input adder in the case of completely parallel implementation of (6). So, we exchanged one multiplier for one three-input adder.

3.2. Algorithm for N = 3

Let_X3×1=^{[_x0,_x1,_x2]T}and_H3×1=^{[_h0,_h1,_h2]T}be 3-dimensional data vectors being convolved and_Y5×1=^{[_y0,_y1,_y2,_y3,_y4]T}be an input vector represented linear convolution for N = 3. The problem is to calculate the product

_Y5×1=_{H5×3 X3×1},

where

_H5×3=_{h0h1h0h2h1h0h2h1h2}.

Direct computation of (9) takes nine multiplications and five addition. Because the matrix_H5×3 also possesses uncommon structure, the number of multiplications in the calculation of the 3-point linear convolution can be reduced too [1,4,21].

An algorithm for computation 3-point linear convolution with reduced multiplicative complexity can be written with the help of following matrix-vector calculating procedure:

_Y5×1(3)=_{A5×6(3) D6(3) A6×3(3) X3×1},

where

_A6×3(3)=111111111, _A5×6(3)=1−1−11−11−11−1−111,

_D6(3)=diag(_h0,_h1,_h2,_h0+_h1,_h0+_h2,_h1+_h2).

Figure 2 shows a signal flow graph of the proposed algorithm for the implementation of 3-point linear convolution. As it can be seen, the calculation of 3-point linear convolution requires only 6 multiplications and 10 additions. Thus, the proposed algorithm saves six multiplications at the cost of six extra additions compared to the ordinary matrix-vector multiplication method.

In terms of arithmetic units, a fully parallel hardware implementation of the processor unit for calculating a 3-point linear convolution (11) will require only six multipliers, four two-input adders, one three-input adder, and one four-input adder instead of 12 multipliers, 2 two-input adders, and 1 two-input adder in the case of fully parallel implementation of (9).

3.3. Algorithm for N = 4

Let_X4×1=^{[_x0,_x1,_x2,_x3]T}and_H4×1=^{[_h0,_h1,_h2,_h3]T}be 4-dimensional data vectors being convolved and_Y7×1=^{[_y0,_y1,_y2,_y3,_y4,_y5,_y6]T}be an input vector represented linear convolution for N = 4.

The problem is to calculate the product:

_Y7×1=_{H7×4 X4×1},

where

_H7×4=_{h0h1h0h2h1h0h3h2h1h0h3h2h1h3h2h3}.

Direct computation of (14) takes 16 multiplications and 9 addition. Due to the specific structure of matrix_H7×4 , the number of multiplication operations in the calculation of (14) can be significantly reduced.

An algorithm for computation 4-point linear convolution with reduced multiplicative complexity can be written using the following matrix-vector calculating procedure:

_Y7×1=_{A7×8(4) A8×9(4) D9(4) A9×8(4) A8×4(4) X4×1},

where

_A8×4(4)=111111−11−1111, _A9×8(4)=1⊕_H2⊕1111⊕_I3,

where_INis anidentityN×Nmatrix,_H2is the(2 × 2) Hadamard matrix, and sign,”⊕” denotes direct sum of two matrices [23,24,25].

_D9(4)=diag(_{s0(4) s1(4)},…,_s8(4)),

_s0(4)=_h0,_s1(4)=(_h0+_h1+_h2+_h3)/4,_s2(4)=(_h0−_h1+_h2−_h3)/4,

_s3(4)=(_h0−_h1−_h2+_h3)/2,_s4(4)=(_h0+_h1−_h2−_h3)/2,_s5(4)=(_h0−_h2)/2,_s6(4)=_h3,

_s7(4)=_h2,_s8(4)=_h3,

_A8×9(4)=1⊕_H2⊕−11−11⊕_I3, _A7×8(4)=111−1−11−1−11−1−111111.

Figure 3 shows a signal flow graph of the proposed algorithm for the implementation of 4-point linear convolution. So, the proposed algorithm saves 7 multiplications at the cost of 11 extra additions compared to the ordinary matrix-vector multiplication method.

In terms of arithmetic units, a fully parallel hardware implementation of the processor unit for calculating a 4-point linear convolution will require only 9 multipliers, 13 two-input adders, 2 three-input adders, and 1 four-input adder instead of 16 multipliers, 2 two-input adders, 2 three-input adders, and 1 four-input adder in the case of fully parallel implementation of (14).

3.4. Algorithm for N = 5

Let_X5×1=^{[_x0,_x1,_x2,_x3,_x4]T}and_H5×1=^{[_h0,_h1,_h2,_h3,_h4]T}be 5-dimensional data vectors being convolved and_Y9×1=^{[_y0,_y1,_y2,_y3,_y4,_y5,_y6,_y7,_y8]T}be an input vector representing a linear convolution for N = 5.

The problem is to calculate the product:

_Y9×1=_{H9×5 X5×1},

where

_H9×5=_h0h1h0⋮_h1⋱_h4⋮⋱_h0h4⋱_h1⋱⋮_h4.

Direct computation of (16) takes 25 multiplications and 16 addition. Due to the specific structure of matrix_H9×5 , the number of multiplication operations in the calculation of (16) can be significantly reduced.

Thus, an algorithm for computation 5-point linear convolution with reduced multiplicative complexity can be written using the following matrix-vector calculating procedure:

_Y9×1=_{A9×11(5) A11×13(5) A13×16(5) D16(5) A16×15(5) A15×11(5) A11×5(5) X5×1},

where

_A11×5(5)=11₀₃11−11−11−11−1111111₀₃11,_A15×11(5)=111111111111−11−1_{011×404×7I4}

and_0M×Nis a null matrix of orderM×N [23,24,25],

_D16(5)=diag(_{s0(5) s1(5)},…, _s15(5)),

_s0(5)=_h0,_s1(5)=_h1,_s2(5)=_h0,_s3(5)=(_h0−_h2+_h3−_h4)/4,_s4(5)=(_h1−_h2+h₃−_h4)/4,

_s5(5)=(3_h2−2_h1+2h₀−2_h3+3_h4)/5,_s6(5)=(−_h0+_h1−h₂+_h3),_s7(5)=(−_h0+_h1−h₂+_h3),

_s8(5)=(3_h0−2_h1+3h₂−2_h3−2_h4)/5,_s9(5)=−_h2+_h3,_s10(5)=_h1−_h2,

_s11(5)=(−_h0−_h1+4h₂−_h3−_h4)/5,_s12(5)=(_h0+_h1+h₂+_h3+_h4)/5,_s13(5)=_h4,_s14(5)=_h3,

_s15(5)=_h4,

_A16×15(5)=111111111111−1_{012×404×11I4},

_A13×16(5)=_I3⊕1111₀₃1111₀₃1111⊕_I4,

_A11×13(5)=_I3⊕1−11−1−1−1−111111111−11⊕_I3,

_A9×11(5)=1_02×6111−1−1_03×51−11−11_02×6−1−11_02×5111.

Figure 4 shows a data flow diagram of the proposed algorithm for the implementation of 5-point linear convolution. The algorithm saves 9 multiplications at the cost of 22 extra additions compared to the ordinary matrix-vector multiplication method.

In terms of arithmetic units, a fully parallel hardware implementation of the processor unit for calculating a 5-point linear convolution will require 16 multipliers, 20 two-input adders, 2 three-input adders, 1 four-input adder, and 1 five-input adder instead of 25 multipliers, 2 two-input adders, 2 three-input adders, 2 four-input adders, and 1 five-input adder in the case of fully parallel implementation of expression (16).

3.5. Algorithm for N = 6

Let_X6×1=^{[_x0,_x1,_x2,_x3,_x4,_x5]T}and_H6×1=^{[_h0,_h1,_h2,_h3,_h4,_h5]T}be 6-dimensional data vectors being convolved and_Y11×1=^{[_y0,_y1,_y2,_y3,_y4,_y5,_y6,_y7,_y8,_y9,_y10]T}be an input vector representing a linear convolution for N = 6.

The problem is to calculate the product:

_Y11×1=_{H11×6 X6×1},

where

_H11×6=_h0h1h0⋮_h1⋱_h5⋮⋱_h0h5⋱_h1⋱⋮_h5.

Direct computation of (18) takes 36 multiplications and 25 addition. We proposed an algorithm that takes only 16 multiplications and 44 additions. It saves 20 multiplications at the cost of 19 extra additions compared to the ordinary matrix-vector multiplication method.

Proposed algorithm for computation 6-point linear convolution can be written with the help of following matrix-vector calculating procedure:

Y₁₁₍₆₎=_{A11(6) A˘11(6) A11×14(6) A^14(6) A˘14(6) A14×16(6) D16(6) A16×14(6) A˘14(6) A14(6) A14×6(6) X6×1},

where

_A14×6(6)=11111_05×3031111₀₃111_02×3111,_A14(6)=_I3⊕(_H2⊗_I3)⊕_I5,

_A˘14(6)=_I3⊕_I2⊗1111−11−1⊕_I5, _A16×14(6)=_I3⊕_I2⊗11111⊕_I5,

_D16=diag(_s0(6),_s1(6),…,_s15(6)),

_s0(6)=6_h0,_s1(6)=6_h1,_s2(6)=6_h0,_s3(6)=_h0+_h3+_h4+_h1+_h2+_h5,_s4(6)=3(_h4+_h1−_h0−_h3),_s5(6)=3(_h2+_h5−_h0−_h3),_s6(6)=3(_h0+_h3)−(_h0+_h3+_h4+_h1+_h2+_h5),_s7(6)=_h0−_h3+_h4−_h1+_h2−_h5),_s8(6)=3(_h4−_h1−_h0+_h3),_s9(6)=3(_h2−_h5−_h0+_h3),_s10(6)=3(_h0+_h3)−(_h0−_h3+_h4−_h1+_h2−_h5),_s11(6)=6_h5,_s12(6)=6(−_h4+_h5),_s13(6)=6(_h3−_h4),_s14(6)=6_h4,_s15(6)=6_h5,

_A14×16(6)=_I4⊕1111_02×401×31_01×302×41111⊕_I5,_A^14(6)=_I4⊕11111⊕_I5,

_A11×14(6)=111⊕(_H2⊗_I3)⊕111111, _A˘11(6)=_{I303×40304×3}1111_{04×304×304×3I4},

and sign “⊗” denotes tensor or Kronecker product of two matrices [23,24,25].

Figure 5 shows a data flow diagram of the proposed algorithm for the implementation of 6-point linear convolution.

In terms of arithmetic units, a fully parallel hardware implementation of the processor unit for calculating a 6-point linear convolution will require 16 multipliers, 32 two-input adders, and 5 three-input adders, instead of 36 multipliers, 2 two-input adders, 2 three-input adders, 2 four-input adders, 2 five-input adders, and 1 six-input adder in the case of completely parallel implementation of expression (18).

3.6. Algorithm for N = 7

Let_X7×1=^{[_x0,_x1,_x2,_x3,_x4._x5,_x6]T}and_H7×1=^{[_h0,_h1,_h2,_h3,_h4,_h5,_h6]T}be 7-dimensional data vectors being convolved and_Y13×1=^{[_y0,_y1,_y2,_y3,_y4,_y5,_y6,_y7,_y8,_y9,_y10,_y11,_y12]T}be an input vector representing a linear convolution for N = 7.

The problem is to calculate the product:

_Y13×1=_{H13×7 X7×1},

where

_H13×7=_h0h1h0⋮_h1⋱_h6⋮⋱_h0h6⋱_h1⋱⋮_h6.

Direct computation of (20) takes 49 multiplications and 36 addition. We developed an algorithm that contains only 26 multiplications and 79 additions. It saves 23 multiplications at the cost of 43 extra additions compared to the ordinary matrix-vector multiplication method.

The proposed algorithm for computation 7-point linear convolution with reduced multiplicative complexity can be written using the following matrix-vector calculating procedure:

_Y13×1(7)=_{A13×26(7) D26(7) A26×7(7) X7×1},

where

_A13×26(7)=_{A13×15(7) A15×20(7) A20×21(7) A21×22(7) A22×25(7) A25×21(7) A21(7) A21×26(7)},

_A26×7(7)=_{A26(7) A26×28(7) A28×21(7) A21×18(7) A18×7(7)}

and

_A18×7(7)=111111111_06×3I404×2−_14×102×41−11−1_06×41111111,_A21×18(7)=_I5⊕1111−11111−1111−1111⊕_I5,

−_14×1=^{[−1,−1,−1,−1]T},

_A28×21(7)=_I7⊕11_03×40311111_04×3111111_06×31_06×31111−11_04×30411⊕_I4,

_A26×28(7)=_I10⊕1−1111−111₀₆11111111⊕_I8

_A26(7)=_I5⊕11111111₀₆₀₇11111−111⊕_I8,

_D26(7)=diag(_s0(7),_s1(7),…,_s25(7)),

_s0(7)=_h0, _s1(7)=_h2−_h1,_s2(7)=_h0−_h1,_s3(7)=_h1,_s4(7)=_h0,_s5(7)=(_h6+_h5+_h4+_h3+2_h2+_h1+_h0)/7,_s6(7)=(−_h6−2_h5+3_h4−_h3−2_h2+_h1+2_h0)/2,_s7(7)=(2_h4−_h3−2_h2+_h1)/2,_s8(7)=(−_h6+_h5+2_h4−_h3−2_h2+3_h1−_h0)/2,_s9(7)=(10_h6+3_h5−11_h4+10_h3+3_h2−11_h1−4_h0)/14,_s10(7)=(−2_h6−2_h5−2_h4+12_h3+5_h2−9_h1−2_h0)/14,_s11(7)=(2_h6+3_h5−_h4−2_h3+3_h2−_h1)/6,_s12(7)=(3_h6−11_h5−4_h4+10_h3+3_h2−11_h1+10_h0)/14,_s13(7)=(−2_h3+3_h2−_h1)/6,_s14(7)=(3_h6−_h5−2_h3+3_h2−_h1−2_h0)/6,_s15(7)=(−_h6+_h4−_h3+_h1)/6,_s16(7)=(−_h3+_h1)/6,_s17(7)=(_h5−_h3+_h1−_h0)/6,_s18(7)=2_h6−_h5−2_h4+3_h3−2_h2−2_h1+_h0,_s19(7)=2_h3−_h2−2_h1+_h0,_s20(7)=−_h6−2_h5+_h4+2_h3−_h2−2_h1+3_h0,_s21(7)=_h6,_s22(7)=_h4−_h5,_s23(7)=_h6−_h5,_s24(7)=_h5,_s25(7)=_h6,

_A21×26(7)=_I6⊕111111_03×602×91111−11−11−11_03×602×9−11−11⊕_I5,

_A21(7)=_I7⊕11111₀₄₀₄11111⊕_I6, _A25×21(7)=_I7⊕11111111_06×406×41111111⊕_I6,

_A22×25(7)=_I6⊕111₀₃1−111111_05×706×711₀₃11−111111⊕_I5,

_A21×22(7)=_I6⊕1−1−1−11111_0505×6111111⊕_I5, _A20×21(7)=_I10⊕11111⊕_I6,

_A15×20(7)=111111⊕_I5⊕1−1−1−111_04×503×5111111,

_A13×15(7)=_{I303×403×50304×3}11111111−1−1−1−11−11−111111_{0303×40303×5I3},

Figure 6 shows a data flow diagram of the proposed algorithm for the implementation of 7-point linear convolution.

In terms of arithmetic units, a fully parallel hardware implementation of the processor unit for calculating a 7-point linear convolution will require 27 multipliers, 49 two-input adders, 7 three-input adders, and 5 four-input adders, instead of 49 multipliers, 2 two-input adders, 2 three-input adders, 2 four-input adders, 2 five-input adders, 2 six-input adders, and 1 seven-input adder in the case of completely parallel implementation of expression (20).

3.7. Algorithm for N = 8

Let_X8×1=^{[_x0,_x1,_x2,_x3,_x4._x5,_x6,_x7]T}and_H8×1=^{[_h0,_h1,_h2,_h3,_h4,_h5,_h6,_h7]T}be 8-dimensional data vectors being convolved and_Y15×1=^{[_y0,_y1,_y2,_y3,_y4,_y5,_y6,_y7,_y8,_y9,_y10,_y11,_y12,_y13,_y14]T}be an input vector representing a linear convolution for N = 8.

The problem is to calculate the product:

_Y15×1=_{H15×8 X8×1},

where

_H15×8=_h0h1h0⋮_h1⋱_h7⋮⋱_h0h7⋱_h1⋱⋮_h7.

Direct computation of (22) takes 64 multiplications and 49 addition. We developed an algorithm that contains only 27 multiplications and 67 additions. Thus, the proposed algorithm saves 22 multiplications at the cost of 18 extra additions compared to the ordinary matrix-vector multiplication method.

Proposed algorithm for computation 8-point linear convolution with reduced multiplicative complexity can be written using the following matrix-vector calculating procedure:

_Y15×1=_{A15(8) A^15(8) A˘15(8) A15×17(8) A17×27(8) D27(8) A27×17(8) A17×15(8) A15×8(8) X8×1}

where

_A15×8(8)=_{I303×5I4I4I4}−_I404I4, _A17×15(8)=_I3⊕(_H2⊗_I2)⊕_I20202I2I2I2⊕_I4,

_A27×17(8)=111111⊕_H2⊕_I4⊗1111⊕11111−1111111,

_D27(8)=diag(_s0(8),_s1(8),...,_s27(8)),

_s0(8)=_h0,_s1(8)=_h2−_h1,_s2(8)=_h0−_h1,_s3(8)=_h1,_s4(8)=_h0,_s5(8)=18(_h0+_h1+_h2+_h3+_h4+_h5+_h6+_h7),_s6(8)=18(_h0−_h1+_h2−_h3+_h4−_h5+_h6−_h7),_s7(8)=14(−_h0+_h1+_h2−_h3−_h4+_h5+_h6−_h7),_s8(8)=14(−_h0−_h1+_h2+_h3−_h4−_h5+_h6+_h7),_s9(8)=14(_h0−_h2+_h4−_h6),_s10(8)=12(_h0−_h1−_h2+_h3−_h4+_h5+_h6−_h7),_s11(8)=12(_h0+_h1−_h2+_h3−_h4−_h5+_h6−_h7),_s12(8)=12(−_h0+_h2+_h4−_h6),_s13(8)=12(_h0−_h1+_h2−_h3−_h4+_h5−_h6+_h7),_s14(8)=12(_h0−_h1+_h2+_h3−_h4+_h5−_h6−_h7),_s15(8)=12(−_h0−_h2+_h4+_h6),_s16(8)=12(−_h0+_h1+_h4−_h5),_s17(8)=12(−_h0−_h3+_h4+_h7),_s18(8)=12(_h0−_h4),_s19(8)=_h4−_h6,_s20(8)=12(_h5+_h6),_s21(8)=12(_h6−_h5),_s22(8)=_h5−_h7,_s23(8)=_h7,_s24(8)=_h6,_s25(8)=_h7,_s26(8)=_h7,

_A17×27(8)=111111⊕_H2⊕_I4⊗1111⊕11111−1−1⊕111,

_A15×17(8)=_I3⊕(_H2⊗_I2)⊕_02I2I2I202I2⊕_I4, _A˘15(8)=_I3⊕(_H2⊗_I4)⊕_I4,

_A^15(8)=_I3⊕_{05×3I5I303×5}⊕_I4, _A15(8)=_I703×4⊕(−_I4)−1−1−1⊕_I4I8.

Figure 7 shows a data flow diagram of the proposed algorithm for the implementation of 8-point linear convolution.

In terms of arithmetic units, a fully parallel hardware implementation of the processor unit for calculating a 8-point linear convolution will require 27 multipliers, 57 two-input adders, 4 three-input adders and 1 four-input adder, instead of 64 multipliers, 2 two-input adders, 2 three-input adders, 2 four-input adders, 2 five-input adders, 2 six-input adders, 2 seven-point adders and 1 eight-input adder in the case of completely parallel implementation of expression (22).

4. Implementation Complexity

Since the lengths of the input sequences are relatively small, and the data flow graphs representing the organization of the computation process are fairly simple, it is easy to estimate the implementation complexity of the proposed solutions. Table 1 shows estimates of the number of arithmetic blocks for the fully parallel implementation of the short lengths linear convolution algorithms. Since a parallel N-input adder consists of N-1 two-input adders, we give integrated estimates of the implementing costs of the sets of adders for each proposed solution expressed as the sums of two-input adders. The penultimate column of the Table 1 shows the percentage reduction in the number of multipliers, while the last column shows the percentage increase in the number of adders. As you can see, the implementation of the proposed algorithms requires fewer multipliers than the implementation based on naive methods of performing the linear convolution operations.

It should be noted that our solutions are primarily focused on efficient implementation in application specific integrated circuits (ASICs). In low-power designing low-power digital circuits, optimization must be performed both at the algorithmic level and at the logic level. From the point of view of designing an ASIC-chip that implements fast linear convolution, you should pay attention to the fact that the hardwired multiplier is a very resource-intensive arithmetic unit. The multiplier is also the most energy-intensive arithmetic unit, occupying a large crystal area [26] and dissipating a lot of energy [27]. Reducing the number of multipliers is especially important in the design of specialized fully parallel ASIC-based processors because minimizing the number of necessary multipliers reduces power dissipation and lowers the cost implementation of the entire system being implemented. It is proved that the implementation complexity of a hardwired multiplier grows quadratically with operand size, while the hardware complexity of a binary adder increases linearly with operand size [28]. Therefore, a reduction in the number of multipliers, even at the cost of a small increase in the number of adders, has a significant role in the ASIC-based implementation of the algorithm. Thus, it can be argued categorically that algorithmic solutions that require fewer hardware multipliers in an ASIC-based implementation are better than those that require more embedded multipliers.

This statement is also true for field-programmable gate array (FPGA)-based implementations. Most modern high-performance FPGAs contain a number of built-in multipliers. This means that instead of implementing the multipliers with a help of a set of conventional logic gates, you can use the hardwired multipliers embedded in the FPGA. Thus, all multiplications contained in a fully parallel algorithm can be efficiently implemented using these embedded multipliers; however, their number may not be enough to meet the requirements of a fully parallel implementation of the algorithm. So, the developer uses the embedded multipliers to implement the multiplication operations until all of the multipliers built into the chip have been used. If the embedded multipliers available in the FPGA run out, the developer will be forced to use ordinary logic gates instead. This will lead to significant difficulties in the design and implementation of the computing unit. Therefore, the problem of reducing the number of multiplications in fully parallel hardware-oriented algorithms is critical. It is clear that you can go the other way—use a more complex FPGA chip from the same or another family, which contains a larger number of embedded multipliers; however, it should be remembered that the hardwired multiplier is a very resource-intensive unit. The multiplier is the most resource-intensive and energy-consuming arithmetic unit, occupying a large area of the chip and dissipating a lot of power; therefore, the use of complex and resource-intensive FPGAs containing a large number of multipliers without a special need is impractical.

Table 2 shows FPGA devices of the Spartan-3 family, in which the number of hardwired multipliers allows one to implement the linear convolution operation in a single chip. So, for example, a 4-point convolution implemented using our proposed algorithm can be implemented using a single Spartan XC3S200 device, while a 4-point convolution implemented using a naive method requires a more voluminous Spartan XC3S400 device. A 5-point convolution implemented using our proposed algorithm can be implemented using a single Spartan XC3S200A chip, while a 5-point convolution implemented using a naive method requires a more voluminous Spartan XC3S1500A chip, and so on.

Thus, the hardware implementation of our algorithms requires fewer hardware multipliers than the implementation of naive calculation methods, all other things being equal. Taking into account the previously listed arguments, this proves their effectiveness. 5. Conclusions

In this paper, we analyzed possibilities to reduce the multiplicative complexity of calculating the linear convolutions for small length input sequences. We also synthesized new algorithms for implementing these operations for N = 3, 4, 5, 6, 7, and 8. Using these algorithms reduces the computational complexity of linear convolution, thus reducing its hardware implementation complexity too. In addition, as can be seen from Figure 1, Figure 2, Figure 3, Figure 4, Figure 5, Figure 6 and Figure 7, the proposed algorithms have a pronounced parallel modular structure. This simplifies the mapping of the algorithms into an ASIC structure and unifies its implementation in FPGAs. Thus, the acceleration of computations during the implementation of these algorithms can also be achieved due to the parallelization of the computation processes.

Figure 1. The signal flow graph of the proposed algorithm for computation of 2-point linear convolution.

Figure 2. The signal flow graph of the proposed algorithm for computation of 3-point linear convolution.

Figure 3. The signal flow graph of the proposed algorithm for computation of 4-point linear convolution.

Figure 4. The signal flow graph of the proposed algorithm for computation of 5-point linear convolution.

Figure 5. The signal flow graph of the proposed algorithm for computation of 6-point linear convolution.

Figure 6. The signal flow graph of the proposed algorithm for computation of 7-point linear convolution.

Figure 7. The signal flow graph of the proposed algorithm for computation of 8-point linear convolution.

Length N	Number of Arithmetical Units (Multipliers — “×” and Adders — “+”)
	Naïve Method		Proposed Solutions		Percentage Stimate
	“×”	“+”	“×”	“+”	“×"	“+”
2	4	1	3	3	25%	66.7%
3	9	4	6	10	33.3%	60%
4	16	9	9	20	43.8%	55%
5	25	16	16	38	57.9%	66%
6	36	25	16	44	55.6%	43.2%
7	49	36	26	79	47%	54.4%
8	64	49	27	67	58%	26.9%

Length N	Features of the Implementation in Spartan-3 Family Devices
	Naïve Method	Proposed Solutions
	Type of Device	Type of Device
2	XC3S50	XC3S50AN
3	XC3S200	XC3S200
4	XC3S400	XC3S200
5	XC3S1400AN	XC3S200AN
6	XC3S2000	XC3S200
7	XC3S4000	XC3S1400AN

Author Contributions

Conceptualization, A.C.; methodology, A.C.; validation, J.P.P.; formal analysis, A.C. and J.P.P.; writing-original draft preparation, A.C.; writing-review and editing, J.P.P.; visualization, A.C. and J.P.P.; supervision, A.C. and J.P.P. Both authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Conflicts of Interest

The authors declare no conflict of interest.

References

1. Blahut, R.E. Fast Algorithms for Signal Processing; Cambridge University Press: Cambridge, UK, 2010.

2. Tolimieri, R.; An, M.; Lu, C. Algorithms for Discrete Fourier Transform and Convolution; Springer Science & Business Media: New York, NY, USA, 1989.

3. McClellan, J.H.; Rader, C.M. Number Theory in Digital Signal Processing; Prentice-Hall: Englewood Cliffs, NJ, USA, 1979.

4. Berg, L.; Nussbaumer, H. Fast Fourier Transform and Convolution Algorithms. Z. fur Angew. Math. Und Mech. 1982, 62, 282.

5. Burrus, C.S.; Parks, T. Convolution Algorithms; John Wiley and Sons: New York, NY, USA, 1985.

6. Krishna, H. Digital Signal Processing Algorithms: Number Theory, Convolution, Fast Fourier Transforms, and Applications; Routledge: New York, NY, USA, 2017.

7. Bi, G.; Zeng, Y. Transforms and Fast Algorithms for Signal Analysis and Representations; Springer Science & Business Media: New York, NY, USA, 2004.

8. Mallat, S.G. A theory for multiresolution signal decomposition: The wavelet representation. IEEE Trans. Pattern Anal. Mach. Intell. 1989, 11, 674-693.

9. Parhi, K.K. VLSI Digital Signal Processing Systems: Design and Implementation; John Wiley & Sons: New York, NY, USA, 2007.

10. Chan, Y.H.; Siu, W.C. General approach for the realization of DCT/IDCT using convolutions. Signal Process. 1994, 37, 357-363.

11. Huang, T.S. Two-dimensional digital signal processing II. Transforms and median filters. In Two-Dimensional Digital Signal Processing II. Transforms and Median Filters; Topics in Applied Physics; Springer: Berlin/Heidelberg, Germany, 1981; Volume 43.

12. Pratt, W.K. Digital Image Processing; John Wiley & Sons: Hoboken, NJ, USA, 2007.

13. Wolberg, G.; Massalin, H. A fast algorithm for digital image scaling. In Communicating with Virtual Worlds; Springer: Tokyo, Japan, 1993; pp. 526-539.

14. Krizhevsky, A.; Sutskever, I.; Hinton, G.E. Imagenet classification with deep convolutional neural networks. Commun. ACM 2017, 60, 84-90.

15. Mathieu, M.; Henaff, M.; LeCun, Y. Fast training of convolutional networks through ffts. arXiv 2013, arXiv:1312.5851.

16. Lin, S.; Liu, N.; Nazemi, M.; Li, H.; Ding, C.; Wang, Y.; Pedram, M. FFT-based deep learning deployment in embedded systems. In Proceedings of the 2018 Design, Automation & Test in Europe Conference & Exhibition (DATE), Dresden, Germany, 19-23 March 2018; pp. 1045-1050.

17. Abtahi, T.; Shea, C.; Kulkarni, A.; Mohsenin, T. Accelerating convolutional neural network with fft on embedded hardware. IEEE Trans. Very Large Scale Integr. (VLSI) Syst. 2018, 26, 1737-1749.

18. Lavin, A.; Gray, S. Fast algorithms for convolutional neural networks. In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, Las Vegas, NV, USA, 27-30 June 2016; pp. 4013-4021.

19. Wang, X.; Wang, C.; Zhou, X. Work-in-progress: WinoNN: Optimising FPGA-based neural network accelerators using fast winograd algorithm. In Proceedings of the 2018 International Conference on Hardware/Software Codesign and System Synthesis (CODES+ ISSS), Turin, Italy, 30 September-5 October 2018; pp. 1-2.

20. Cariow, A.; Cariowa, G. Minimal Filtering Algorithms for Convolutional Neural Networks. arXiv 2020, arXiv:2004.05607.

21. Wang, Y.; Parhi, K. Explicit Cook-Toom algorithm for linear convolution. In Proceedings of the 2000 IEEE International Conference on Acoustics, Speech, and Signal Processing, Istanbul, Turkey, 5-9 June 2000; Volume 6, pp. 3279-3282.

22. Ju, C.; Solomonik, E. Derivation and Analysis of Fast Bilinear Algorithms for Convolution. arXiv 2019, arXiv:1910.13367.

23. Cariow, A. Strategies for the Synthesis of Fast Algorithms for the Computation of the Matrix-vector Products. J. Signal Process. Theory Appl. 2014, 3, 1-19.

24. Regalia, P.A.; Sanjit, M.K. Kronecker products, unitary matrices and signal processing applications. SIAM Rev. 1989, 31, 586-613.

25. Granata, J.; Conner, M.; Tolimieri, R. The tensor product: A mathematical programming language for FFTs and other fast DSP operations. IEEE Signal Process. Mag. 1992, 9, 40-48.

26. Wong, H.; Betz, V.; Rose, J. Comparing FPGA vs. custom CMOS and the impact on processor microarchitecture. In Proceedings of the 19th ACM/SIGDA International Symposium on Field Programmable Gate Arrays, Monterey, CA, USA, 27 February-1 March 2011; pp. 5-14.

27. Jhamb, M.; Lohani, H. Design, implementation and performance comparison of multiplier topologies in power-delay space. Eng. Sci. Technol. Int. J. 2016, 19, 355-363.

28. Oudjida, A.K.; Chaillet, N.; Berrandjia, M.L.; Liacha, A. A new high radix-2r (r ≥ 8) multibit recoding algorithm for large operand size (N ≥ 32) multipliers. J. Low Power Electron. 2013, 9, 50-62.