An Algorithm to Compute the Nearest Point in the Lattice $A_{n}^*$

ROBBY G. MCKILLIAM ET AL. , NE AREST POINT IN T HE LA TTICE A ∗ N 1 An Algorithm to Compute the Nearest Point in the Lattice A ∗ n Robby G. McKilli am, I. V aughan L. Clarkson and Barry G. Quinn Abstract —The lattice A ∗ n is an important lattice b ecause of its cov ering p roperties in low dimensions. Clarkson [1] described an algorithm to comput e the nearest lattice point in A ∗ n that requires O ( n log n ) arith metic operations. In this paper , we describe a new algorithm. While the complexity is still O ( n log n ) , it is significantly simpler to describe and verify . In practice, we find that the new algorithm also ru ns f aster . Index T erms —Lattice theory , nearest point algorithm, quanti- zation, channel coding, frequency estimation, direction-of-arriva l estimation, synchronization I . I N T RO D U C T I O N T HE study of point lattices is of grea t impor tance in se veral areas of numbe r theory , par ticularly the studies o f quadra tic fo rms, the geometry of numbers and s imultaneo us Diophantin e approxim ation, and also to the pr actical engine er- ing pro blems of quantizatio n and channel cod ing. They are also imp ortant in studying the sph ere pack ing problem an d the k issing num ber problem [1], [2]. A lattice , L , is a set of points in R n such that L = { x ∈ R n | x = Bw , w ∈ Z n } where B is term ed the gene rator ma trix . The lattice A ∗ n is an inter esting lattice due to its covering proper ties in low dimension s. It g iv es th e thinn est covering in all dimen sions up to 8 [2]. A ∗ n has also fou nd app lication in a number of estimation p roblems in cluding perio d estimatio n from sparse timing data [3], fr equency estimation [4] and direction of arriv al estimation [5 ]. The near est lattice po int prob lem is: Given y ∈ R n and some lattice L whose lattice points lie in R n , find the lattice point x ∈ L such that the E uclidean distance between y and x is min imized. I f th e la ttice is u sed f or vector qu antization then the near est lattice point cor respond s to the m inimum distortion point. If the lattice is used as a cod e for a Gau ssian channel, then the nearest lattice poin t correspo nds to maximu m likelihood decod ing [6 ]. Conway and Sloan e [6 ] appear to have been the first to study the problem of computin g the nearest lattice point in A ∗ n . By decomposing A ∗ n into a unio n of translations of its dual lattice A n , they d iscovered an algor ithm for com puting the nearest lattice p oint to a g iv en p oint in O ( n 2 log n ) arithmetic Robby McKillia m is partly supported by a scholarship from the W irele ss T echnologie s Laboratory , CSIR O ICT Centre, Sydney , Australia Robby McKilliam and V augh an Clarkson are with the School of Infor- mation T echnol ogy & Electri cal E nginee ring, T he Univ ersity of Queensland, Qld., 4072, Australia Barry Quinn is with the D epartme nt of Statistic s, Macquari e Univ ersity , Sydney , NSW , 2109, Australia operation s. Later [7], they we re able to imp rove th e execution time of the algo rithm to O ( n 2 ) op erations. Clarkson [1] further improved upo n th e work of Conway and Sloan e and described an algorithm to co mpute the n earest lattice po int that req uires only O ( n log n ) arithmetic opera- tions. In this paper we describe an algo rithm that is s imilar to Clarkson’ s algo rithm. L ike Clarkson ’ s algor ithm, our algo - rithm requires O ( n log n ) a rithmetic op erations. Howev er , ou r algorithm and its de riv ation are simp ler . The new alg orithm, although o f the same ord er of co mplexity , is com putationally superior . W e now describe how the pap er is organ ized. Section II introdu ces s ome prelimin ary results and definitions. In Sec- tion III we der i ve all results ne cessary to prove that the algo - rithm does find th e nearest lattice po int. Section IV describ es the algorithm . A pseudoc ode implemen tation is provided . Is Section V, the arithmetic complexity of the algorithm is shown to b e O ( n log n ) . W e also tabulate some practical compu tation times that show the new algorithm to be com putationally superior to Clarkson’ s o riginal algorith m. I I . P R E L I M I N A RY T H E O RY V ectors an d matrices are written in b old. The i th elem ent in a vecto r is den oted by a subscr ipt: x i . T he transpose of a vector is ind icated by superscrip t T : x T . W e let 1 be a column vector of 1’ s and e i be a co lumn vector of zeros with a 1 in the i th position. The V or onoi re gio n o r near est-neighb or r e gion V ( x ) of a lattice poin t x is the subset o f R n such that, with re spect to a giv en norm, all points in V ( x ) are near er to x than to any other point in the lattice. The V oronoi regions are n dimensio nal polytop es [2 ]. The cubic lattice Z n is the set of n dimensiona l v ectors with integer elemen ts. The V oronoi region s of Z n are hyper cubes of side length 1 . The lattice A ∗ n can b e defin ed as the projectio n of the cu bic lattice Z n +1 onto th e hy perplan e o rthogo nal to 1 . Th is is, A ∗ n =  Qx | x ∈ Z n +1  (1) where Q is the projection ma trix Q =  I − 11 T n + 1  (2) where I is the ( n + 1) × ( n + 1) identity matrix. Let Π be a permu tation matrix. Observe the following elementary pr operties: 1) Π1 = 1 , ROBBY G. MCKILLIAM ET AL. , NE AREST POINT IN T HE LA TTICE A ∗ N 2 2) 1 T Π = 1 T , 3) k Πx k = k x k . Lemma 1: Th e matrices Π and Q com mute, i.e. , ΠQ = QΠ . Pr oo f: Using the pro perties of the per mutation matrix , observe that ΠQ = Π  I − 11 T n + 1  = Π − Π11 T n + 1 = Π − 11 T n + 1 = Π − 11 T Π n + 1 =  I − 11 T n + 1  Π = QΠ . Cor olla ry 1: For all z ∈ R n +1 , k Q z k = k QΠz k . Cor olla ry 2: x ∈ A ∗ n if an d only if Πx ∈ A ∗ n . Pr oo f: Because th e in verse o f a p ermutation matrix is also a permu tation matrix, we need only prove sufficiency . If x ∈ A ∗ n then x = Qk with k ∈ Z n +1 . Theref ore, Πx = ΠQk = QΠk = Qk ′ where k ′ = Πk ∈ Z n +1 and so Πx ∈ A ∗ n . Cor olla ry 3: Th e lattice p oint x is a closest p oint in A ∗ n to y if a nd only if Πx is a closest p oint in A ∗ n to Πy . Pr oo f: As for Coro llary 2, we need o nly show suf ficiency . W e do this b y co ntradiction . Suppose Πx is not clo sest to Πy but the re is instead some Πz ∈ A ∗ n such that k Π ( z − y ) k < k Π ( x − y ) k This im plies that k z − y k < k x − y k which co ntradicts the assumption that x is a closest po int to y in A ∗ n . Hence, in conside ring an a lgorithm to find a closest poin t in A ∗ n to y , it is sufficient to consider a c anonical p ermutation of y . W e will see that it is very co n venient to c onsider the permutatio n in which the (centered) fractiona l parts of y , i.e. , { y i } = y i − ⌊ y i ⌉ , are sorted in descend ing order . Th at is, in the sequ el, excep t wh ere other wise noted, we will assume that { y 1 } ≥ { y 2 } ≥ · · · ≥ { y n +1 } . (3) In the case that two or more { y i } are equal then multip le orderin gs o f y satisfy (3 ). Th e following arguments and the subsequen t algorith m are v alid fo r any ordering of y that satisfies (3 ). I I I . C L O S E S T P O I N T I N A ∗ n Lemma 2: If x = Qk is a closest p oint in A ∗ n to y ∈ R n +1 then there exists some λ ∈ R for which k is a closest point in Z n +1 to y + λ 1 . Pr oo f: Decompose y into orthogon al comp onents Qy and t 1 fo r some t ∈ R . Then k y − x k 2 = k Q ( y − k ) k 2 + t 2 ( n + 1) . (4) Observe th at Q ( y − k ) = y + λ 1 − k where we set λ = 1 T ( k − y ) n + 1 . Suppose k is not a closest point in Z n +1 to y + λ 1 . Supp ose k ′ is closer . Let x ′ = Qk ′ . Th en k y − x ′ k 2 = k Q ( y − k ′ ) k 2 + t 2 ( n + 1 ) ≤ k y + λ 1 − k ′ k 2 + t 2 ( n + 1) < k y + λ 1 − k k 2 + t 2 ( n + 1 ) = k y − x k 2 , contradictin g the assumption tha t x is a closest p oint in A ∗ n to y . Now con sider the f unction f : R 7→ Z n +1 defined so that f ( λ ) = ⌊ y + λ 1 ⌉ where ⌊·⌉ ap plied to a vector deno tes the vector in wh ich each element is round ed to a nearest integer 1 . T hat is, f ( λ ) gives a nearest point in Z n +1 to y + λ 1 as a fu nction of λ . Observe that f ( λ + 1) = f ( λ ) + 1 . Hen ce, Qf ( λ + 1 ) = Qf ( λ ) . (5) Lemma 2 implies there exists some λ ∈ R such that x = Qf ( λ ) is a closest poin t to y . Furth ermore, we see from (5 ) that λ can be fo und with in an interval of len gth 1. Hence, if we de fine the set S = { f ( λ ) | λ ∈ [0 , 1) } then Q S co ntains a closest point in A ∗ n to y . If the fr actional parts of y are sorted as in (3) , it is clear that S contain s at mo st n + 2 vectors, i.e. , S ⊆  ⌊ y ⌉ , ⌊ y ⌉ + e 1 , ⌊ y ⌉ + e 1 + e 2 , . . . , ⌊ y ⌉ + e 1 + · · · + e n +1  . (6) It can be seen that the last vector listed in the set is simp ly ⌊ y ⌉ + 1 a nd so, on ce multiplied by Q , th e fir st and th e last vector are identical. An algorithm im mediately suggests itself: test each of the n + 1 distinct vectors and find the clo sest o ne to y . Inde ed, this is exactly th e prin ciple of the alg orithm we propo se here . It only remain s to show that th is can be do ne in O ( n log n ) arithmetic oper ations. I V . A L G O R I T H M W e lab el the e lements of S acco rding to the order g iv en in (6). That is, we set u 0 = ⌊ y ⌉ and, f or i = 1 , . . . , n , u i = u i − 1 + e i . (7) Let z i = y − u i . Clearly , z 0 = { y } . Following (4), the squared distance betwee n Qu i and y is k y − Qu i k 2 = d i + t 2 ( n + 1 ) (8) where we define d i as d i = k Qz i k 2 =     z i − z T i 1 n + 1 1     2 = z T i z i − ( z T i 1 ) 2 n + 1 . (9) 1 The direction of rounding for half-inte gers is not important. Howe ver , the authors hav e chosen to round up half-inte gers in their own implementat ion. ROBBY G. MCKILLIAM ET AL. , NE AREST POINT IN T HE LA TTICE A ∗ N 3 W e kn ow th at th e near est point to y is that Qu i which minimizes ( 8). Since the ter m t 2 ( n + 1) is independ ent o f the index i , we ca n igno re it. That is, it is sufficient to minimize d i , i = 0 , . . . , n . W e n ow show that d i can be calculated inexpen si vely in a recursive fashion. W e defin e two new quan tities, α i = z T i 1 and β i = z T i z i . Fro m (7), α i = z T i 1 = ( z i − 1 − e i ) T 1 = α i − 1 − 1 (10) and β i = z T i z i = ( z i − 1 − e i ) T ( z i − 1 − e i ) = β i − 1 − 2 { y i } + 1 . (11) Input : y ∈ R n +1 z = y − ⌊ y ⌉ 1 α = z T 1 2 β = z T z 3 s = dsortindices ( z ) 4 D = β − α 2 n +1 5 m = 0 6 for i = 1 to n do 7 α = α − 1 8 β = β − 2 z s i + 1 9 if β − α 2 n +1 < D then 10 D = β − α 2 n +1 11 m = i 12 k = ⌊ y ⌉ 13 for i = 1 to m do 14 k s i = k s i + 1 15 x = k − 1 T k n +1 1 16 return x 17 Algorithm 1 : Algorithm to find a nearest l attice point in A ∗ n to y ∈ R n +1 Algorithm 1 n ow follows. The main loop beginnin g at line 7 calculates the α i and β i recursively . There is no need to retain their previous values, so the subscripts are dropp ed. The variable D maintains th e minimum value of th e (implicitly calculated values of) d i so far encoun tered, and m the co rrespon ding index. V . C O M P U TA T I O NA L C O M P L E X I T Y Each line o f the main loo p requ ires O (1) arithm etic com pu- tations so the loop (an d that on line 1 4) req uires O ( n ) in total. On line 4 the fun ction dsortindices ( z ) return s the vector s such that z s 1 ≥ z s 2 ≥ · · · ≥ z s n +1 . This sorting o peration requires O ( n log n ) arithmetic opera tions. The vector operations on lines 1 –3, 13 an d 16 all r equire O ( n ) op erations. It can b e seen, then, that the computatio nal cost of the algorithm is dom- inated by the sorting operatio n a nd is therefore O ( n lo g n ) . Clarkson’ s orig inal algo rithm re quired two sorts o f n + 1 elements. Th e new algo rithm requ ires only a single sort. See- ing as the sort domin ates the comp lexity of b oth algorith ms, we might expect our algorithm to require appro ximately half the arith metic op erations of Clarkson’ s original algorith m. This ap pears to be the case for small n . T able I shows the practical comp utational perfor mance of Clarkso n’ s algorith m versus our n ew algorith m. It is evident that the new algo rithm is comp utationally superior, p articularly for small n . It ap pears that the compu tational performance of the algorithms conv erge for large n . The computer used for these trials is an Intel C ore2 runnin g at 2.13Gh z. T ABLE I C O M P U TA T I O N T I M E I N S E C O N D S F OR 10 5 T R I A L S Algorith m n=20 n=50 n=100 n=500 Clarkson 4.57 6.97 11.11 47.81 Ne w 2.05 3.86 7.125 35.44 As a final note, the algorithm proposed here can be extended to o ther lattices for which Lemmata 1 and 2 hold. Potential candidates are the Coxeter lattices [8], [9]. R E F E R E N C E S [1] I. V . L. Clarkson, “ An algorithm to compute a nearest point in the lattice A ∗ n , ” in Applied Algebr a, Algebr aic A lgorithms and Err or-Corr ecti ng Codes , Marc Fossorier , Hideki Imai, Shu Lin, and Alain Poli, Eds., vol. 1719 of Lectur e Notes in Computer Science , pp. 104–120. Springe r , 1999. [2] J. H. Conway and N. J. A. Sloane, Sphere pac kings, lattices and gr oups , Springer , 3rd edition, 1998. [3] I. V . L. Clarkson, “ Approximate maximum-lik elihood period estimation from sparse, noisy timing data, ” IEEE T rans. Signal Pr ocess. , vol. 56, no. 5, pp. 1779–1787, May 2008. [4] I. V . L. Clarkson, “Frequen cy estimation, phase unwrapping and the nearest lat tice point problem, ” Proc. Inte rnat. Conf. Acoust. Speec h Signal Pr ocess. , vol. 3, pp. 1609–1612, 1999. [5] B. G. Quinn, “Estimating the m ode of a phase distribut ion, ” Asilomar Confer ence on Signals, Systems and Computers , pp. 587–591, Nov 2007. [6] J. H. Conway and N. J. A. Sloane, “Fast quanti zing and decoding and algorit hms for lattice quantizers and codes, ” IEEE T rans. Inform. Th. , vol. 28, no. 2, pp. 227–232, 1982. [7] J. H. Conway and N. J. A. Sloane, “Soft decoding technique s for codes and lattic es, including the Golay code and the Leech lattice, ” IEEE T rans. Inform. Th. , vol. 32, no. 1, pp. 41–50, 1986. [8] H.S.M. Coxe ter , “Extreme forms, ” Canad. J. Math. , vol. 3, pp. 391–441, 1951. [9] J. Martine t, P erfect lattices in E uclide an spaces , Springer , 2003.