UFDC Home  |  Search all Groups  |  UF Institutional Repository  |  UF Institutional Repository  |  UF Theses & Dissertations  |  Vendor Digitized Files |   Help

# The generalized matrix product and fast Fourier transform for permutohedral aggregates

## Material Information

Title:
The generalized matrix product and fast Fourier transform for permutohedral aggregates
Physical Description:
ix, 71 leaves : ; 29 cm.
Language:
English
Creator:
Zapata, Jaime L
Publication Date:

## Subjects

Subjects / Keywords:
Mathematics thesis, Ph.D   ( lcsh )
Dissertations, Academic -- Mathematics -- UF   ( lcsh )
Genre:
bibliography   ( marcgt )
theses   ( marcgt )
non-fiction   ( marcgt )

## Notes

Thesis:
Thesis (Ph.D.)--University of Florida, 2000.
Bibliography:
Includes bibliographical references (leaves 69-70).
General Note:
Printout.
General Note:
Vita.
Statement of Responsibility:
by Jaime L. Zapata.

## Record Information

Source Institution:
University of Florida
Rights Management:
All applicable rights reserved by the source institution and holding location.
Resource Identifier:
aleph - 023016300
oclc - 45080607
System ID:
AA00020424:00001

Full Text

THE GENERALIZED MATRIX PRODUCT
AND FAST FOURIER TRANSFORM
FOR PERMUTOHEDRAL AGGREGATES

By

JAIME L. ZAPATA

A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY

UNIVERSITY OF FLORIDA

2000

To my sister, Luz Helena. In her presence I get to be everything that is good,

powerful, and worthy.

ACKNOWLEDGEMENTS

First of all, I would like to thank Dr. Ritter, friend and mentor, for introducing

me to image algebra, for giving me the opportunity to do independent research, and
for serving as chair of my doctoral committee. Thanks to the other members of my

committee: Dr. David Wilson, Dr. Joseph Wilson, Dr. James Keesling, and Dr.

Murali Rao. I am indebted to them for their support and encouragement and for
everything they taught me. Special thanks to Dr. Andrew Vince for carefully going

over everything I had written and for helping me to simplify the presentation of ideas.
Thanks to the faculty of the Department of Mathematics at the University of
Florida. Nearly everyone was a teacher or friend at one point or another. Special
thanks to Dr. Beverly Brechner for giving me a deep appreciation of topology; to
Dr. Neil White for listening at a time when I was talking about the protein folding

problem; to Dr. Gang Bao for his encouragement and mentorship; and to Mrs.

Townsend-I enjoyed working with you.
More thanks to Sandy Gagnon, Sharon Easter, Sonia Simmons, and Julia

Porchiazzo, all with the staff at the Department of Mathematics. Very special thanks

to Gretchen Garrett. All of you helped me with the little things and the big things

that make a graduate career possible.
Thanks to my parents for creating an open-ended possibility that bears my
name, for their unconditional love and support, and for being such powerful witnesses.

Thanks to Lily, Rosana, and Luz Helena for witnessing me as a brother and a scholar,

among other things.

Finally, I would like to thank Haylee for coming into my life at the right time,

for making me laugh, for her love and help during the difficult times, and for having

the courage to say goodbye.

ACKNOWLEDGMENTS ............................. iii

LIST OF FIGURES ................................ vii

ABSTRACT .................................... viii

CHAPTERS
1 INTRODUCTION ............................. 1

1.1 Hexagonal and Truncated Octahedral Aggregates ........... 4
1.2 Overview ............................ ... 5
1.3 Summary of Results ......................... 6

2 THE DISCRETE FOURIER TRANSFORM .............. 7

2.1 Introduction ............... ... .. .. .. ..... 7
2.2 Vector Spaces Over Finite Abelian Groups .............. 7
2.3 The Discrete Fourier Transform and Inverse DFT ........ 11
2.4 The Convolution Theorem .......................... 17
2.5 The Cooley-Tukey Factorization .................... 20

3 QUOTIENT LATTICES ......................... 23

3.1 Introduction ............................. 23
3.2 The DFT and Inverse DFT for Quotient Lattices .......... 23
3.3 The Cooley-Tukey Factorization for Quotient Lattices ..... 28
4 PERMUTOHEDRAL AGGREGATES ................... 31

4.1 Introduction ............................. 31
4.2 Geometric Description ........... ..................... 32
4.3 Quotient Lattices Associated with Permutohedral Aggregates . 35
4.4 Aggregates as Sets of Coset Representatives . . . . .... ..41
4.5 The Cooley-Tukey Factorization for Permutohedral Aggregates 48
5 THE GENERALIZED MATRIX PRODUCT .... . . . . .. 50

5.1 Introduction ....................................50
5.2 Digital Indexing Sets . . . . . . . . . . ... .. 51

5.3 The Generalized Matrix Product for Permutohedral Aggregates 52
5.4 Digital Addressing for Permutohedral Aggregates . . . ... ..54

6 DIGITAL ALGORITHMS FOR PERMUTOHEDRAL AGGREGATES 57

6.1 Introduction . . . . . . . . . . . . . ... .. 57
6.2 The Discrete Fourier Transform and Inverse DFT . . . ... ..58
6.3 Cooley-Tukey Factorization Algorithms . . . . . . . 59
6.4 The Fast Fourier Transform and Inverse FFT Algorithms . 64

7 SUMMARY AND SUGGESTIONS FOR FURTHER RESEARCH . 67

REFERENCES . . . . . . . . . . . . . . . . ... .. 69

BIOGRAPHICAL SKETCH . . . . . . . . . . . . ... .. 71

LIST OF FIGURES

4.1 Level 1 hexagonal aggregate . . . . . . . . . . ... .. 33

4.2 Hierarchical structure of a level 3 hexagonal aggregate . . ... .. 34

4.3 Rhombic dodecahedron . . . . . . . . . . .... .. 35

4.4 Truncated octahedron . . . . . . . . . . . .... ..36

4.5 Level 1 truncated octahedral aggregate . . . . . . . ... .. 37

4.6 Level 1 hexagonal aggregate in the spatial domain . . . ... .. 45

4.7 Level 2 hexagonal aggregate in the spatial domain . . . ... .. 46

4.8 Level 2 hexagonal aggregate in the frequency domain . . . ... .. 47

5.1 Indexed level 3 hexagonal aggregate in the spatial domain ..... .. 56

Abstract of a Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy

THE GENERALIZED MATRIX PRODUCT
AND FAST FOURIER TRANSFORM
FOR PERMUTOHEDRAL AGGREGATES

By

Jaime L. Zapata

May 2000

Chairman: Dr. Gerhard X. Ritter
Major Department: Mathematics

Permutohedral aggregates are hierarchical arrangements of cells. Special cases

include hexagonal aggregates and truncated octahedral aggregates in dimensions 2

and 3, respectively. These cells can be efficiently addressed using a scheme known as

generalized balanced ternary, or GBT. The cells of a D-dimensional permutohedral

aggregate are addressed, or indexed, using finite strings of a set of n = 2D+1 1 digits.
The objects of interest in this work are functions whose domains are permutohedral

aggregates.

The generalized matrix product, also known as the p-product, is a generaliza-

tion of standard matrix multiplication.

The present work defines the discrete Fourier transform and the p-product for
permutohedral aggregates, and derives a radix-n, decimation-in-space, fast Fourier

transform (FFT) for permutohedral aggregates in terms of the p-product. Data
reordering (also known as shuffle permutations) is generally associated with FFT

algorithms. However, use of the p-product makes data reordering unnecessary.

CHAPTER 1
INTRODUCTION

The goal of this dissertation is to describe algorithms that efficiently compute
the discrete Fourier transform and the corresponding inverse transform for functions

with discrete domains known as permutohedral aggregates. The Fourier transform is

one of the most fundamental transforms in signal analysis and image processing. This

transform extracts the frequency content of a function. Permutohedral aggregates are

hierarchical arrangements of certain points in D-dimensional Euclidean space. These

aggregates are of interest because they possess several useful properties:

1. they can be addressed efficiently using finite strings of digits;

2. they have a ring structure and addition and multiplication can be done entirely

in terms of bit string operations;

3. in dimension two, aggregates approximate the region of light-sensing cells in

the human visual system;

4. in dimensions two and three, aggregates represent optimal sampling strategies

in the spatial domain.

Algebraically, A-level permutohedral aggregates are sets of coset representa-
tives of quotient lattices of the form L/NAL, where L is a D-dimensional lattice
having a ring structure, N is an endomorphism of L, A is a positive integer, and

N'L = NA(L). These quotient lattices consist of n' elements, where n = 2D+1 1.

The Voronoi cell of a point x in a lattice I is the set {y C RD : |y xa <

ly zJ, for all z E I}. The Voronoi regions or cells associated with permutohe-
dral aggregates are known as permutohedra. From an image processing perspective,
the Voronoi cells associated with aggregate elements play the role of picture elements,
or pixels. Important special cases include hexagonal cells for D = 2, and truncated
octahedral cells for D = 3. The cells of a hexagonal aggregate are centered on a
hexagonal lattice and the cells of a truncated octahedral aggregate are centered on a
body-centered cubic lattice. The simplest (level 1) hexagonal aggregate can be visual-
ized as a hexagonal cell with six neighboring cells, while a level 2 hexagonal aggregate
consists of 49 cells. A central cell in a truncated octahedral aggregate has fourteen
neighbors. See the figures in Chapter 4. The cells in a A-level permutohedral aggre-
gate can be addressed using strings of length A from the set of digits {0, 1,..., n 1}.
This form of addressing, or indexing, is known as generalized balanced ternary (GBT),

and was developed by Gibson and Lucas [8].
The majority of image processing techniques and algorithms can be broadly
classified as either spatial domain or frequency domain techniques. Spatial domain

techniques include morphological transforms, while frequency domain techniques in-
clude highpass and lowpass filtering (see Gonzalez and Woods [9], or Ritter and
Wilson [18]). Frequency domain techniques depend on the discrete Fourier trans-
form (DFT), a linear transform that extracts the frequency content of images with
discrete domains. This dissertation is concerned with the computation of the dis-
crete Fourier transform and the inverse DFT of functions, or images, whose domains
are permutohedral aggregates. The geometric duals of permutohedral aggregates are
function domains in the frequency space.
The fast Fourier transform (FFT) is an efficient algorithm for the computation
of the discrete Fourier transform. Gilbert Strang [19] calls the FFT "the most valuable

numerical algorithm in our lifetime." The complexity of the FFT is O(M log M),
where M is the cardinality of the discrete image domain. Mersereau and Speak [15],

and Dudgeon and Mersereau [7] describe the FFT for periodic functions defined on

multidimensional lattices, and Tolimieri, An, and Lu [20, 21], describe the FFT for
functions defined over Zn, the ring of integers modulo n. The FFT is based on the

Cooley-Tukey factorization of the discrete Fourier transform, which was published in

1965 [5]. However, some unpublished manuscripts by Gauss contained the essentials

of the Cooley-Tukey factorization and the FFT algorithm around 1805 (see Heideman

et al. [11]).
Data reordering, also known as shuffle permutations, is generally associated
with FFT algorithms. Zhu and Ritter [24, 26], derived algorithms for the computa-
tion of the FFT of images defined on rectangular arrays that did not require a data
reordering stage. They accomplished this by using the generalized matrix product, a
generalization of standard matrix multiplication. This product, also known as the

p-product, was first described by Ritter [17]. The present work defines the gener-

alized matrix product for arrays of data indexed by the GBT addressing scheme.

Such a definition is necessary for the description of image processing algorithms for

permutohedral aggregates and their duals.

This dissertation exploits the hierarchical structure of permutohedral aggre-
gates and the GBT addressing scheme to derive pyramidal data processing algorithms

that compute the FFT and the inverse FFT for aggregates. These algorithms have

complexity O(M log M), where M = n' and n = 2D+l 1; they are expressed in
terms of the generalized matrix product for permutohedral aggregates.

1.1 Hexagonal and Truncated Octahedral Aggregates

We include some further remarks regarding hexagonal and truncated octahe-

dral aggregates since properties of these two special cases of permutohedral aggregates

motivated the present study.

Watson and Ahumada, Jr [23] derived a mathematical transform that models

space-frequency coding of images in the primary visual cortex of the brain. They

noted that while retinal ganglion cells represent visual images with a spatial code, the

cells of the primary visual cortex code spatial, frequency, and orientation information

about an image. Furthermore, ganglion cell receptive fields form an approximate

hexagonal lattice near the fovea, the region in the central portion of the retina that

allows human beings to resolve fine details in an image. The input lattice for the

image transform presented by Watson and Ahumada, Jr. is a hexagonal array that

is equivalent to the hexagonal aggregates described here.
Petersen and Middleton [16] provide further motivation for the study of multi-

dimensional functions sampled on hexagonal and body-centered cubic lattices. They

showed that the most efficient sampling lattices are not in general rectangular. Ef-

ficient sampling lattices attempt to use a minimum number of sampling points to

exactly reconstruct bandlimited functions. In particular, sampling two-dimensional

isotropic functions on a hexagonal lattice and sampling three-dimensional isotropic

functions on a body-centered cubic lattice is significantly more efficient than sam-

pling such functions on simple square and cubic lattices, respectively. The defining
property of an isotropic function is that the Fourier transform of the function cuts off

at the same magnitude in all directions. Thus, the support of the Fourier transform

of an isotropic function is a hypersphere or ball. Petersen and Middleton specify
minimum sampling lattices for isotropic functions in dimensions one through eight.

1.2 Overview

Chapter 2 presents a number of standard definitions and results from finite

abelian group harmonic analysis: the discrete Fourier transform, the inverse DFT,

the inversion theorem, the convolution theorem, and the Cooley-Tukey factorization

theorem.

Chapter 3 covers the discrete Fourier transform, the inverse DFT, and the

Cooley-Tukey factorizations for finite abelian groups that arise as quotient lattices

and their geometric duals, respectively.

Permutohedral aggregates and the quotient lattices AN associated with per-

mutohedral aggregates are described in Chapter 4. Permutohedral aggregates are
particular sets of coset representatives of the quotients A\. Chapter 4 also discusses

the subgroups B# in A\, which lead to quotients of the form Qa = A\/BO, and con-

tains statements of the Cooley-Tukey factorizations for the quotient lattices A,\ and

the geometric duals A', respectively.

Chapter 5 defines F-indices-lexicographically ordered finite strings of digits-

and the generalized matrix product over F-indices. The F-indices are used to address

the elements of permutohedral aggregates and their duals.

Two theorems in Chapter 6 are formulations of the DFT and the inverse DFT

for permutohedral aggregates in terms of the generalized matrix product over F-
indices, or digital indexing sets. Chapter 6 also states four main digital algorithms in

terms of the generalized matrix product for permutohedral aggregates: the Cooley-
Tukey (CT) factorization, the inverse CT factorization, the fast Fourier transform,

and the inverse FFT.

Chapter 7 contains suggestions for further research.

1.3 Summary of Results

The main results established in this work are as follows:

1. algebraic formulation of the discrete Fourier transform (DFT) and the inverse

DFT for permutohedral aggregates;

2. definition of the generalized matrix product for permutohedral aggregates;

3. discrete Fourier transform and inverse DFT algorithms in terms of the gener-

alized matrix product for permutohedral aggregates;

4. Cooley-Tukey factorization and inverse CT factorization algorithms in terms of

the generalized matrix product for permutohedral aggregates;

5. a decimation-in-space, fast Fourier transform (FFT) algorithm for permutohe-

dral aggregates;

6. a decimation-in-frequency, inverse FFT algorithm for permutohedral aggre-

gates.

CHAPTER 2
THE DISCRETE FOURIER TRANSFORM

2.1 Introduction

The purpose of this chapter is to present background material from finite
abelian group harmonic analysis. Section 2.2 describes the normed vector space L'(A)
of complex-valued functions on a finite abelian group A, the algebraic dual A, and
the space L2(A). Section 2.3 presents some standard definitions and results for the

discrete Fourier transform and the inverse DFT for finite abelian groups. We recall

that the Fourier transform is an isometry and state the inversion theorem. Section 2.4

defines cyclic convolution for functions in L(A) and presents the convolution theorem

for finite abelian groups. In Section 2.5 we prove the Cooley-Tukey factorization
theorem.
Denote the set of positive integers by Z+ and the field of complex numbers

byC.

2.2 Vector Spaces Over Finite Abelian Groups

Denote by card(S) the cardinality of a set S. Let A be a finite abelian group.
The set L(A) of all complex-valued functions defined on A is a vector space over C
with addition and multiplication given by

(a + b)(r) = a(r) + b(r), for all a,b E L(A),r E A,

(Ca)(r) = Ca(r), for all a E L(A), E C,r E A.

For a,b E L(A), define ab E L(A) by (ab)(r) = a(r)b(r), for all r E A.
Define an inner product on L(A) as follows. For a, b E L(A), let

(a,b) = Z a(r)b(r),
rEA
where b(r) denotes the complex conjugate of b(r). For a E L(A), define the norm of
a by
/ \ 1/2
ail = V(a,a) = ( Ja(r)2 12

Denote by L2(A) the inner product space L(A) together with the norm ].
We now describe a basis of L2(A). For r, s C A, define the Kronecker delta
over indices in A by
3 1 if r =s
ur= 0 if r s.

For r E A, define er E L(A) by er(s) = Sr,, for all s E A. We will also have occasion
to refer to the support of a function, defined as follows. If f is any function with
domain A, R is any ring, and f : A -4 R, then the support of f, written supp(f), is
defined by supp(f) = {d E A f(d) 0}. Thus, supp(er) = r, for all r E A.

Proposition 1 The set C = {er r E A} is an orthonormal basis of L2(A).

Proof. If a E L2(A), then
a = Z a(r)er.
rEA
Hence, the set of vectors {er : r E A} spans the space L 2(A).
If r and s are distinct elements of A, then

(er, e.) = E er(t)e.(t) = e,(s) + e.(r) = J'. + S.r = 0.
tEA
Thus, is a set of (nonzero) orthogonal vectors, which implies that it is a linearly
independent set. Since C also spans L2(A), it is a basis of L2(A).

If r = s, then

(er,e) = er,) = Zer() e r(() = 1.
tEA
Hence, is an orthonormal basis of L2(A). D

Refer to as the standard basis of L2(A). If card(A) = n, then L2(A) is an
n-dimensional normed vector space, since card(E) = card(A). Section 2.3 describes
a second basis of L2(A).
Let U = {z E C: jzj = 1}, the unit circle, and for n E Z+, let U(n) denote
the multiplicative group of all the nth roots of unity. Define the algebraic dual A
of a finite abelian group A of cardinaliy n by A = Hom(A, U(n)), the set of all
homomorphisms from A into U(n). Thus, 4 E A implies 0 : A -+ U(n) such that
4)(r + s) = 4)(r)}4(s), for all r, s E A. Define addition in A as follows: for ,), b E A, let
(0 + 0)(r) = O(r)t(r), for all r E A. Note that kOc(r) = O(r)k, for all k E Z, r E A.
This addition renders A an abelian group. Denote by 0o the zero element or trivial
homomorphism in A. Hence, 00o(r) = 1, for all r E A.
A standard result for finite abelian groups follows.

Theorem 2 If A is a finite abelian group with card(A) = n and A = Hom(A, U(n)),
then A A. In particular, card(A) = n.

The elements of the dual group A enjoy the following properties.

Proposition 3 If A is a finite abelian group and 4) E A, then

1. 0(O) = 1;

2. I1)(r)12 = O(r)O(r) = 1, for all r E A;

3. The additive inverse -4) of 4) satisfies -O(r) = 4)(r)-' = O)(r), for all r E A.

Proof. Any homomorphism 4): A -* U(n) maps the identity element 0 E A to 1,
the identity element of U(n). Hence, 0(0) = 1,for all E A. which proves (1).
If r E A, then ](r)l2 is the norm of O(r) E U(n). The elements of U(n) lie
on the unit circle. Hence, I|O(r)l2 = 1,for all r E A, which proves (2).
If -0 is the additive inverse of c A, then (0 + (-O))(r) = Oo(r) =
1,for all r E A. Since (0 + (-))(r) = ((r))(-(r)) = 1, -O(r) is the multi-
plicative inverse of O(r) E U(n). Hence, -O(r) = 0(r)-1 and 0(r)-1 = O(r), for all
r E A, since -.O(r) E U, for all r E A, which proves (3). 0

Note that E A implies : A -+ U(n) C C, so that 4 L(A), for all ) E A.
Thus, A C L(A).
The set L(A) = {a: A -+ C} of all complex-valued functions on A is a vector

space over C with addition and scalar multiplication given by

(a+ b)(0) =a(0) +b(O), for all a,bE L(A),4 E A,
(-)(O) = Ca(), for all a E L(A), E C, E A.

For a,b E L(A), define a. b E L (A) by (. b)(0) = a()b()), for all E A.
Define an inner product and a norm on L(A) as follows. For a,b 6 E A, the inner
product of a and b is given by

(a,b) = a(4)b4),

and the norm of a is given by
1/2
||a||=^2)= ja(0)|12

In Section 2.3 we define a basis of L(A) with cardinality n. Denote by L2(A)
the inner product space L(A) together with the norm I1 ||. Thus, L2(A) is a normed
vector space over C of dimension n.

In order to define the Fourier transform for finite abelian groups in Section 2.3
we need a bilinear form (.,.): A x A -* U(n).

Proposition 4 If A is a finite abelian group of cardinality n, then the map

(.):A x A -* U(n)

defined by
(r,) =(r),

for all r E A and 0 E A, is a bilinear form, that is,

1. (r + s, b) = (r, 0)(s, ), for all r, s E A, E A;

2. (r, + 0) = (r, 4)(s, 4), for all r E A, EA;

3. (kr, ) = (r,))k = (r,ko), for all r E A,4) A.

In particular, Proposition 4 implies (-r,4) = (r,)-' = (r, 4) = (r,-),
for all r E A, 0 E A.

2.3 The Discrete Fourier Transform and Inverse DFT

Definition 1 (Discrete Fourier Transform) If A is a finite abelian group with
card(A) = n, then the Fourier transform of a E L(A) is the function a E L(A)
defined by
(4) = a(r)(r), (2.1)
rEA
for all ) E A. We shall write TA for the map

TA : L(A) ---+ L(A)

such that A (a) = i.

The complex numbers a(O) are known as the Fourier coefficients of a E L(A).
Computation of the DFT of a E L(A) via equation 2.1 requires no more than n2
complex multiplications (not counting multiplication by the normalization constant).
The factors (r, 0) can be kept in a lookup table. Hence, their computation does not
contribute to the complexity of equation 2.1.

Example. For a E L(A) and s E A, the shift of a by s is the function a. E L(A)
defined by a.(r) = a(r s), for all r E A. Applying the Fourier transform to a. yields
a,(0) = (s, 4) a( () = 0(s) fa(O), for all q E A.

One may easily verify linearity of the Fourier transform, which we state for-
mally.

Proposition 5 If A is a finite abelian group with card(A) = n, then FA is a linear
transform, that is,

.A(a + b) = YA(a) + YA(b), for all a,b E L(A),

and
YA(a) = YA(a), for all a E L(A),C E C.

For an application of the Fourier transform, we determine a basis of L2(A)
by computing the Fourier transforms of the basis elements in E = {er : r E A}. Let
er = A(er), for all r E A.

Lemma 6 If r E A, then
1
(or al ( A.r)
for all 0 E A.

Proof. If r E A then

= -=Eer(s)(S,( )
sEA
1

= 1v

for all E A.

Noteworthy is the function eo, where 0 is the zero element of A. Since

eo() = 0(0) = 1

for all 0 E A, -o is constant on A.

Lemma 7 If r E A, then

Proof. If r = 0, then (r,
from the cardinality of A.

is, (r, 4) # 0. Thus,

Sn, ifr=0
(r ) = 0, otherwise.

0) = (0, 0) = 0(0) = 1, for all 4 E A. The result follows

If r # 0, then there exists 0 E A such that O(r) # 1, that

,EA

OEA

OEA

which implies

D(,) =0.
'kEA

Proposition 8 The set 6 = {fr : r E A} is an orthonormal basis of L2(A).

Proof. Suppose r, s E A.

on A x A, we have

Using the inner product in L2(A) and the bilinear form

(Fr,es) = Fr(;b)?s()
EA
1- E 4>(r)4O(s) (by Lemma 6)
n -
4'EA

E !(rs0) (s
n -
_-
EA

n
OA

n
O6A
1, ifr=s
S 0, otherwise,

the last equality by Lemma 7. Thus, (e,, -',) = Sr,, for all r, s E A, which implies &

is an orthonormal basis of L2(A). 0

Since the cardinality of E is n, L(A) has dimension n. Furthermore, since TA

maps an orthonormal basis to an orthonormal basis, we have the following result.

Theorem 9 The Fourier transform FA : L2(A) -+ L2(A) is an isometry, that is,

IIYA(a)II = \all, for all a E L2(A).

Before defining the inverse Fourier transform, we present a dual to Lemma 7.

Lemma 10 If q5 E A, then

D r ) n, if b= 0
EA' f 0= otherwise,

where O is the trivial homomorphism.

15

Proof. If 4 = o, then (r, ) = 0o(r) = 1, for all r E A. The result follows from the

cardinality of A. If 0 # 0o, then O(s) 5 1, for some s E A, and

rEA rEA

rEA
= (s,4))Z(r,4)
reA

implies
D(r,4) = 0.
rEA
0

The next example yields another basis of L2(A).

Example. If 0, l E A, then 0, E L(A). We compute the inner product of and

4':

( 4 = S (r)O(r)
rEA

= D(r, 0)(r, V)
rEA
= 5(r,4)(r,-4')
rEA

= 5(r,4-4)
rEA
n, if=
= 0, otherwise,

the last equality by Lemma 10. Normalization yields an orthonormal basis

e = A = -
fVn-O 0 v}= n

of L 2(A).

Example. If 4 E A, then 4 E L(A). What are the Fourier coefficients of FA= = (O)?
If E C A, then

=--E )
V rEA
1
rEA
1+
VnrEA
Vfn, ifb=-
= 0, otherwise,

the last equality by Lemma 10.

Definition 2 (Inverse Discrete Fourier Transform) If A is a finite abelian group
with card(A) = n, A is the algebraic dual of A, and a E L(A), then the inverse Fourier
transform of a is the function a = )1(ad) E L(A) given by
1
a(r) = 4)(r,), (2.2)
OEA
for all r E A.

Computation of the inverse DFT of a E L(A) via equation 2.2 requires at
most n2 complex multiplications.
The composition of FA and F;'1 is the identity transformation on L(A), as
stated in the next result.

Theorem 11 (Inversion Theorem) If a E L(A) and a = FA(a), then a = F'(-a).

Proof. Let a E L(A) and b = J-^(ia). If s E A, then

b(s) = ()s,)
sEA

V- (-rEA

OEA r A

OEA tEA

= -a(s) n (by Lemma 7)
G~A rEA

Z:a (r) E(r- s

n
rE-A (

a(s).-n (by Lemma?7)

= a(s),

which implies b = a. E

2.4 The Convolution Theorem

Definition 3 (Convolution) If a, b E L(A), then the (cyclic) convolution c E L(A)

of a and b, written c = a b, is defined by

c(t) = a(r)b(t- r),
rEA

for all t E A.

Example. We compute the convolution of a basis element e, E C with a function

a E L(A) at t E A:

(e, a)(t) = E e,(r)a(t r) = a(t s) = a,(t),
rEA

where a, is the shift of a by s A. Thus, e, a = a,.

Example. If , 0 E A, then , 0 E L(A) and

(0 )(0)

= 0(r)7(0 r)
rEA

rEA

rEA

rEAi
f n, if -0=0

0, otherwise,

the last equality by Lemma 10.

The next example shows that we may regard the Fourier transform as a con-

volution.

Example. If a E L(A) and 4 E A, then the properties of the bilinear form on A x A

imply

(a (-0))(0)

= a(r)(-0(0 r))
rEA
E a(r)(-r,-*)
rGA

rEA
= /n. -(<-).

Thus,

a40) =

(a (-k))(0).

The Fourier transform relates convolution and pointwise multiplication. This

relation is known as the convolution theorem.

Theorem 12 (Convolution Theorem) If A is an abelian group with card(A) =n,
a, b E L2(A) and c = a* b, then F= V/ni-i a. b.

Proof. If a, b, c E L(A), with c = a b, and E A, then

(a.6)() = b(1.(

( 1 a(r) (r ) ( Zbss)(
i/n rA ) ( / E
rEA SsEAA

Za~r) (Zb(s)(sc~) (r,4,

n rEA sEA

n
rEA sEA
= 1TZ a(r)b(s)(s + r,,),
n
rEA sEA
using bilinearity. Letting t = s + r and interchanging the summations,

(a. b)() = a(r)b(t r)(t, 0)
rEA tEA

= (ca(r)b(t,-r)) (t,)
n
tEA

tEA
1 -
v"i
V- ,/ CIO).

0

Example. Since a, = e, a from a previous example, the Fourier transform of a, is

a, = v/n a, by the convolution theorem.

2.5 The Cooley-Tukey Factorization

The Cooley-Tukey (CT) factorization is an algorithm that computes the dis-
crete Fourier transform of a function a E L(A), where A is a finite abelian group. An
immediate corollary computes the inverse DFT of a function a E L(A).
For the spatial domain, consider the short exact sequence

0 --- B --- A --- AIB ---- U,

where A is a finite abelian group and B is a subgroup of A. Note that B does not
necessarily split in A, that is, it is not generally true that there exists a subgroup C
of A such that A = B E C (internal direct sum). However, if R C A is a complete
and nonredundant set of coset representatives of the quotient group A/B such that
B n R = {0}, then A = B + R, and for every r E A there exist unique r3 E B and
r, E R such that r = r# + r. We refer to such a set as a complete system of coset
representatives (or residues) of A/B.
Let B. = { E A : (r,>) = 1, for all r E B}. It is easy to see that B. is a
subgroup of A. In the frequency domain consider the short exact sequence

Oo -- B. -- A -* A/B. -- o

where Oo is the trivial homomorphism in A. If R. C A is a complete system of
coset representatives of A/B., then A = B. + R. and for every 0 E A there exist
unique E5 C B. and 0# E R. such that 4 = (f + 0,. One can show that A/B. B
and B. A/B. Let n = card(A), np = card(B), and n, = card(A/B). We have
n = nnp = card(A), card(R) = n, = card(B.), and card(A/B.) = np = card(R.).
An et al. [1], give a detailed discussion of the next result, although their approach
depends on an isomorphism from A onto A.

Theorem 13 (Cooley-Tukey Factorization) Let A be a finite abelian group, B
a subgroup of A, and R C A a complete system of coset representatives of A/B. If
a E L(A) and a = FA(a), then

a()= ,),a(r)(r# O) *(r,,)) (r ) (2.3)
'rBER
for all 4 = 0. + # c- B. + R. = A, where R. C A is a complete system of coset
representatives of A/B., and r = rp + r. E B + R = A.

Proof. Let 4> = 0,, + o E A with 4O E B. and O E R.. If rp E B, then (r, 4o) = 1
by the definition of B.. Hence,

(r,4) = (r +r,. + Op)
= (r/3, o)(r/3, )/)(ro,(fa)(ro,(f>0)
= (r/5, )(ro, )(ro, a),

for all r = rp + ra E A and 4 = + Op E A. The Fourier transform

a() = a(r)(r,'), (2.4)
V rEA
where n = nnn is the cardinality of A, can therefore be written as in the statement
of the theorem (equation 2.3). 0

Corollary 14 Let A be a finite abelian group, B a subgroup of A, and R. C A a
complete system of coset representatives of A/B.. If a E L(A) and a = ,'J(a), then

a^ 1 =O 1 f, ) ).)
v OPAvER 7 76 EB
for all r = ro + r. E B + R = A, where R C A is a complete system of coset
representatives of A/B, and 4 = 0, + 4p E B. + R. = A.

The numbers (r,, 0) and (r,, ) in Theorem 13 and Corollary 14, respec-
tively, are known as twiddle factors. The Cooley-Tukey factorization is a three-stage
algorithm. The first stage is the inner summation in equation 2.3 which requires no
more than npn,, = n complex multiplications for each Op, for a total of n np complex
multiplications. The second stage is an entry-by-entry multiplication by the twiddle
factors consisting of na"n complex multiplications. The last stage is the outer sum-
mation in equation 2.3 which requires no more than n,"np complex multiplications for
each 0,, for a total of n na complex multiplications. The total number of complex
multiplications is bounded by n n + n n, a + n 6n. Consider the computational
savings for the case n = 104, n = 100 = n,. Direct computation using equation 2.4
involves n2 = 108 complex multiplications. The Cooley-Tukey factorization algorithm
requires no more that 2 106 + 104 complex multiplications.

CHAPTER 3
QUOTIENT LATTICES

3.1 Introduction

In this chapter we consider the discrete Fourier transform for functions defined

on finite abelian groups that arise as quotients of D-dimensional lattices. Instead of

working with algebraic duals in the frequency domain, we will take a more intuitive

approach and work with the equivalent geometric duals. This approach differs from

standard presentations of the discrete Fourier transform that emphasize the algebraic

structure of the underlying indexing sets for function domains; our emphasis is on
the geometric structure of the underlying indexing sets or coset representatives.

Section 3.2 presents the discrete Fourier transform and the inverse DFT for
quotient lattices. Since quotient lattices have the structure of finite abelian groups,

all the results of Chapter 2 are immediately applicable. In particular, the inversion

and convolution theorems apply to the quotients under consideration in this chapter.

Section 3.3 states the Cooley-Tukey factorization for quotient lattices in the spatial

and frequency domains.

Denote the field of real numbers by R, D-dimensional Euclidean space by RD,

and the set of integers by Z. Denote the inner product of x, y E RD by (x, y).

3.2 The DFT and Inverse DFT for Quotient Lattices

A D-dimensional lattice is the set of all linear integer combinations of D

linearly independent vectors in RD. Let L be a D-dimensional lattice. The geometric

dual (or reciprocal lattice) L* of L is defined by

L = {s E RD: (r,s) E Z, for all r E L}.

Note that L* is also a D-dimensional lattice. Denote a sublattice relationship by
L1 -< Lo for two D-dimensional lattices L0 and L1. Note that "-<" is also an abelian
subgroup relationship.

Proposition 15 (Lattice Duality Relationship) If Lo and L1 are D-dimensional
latttices, then
L, -< Lo ==> L* -< L\

Remark 1 From a signal processing perspective, we regard Lo in Proposition 15 as

a spatial domain sampling lattice and L, as a spatial domain periodicity lattice.

Similarly, we regard L\ as a frequency domain sampling lattice and LQ as a frequency
domain periodicity lattice.

The sublattice relationship also applies to sublattices of quotient lattices of
the form Lo/Li. Quotient sublattices arise from sequences of the form L, -< L2 -< L0.

Corollary 16 If Lo, L1, and L2 are D-dimensional latttices such that L, -< L2 -< Lo,
then

1. L* -< L -< L\;

2. L2/L1 -< Lo/L, and L2/L; -< L*/Lo;

3. (Lo/L1)/(L2/Lj) Lo/L2 and (L/L*)/(L*/L ) LI/L;.

Define quotient lattices for the spatial and frequency domains by A = Lo/LI
and A' = L /L respectively. If n is the index of L1 in Lo, then A is a finite abelian

group with cardinality n, and nr E L1,for all r E L0. By Proposition 18 (below), the
index of L; in L\ is also n. Hence, ns E L, for all s E L\. The spatial and frequency
domain short exact sequences corresponding to A and A' are

0 ---+ L, -- Lo -- A --0

and
0 -*-+ L L A' 0,

respectively.
Recall that U(n) denotes the multiplicative group of all the nth roots of unity.
The next result defines a bilinear form on A x A'.

Proposition 17 Let Lo and L1 be D-dimensional lattices with L, -< Lo. If A =
Lo/LI1, A' = L*/Lo, and n = card(A), then the map

(, A x A' -+ U(n)

defined by

(:) = exp[-27ri(r,s)],

where i = V/2- and F and represent the costs containing r and s, respectively, is
a bilinear form.

Proof. Since A = Lo/L1 and nr E L1, for all r E Lo, we have (nr, s) E Z, for all
s L*, that is, (r,s) E IZ, for all r E L0 and s E L which implies (F,) =
(r + L1,s + L ) = exp[-27ri(r,s)] E U(n). To show that (r,s) is independent of the
choice of coset representatives, let rl, r2 E L0 with ri = r2 mod L1. Then r, r2 = u,
for some u E L1. Now, if s E L*, then

exp[-27ri(ri,s)] = exp[-27ri(r2, s)] 4== exp[-27ri(r1 r2,s)] = 1

4= exp[-27ri(u,s)] = 1,

which is true for all s E L\ since u E L1. A similar argument shows that (F,3) is
independent of the choice of coset representative of the equivalence class X. Thus,
(, .) is well defined.
The identities

(OF + 'F2, ) = (TI,5)(F2, ), for all F ,T2 E A, j E A',

(F, 3 +s2) = (F, 1)(F,52), for all T E A, 31,32 E A',

and
(cj) = (r,3)c = ( ,cs), for all c E Z, r E A, s E A'

are all easily verified. ]

Note that (r,s) = 1 if and only if (r, s) E Z, where r is any coset representative
of F E A and s is any coset representative of 5 E A'.
The quotient lattice A' is the geometric analogue of the dual group A =
Hom(A, U(n)). Two abelian group isomorphisms follow.

Proposition 18 If Lo and L, are D-dimensional lattices such that L, -< Lo, A =
Lo/L1, A' = L*/L;, and n = card(A) (as in Proposition 17), then (1) A' A, and
(2) A' ^ A.

Proof. Part (1). For every = s + L; E A', define a map < : A -+ U(n), where
n = card(A), by 0,(F) = (F, ), for all F = r + Li E A. The linearity of (-,3) implies
A,, that is, 4> is a homomorphism. Define 'P : A' -4 A by 'P(3) = -^, for all
E A'. The linearity of (F, .) implies %F is a homomorphism. Since card(A) = n =
card(A) (Theorem 2), we need only show T is injective. If i,S2 E A', then

T(1) = 'I(52) == (T,;i) = (f,S2), for all F E A

==# (F(,31 52) = 1

=== exp[-27ri(r,si S2)] = 1, for all r E Lo,

= S1 = 52

== 5I =S2-

Hence, T is an isomorphism.
Part (2). By Theorem 2, A A; by the first part of the proposition, A' A.
Hence, A' A. 0

Denote by L(A) the n-dimensional space of all complex-valued functions on
A and by L(A') the n-dimensional space of all complex-valued functions on A'. The
geometric analogue of the Fourier transform (Definition 1) of a function a E L(A) is
the function a' = .A(a), where FA : L(A) -+ L(A'), defined by

a,(P) = -= a(f)(T,-) = -= a(f)exp[-27ri(r,s)], for all 5 6 A'. (3.1)
v EA V 7EA
The inverse Fourier transform of a function a' E L(A') (the geometric analogue of
Definition 2) is the function a = Jc'(a') E L(A) defined by

a(T) = a'(3)(T, ) = n a'(3) exp[27ri(r,s)], for all F E A. (3.2)
E A' 3EA'
In the sequel we will work with sets of coset representatives of A whenever A
is a quotient lattice, that is, whenever A = Lo/L1 for D-dimensional lattices L0 and
LI. Similarly for the dual quotient A' = LI/L In this context, let G C L0 be a
complete system of residues for A and let G' C L\ be a complete system of residues
for A'. Assume 0 E G and 0 E G', where 0 is the zero vector in RD. Define a map

(, ) on Lo x L; by

(r,s) = exp[-27ri(r,s)], for all r E Lo, s E L*. (3.3)

Since Lo and L1 are Z-modules, and from the proof of Proposition 17, we have

(., .): Lo x L* U(n),

where n = card(A), is a bilinear form. Now, define a: G -4 C by a(r) = a(f), for all
r E G, and a': G' -+ C by a'(s) = a'(3), for all s E G'. Note that G x G' C Lo x L,.
In view of Proposition 17, compute the DFT of a E L(A) as
1 1
a'(s) = E a(r)(r, s) = E a(r)exp[-2ri(r, s)], for all s E G', (3.4)
SrEG V rEG
and the inverse DFT of a' E L(A') as

a(r) = 1: > a'(s)(r,s) = -= j a'(s)exp[27ri(r,s)], for all r E G. (3.5)
BEG' sEG'

3.3 The Coolev-Tukey Factorization for Quotient Lattices

The purpose of this section is to state the Cooley-Tukey (CT) factorization
(Theorem 13) and the inverse CT factorization (Corollary 14) for the pair of sequences
{ B -+ A -+ A/B, B' -+ A' -* A'/B'}, respectively, where A and A' arise as quotients
of D-dimensional lattices, B is a subgroup of A, and B' is a subgroup of A'. If
n = card(B) and n, = card(A/B), then n = nono = card(A).
In the spatial domain, let L0, L1, and L2 be D-dimensional lattices such that
L, -< L2 -< L0. Consider the short exact sequence

0 --* B ---+ A --+ A/B ---0,

where A = Lo/L1 and B = L2/Li -< A. If G C Lo is a complete system of residues for
A, Kp C L2 a complete system of residues for B, and Gc, C Lo such that G = KO+G,

and KfnGl = {0}, then for every r E G there exist unique rj3 E KG and ro E G, such
that r = ro + r,. Furthermore, card(G) = n, card(Kp) = n6, and card(Ga) = n,.

Example. Let L0 = Z, L2 = 10Z, and L1 = 1000Z. Then A = Z/1000Z and
B = 10Z/1000Z. The set G = {0, 1,..., 999} C Z is a complete system of residues
for A. Furthermore, if K# = {0,10,20,... ,990} C 10Z and G, = {0,1,...,9}, then
K# is a complete system of residues for B, Kp n G, = {0} and G = K# + Ga,.

In the frequency domain, we have Lo -< L2 -< LI (Corollary 16). Consider the
short exact sequence
0 --+ B' --- A' ---+ A'/B' -- ,

where A' = LI/L* and B' = L'/L -< A'. Note that B' # L/L*. (To understand the
relationship between B and B', let F = r + L1 E B = L2/L1 and = s + L; E B' =
L2/L;. Since r E L2 and s E L, we have (r,s) E Z. Hence, (r,s) = exp[-27ri(r,s)] =
1. Thus, T E B and E B' implies (r,s) = 1.) If G' C LI is a complete system of
residues for A', K" C L a complete system of residues for B', and G9 C L, such
that K"f n G6 = {0} and G' = K" + G', then for every s C G' there exist unique
sa E K' and sp E GO such that s = sa + sp. Corollary 16 and Proposition 18 imply
B' A/B and A'/B' B. Since A' A, we have card(G') = n,card(K') = n,
and card(GC) = n = card(A'/B').

Example. If Lo = Z, L2 = 10Z, and L1 = 1000Z, then L; = Z,L = 0.1Z,
and L\ = 0.001Z. Furthermore, A' = 0.001Z/Z and B' = 0.1Z/Z. The set G' =
{0,0.001,0.002,... ,0.999} C 0.001Z is a complete system of residues for A' and
K' = {0,0.1,0.2,... ,0.9} C 0.1Z is a complete system of residues for B'. If G#
= {0,0.001,0.002,..., 0.099} C 0.001Z, then K", Gn = {0} and G' = K" + G'p.

If a E L(A) and a' = 17A(a) E L(A'), then the spatial domain CT factoriza-
tion for quotient lattices is determined by

a'() = ( a a(r)(r, sf) (r.,so)) (r.,s.), (3.6)
raEG .( rpEKe

for all s = a+s + s E K + G = G', where r = r# +r, E K +G,, = G. If a' E L(A')
and a = .F'(a') E L(A), then the inverse or frequency domain CT factorization for
quotient lattices is determined by

a(r)- l 1 E a'(s)(r, sa)) (r.,so)) (rs), (3.7)

for all r = ro + r,, E Kc + G,, = G, where s = so, + so E K + G' = G'.

CHAPTER 4
PERMUTOHEDRAL AGGREGATES

4.1 Introduction

The purpose of this chapter is to present a brief introduction to certain hier-

archical arrangements of points (or cells) known as permutohedral aggregates. These

aggregates and their geometric duals are the objects of interest in this dissertation.

This chapter applies the material from Chapter 3 to quotient lattices that are as-

sociated with permutohedral aggregates. Specifically, we describe the Cooley-Tukey

factorization algorithms in the spatial and frequency domains, respectively, for these

quotients.

Section 4.2 gives a geometric description of permutohedral aggregates. In

dimensions two and three they are hexagonal and truncated octahedral aggregates,

respectively. See Conway and Sloane [4] for the necessary background. Section 4.3

discusses quotient lattices of the form A\ = L/N"L, where L is a D-dimensional

lattice that has a ring structure, N is an endomorphism of L, A is a positive integer,

and N'L = NA(L). A A-level permutohedral aggregate is a particular set of coset

representatives of the quotient AA. These representatives are discussed in greater

detail in Section 4.4. Section 4.5 contains statements of the Cooley-Tukey (CT)

factorization and inverse CT factorization in terms of coset representatives of A\ and

the geometric duals A'. Although the quotient lattices A,\ have a ring structure, the

factorization algorithms discussed here do not depend on the multiplicative structure
of these quotients.

4.2 Geometric Description

Permutohedral aggregates can be defined constructively starting with a poly-
tope, known as a permutohedron, associated with the classical D-dimensional root
lattice
LD = {(X0,X1,...,XD) E ZD+l1 :X o + *" +zrD = 0},

which lies in the hyperplane E xi = 0 in RD+l. A D-dimensional permutohedron is
the Voronoi cell of the dual lattice LB; its vertices are the (D + 1)! points obtained
by permuting the coordinates of

(G(-D), |(-D + 2), (-D + 4),..., (D 2), tD).

Permutohedra tile D-dimensional Euclidean space.
In dimension one, L, and L\ can be represented by Z in R. The Voronoi
cells are integer translates of the interval [-1/2,1/2]. A level 1 aggregate consists of
the points G = {-1,0,1}. For A E Z+, a A-level aggregate is the set of 3x integers
G\ A ^=0 3kG. We also refer to the set of Voronoi cells associated with GA as a
A-level aggregate.
In dimension two, L2 and L2 are hexagonal lattices. The hexagonal lattice is
an optimal packing strategy for disks in the plane. The Voronoi cells of a hexago-
nal lattice are regular hexagons. A level 1 hexagonal aggregate consists of a central
hexagon together with six neighboring hexagonal tiles as shown in Figure 4.1. A
level 2 hexagonal aggregate consists of seven neighboring level 1 aggregates. In this
manner, higher-level aggregates are constructed from lower-level aggregates. The hi-
erarchical structure of a level 3 hexagonal aggregate consisting of 73 = 343 hexagonal
cells is shown in Figure 4.2.
In dimension three, L3 is a face-centered cubic (fcc) lattice and LZ is a body-
centered cubic (bcc) lattice. The fcc lattice is an optimal packing strategy for spheres

Figure 4.1: Level 1 hexagonal aggregate.

Figure 4.2: Hierarchical structure of a level 3 hexagonal aggregate.

of equal radius in three-space. The Voronoi cells of an fcc lattice are rhombic do-
decahedra (Figure 4.3). The Voronoi cells of a bcc lattice are truncated octahedra

Figure 4.3: Rhombic dodecahedron.

(Figure 4.4). The bcc lattice can also be defined as the set of vectors (x,y,z) E R3
where the coordinates x, y, and z are all even or all odd integers. Figure 4.5 shows a

level 1 truncated octahedral aggregate, which consists of a central truncated octahe-
dron and 14 neighboring cells. Incidentally, beehive chambers are hexagonal prisms
capped by rhombic dodecahedral ends (see Kaproff [12]).
In the next section we give an alternate description of permutohedral aggre-
gates.

4.3 Quotient Lattices Associated with Permutohedral Aggregates

Let p represent the polynomial p(x) = 1 + x + ... + xD, where D E Z+,
and let A = Z[x]/(p), a commutative ring with identity T = 1 + (p). Let = =

Figure 4.4: Truncated octahedron.

x + (p) and let r = 2 C. The quotient ring A is a free abelian group of rank D
with ordered basis X = (1, ,... ,i"D-), where 1 is the identity element. We have:
D = -(1 + + ... + CD-1) and CD+1 = 1. The quotient ring A can be realized as a

D-dimensional lattice in RD by mapping the basis elements to any set of D linearly
independent vectors in RD. We will simply identify the basis elements in X with a
chosen set of D linearly independent vectors in RD and regard the elements of X as
labels or names for the assigned vectors. The element r determines an endomorphism
N, on A defined by N, : p + rp, for all p E A. For example, if D = 2, then
N,r(C) = r7 = (2 ) = 2-_ C2 = 2C + 1 + = 1 + 3. Let NA = Nr(A). Note that
N7A = rA = (r) is a principal ideal of A. Furthermore, rX = (r, T, ..., rTD-l) is a

Figure 4.5: Level 1 truncated octahedral aggregate.

basis of rA. With respect to the basis X,
S2 00 ... 01
-1 2 0 ... 0 1
0 -1 2 ... 0 1
Nr= : : (4.1)
0 0 0 .** 21
0 0 0 ... -1 3
Furthermore, det N, = 2D+1 1 (see Vince [22]). More generally, for A E Z+,
NAA = NA(A) = TrA = (r7), is an ideal of A, and det N,' = n. Kitto, Vince, and
Wilson [13] prove the following isomorphism of rings.

Theorem 19 If D, A E Z+, then A/NAA A Z,, (ring isomorphism), where n =
2D+1 1. In particular, A/NAA is cyclic and card(A/NAA) = nA.

Let {uj : j = 0,1,..., D 1} be a set of D linearly independent vectors in RD
and denote by U the ordered basis (uo, U1,... u D-1) of RD. Identify Ci with uj, for
j = 0,1,..., D 1, and denote by L the D-dimensional lattice generated by U. Let
N be the endomorphism of L such that the matrix of N with respect to the basis U is
given by the matrix in equation 4.1. Thus, NL = N(L) is a D-dimensional sublattice
of L and the index of NL in L is n = det N = 2D+1 1. More generally, if A E Z+
and NAL = NA(L), then NAL is a D-dimensional sublattice of L that is generated
by the set of vectors {NA(uj) : j = 0,1,..., D 1}. Furthermore, the index of N'L
in L is det NA = (det N) = n.
We turn our attention to quotient lattices derived from L and the sublattices
NAL. Let A = L/NL, that is, A is the cokernel of N E End(L), the set of endo-
morphisms of L. We have card(A) = n = 2D+1 1. More generally, for A E Z+, let
A. = L/N'L, that is, A, is the cokernel of N' C End(L). We have card(A,) = nA.

Remark 2 A A-level permutohedral aggregate (in the spatial domain) is a particular
set of coset representatives of A,\ as a quotient lattice. For permutohedral aggregates,

the uj E U are chosen such that the set of D+ 1 points {u, :j = 0,1,...,D}, where
UD = ]jDo1 u3, is the set of vertices of a regular D-simplex with barycenter at the
origin. See Lucas [14] for a particular choice of coordinates for these vectors.

By the lattice duality relationship (Proposition 15), NL -< L implies L* -<
(NL)*. Hence, the geometric counterpart of A in the frequency domain is the quotient
lattice A' = (NL)*/L*. Similarly, for A E Z+, the geometric counterpart of AA in the
frequency domain is A\ = (NAL)*/L*. We can regard the (injective) endomorphisms
NA as L-invariant, invertible, linear operators on RD (and therefore bounded). Recall
that every bounded linear operator on a Hilbert space has a unique adjoint.

Theorem 20 If H is a Hilbert space and M is a bounded linear operator on H,
then there exists a unique linear operator M*, known as the adjoint of M, satisfying
(M(x),y) = (x, M*(y)), for all x,y E H.

For the elementary results from linear operator theory in this section and for a
proof of Theorem 20, which uses the Riesz representation theorem for a Hilbert space,
see DeVito [6]. We shall use the following results for a linear operator M: if M is
an invertible linear operator and M* is the adjoint of M, then M* is invertible and
(M*)-1 = (M-1)*; if M is given by a real-valued square matrix with respect to a
particular basis, then M* = Mt, the transpose matrix with respect to the same basis;
if A E Z+, then (MA)* = (M*)A; finally, det M* = det M. The next proposition will
allow us to express A' and A' in terms of N*, the unique adjoint of the L-invariant,
linear operator N. In general, for I a D-dimensional lattice and M an I-invariant
linear operator on RD, let MI = M(I). Note that if M is invertible, then MI is a
D-dimensional lattice and MI -< I (sublattice and subgroup relationship).

Proposition 21 If I is a D-dimensional lattice, I* is the dual lattice of I, M is an
I-invariant, invertible linear operator on RD, and M* is the unique adjoint of M,
then (1) (MI)* = (M*)-1I*, and (2) if I' = (MI)*, then M" is I'-invariant.

Proof. Part (1). We first show (MI)* C (M*)-1I*. Let u E (MI). Then
(MA(r),u) = (r, M*(u)) E Z, for all r e I. Hence, M*(u) E I*, that is, there exists
s E I* such that M*(u) = s, or u = (M")-(s). Thus, u E (M*)-1I* and (MI)* C
(M*)-'I*. To show (M*)-1I* C (MI)*, let s E I*. We must show (M*)-'(s) E
(MI)*, that is, (M(r), (AM)-'(s)) E Z, for all r E I. Let r E I. Then M(r) E MI
and (A4(r),(M*A)-'(s)) = (M(r),(M-)')(s)) = (M-'(MA(r)),s) = (r, s) E Z, since
r E I and s E 1*. Hence, (M )-'1(s) E (MAI) and (M*)-'I* C (MI)*.
Part (2). The hypotheses imply MI -< I. Let I' = (MI)* = (M,)-I*,
using the first part of the proposition. Then I* = M*I'. Now, MI -< I implies
I* -< (MI)* (Proposition 15), that is, M*I' -< I'. Hence, MA E End(F'). 0

For D, A E Z+ fixed, let L' = (NAL)*. Part (1) of the proposition implies
L'\ = ((NA')*)-L* = (N*)-L*; part (2) implies (NA)* = (N*)A is L'-invariant.
Hence, L* = (N*)AL' and A' = L/1(N))'L'. Just as AA is the cokernel of N' E
End(L), A' is the cokernel of (N*)' E End(LZ\). The corresponding short exact
sequences-for A\ and A', respectively-are

0 NL -4 L +L/NL 0,

in the spatial domain, and

0 -4 (N)'L' -+ L L/(N)L', -40,

in the frequency domain.

Remark 3 Note that N* is L')-invariant. To see this, consider the chain of lattices
NAL -< N-IL -< ... -< NL -< L. In particular, N'L -< N-L implies (N-IL)" .<
(NAL)*, or (N*)'-L* -< (N*)-AL", by the first part of Proposition 21. Hence, if
L' = (N*)-AL*, then N'L'\ -< L".

Denote by L(A\) the nA-dimensional space of complex-valued functions on
A\, and denote by L(A') the n'-dimensional space of complex-valued functions on
A\. Furthermore, denote the discrete Fourier transform (equation 3.1) of a function
a E L(AA) by a' = .Fx(a) E L(A'), and denote the inverse DFT (equation 3.2) of a
function a' E L(A') by a = .T'(a') E L(AA,).

4.4 Aggregates as Sets of Coset Representatives

Let I be a D-dimensional lattice in RD and M E End(I) such that I det M\ >
1. A set R C I is called a residue system for M if 0 E R and I = U{r + MI : r E R}
(see Bandt [3]). Thus, R is a set of coset representatives of I/MI and card(R) =
det M. A finite address of a lattice point r E I is a finite sequence rA-l "'riro
(A E Z+) such that rk E R, for all k, and r = E\= M k(rk). Vince [22] shows that
the finite address of a lattice point, if it exists, is unique. A A-level aggregate (not
necessarily a permutohedral aggregate) is the set of lattice points having addresses of
length no greater than A. We also refer to the corresponding set of Voronoi cells as a
A-level aggregate. The following lemma states that A-level aggregates are particular
sets of coset representatives of I/MAI.

Lemma 22 Let D, A E Z+. If I is a D-dimensional lattice, M E End(I), I det MI >
1, and R is a residue system for M, then
A-I
R, = Mk R)
k=O

is a set of coset representatives of IIM/I.

Proof. Since R is a set of coset representatives of I/MI and MI fl R = {0}, we
have I = MI + R. Now,

I = .MI+R

= M(M/I+R))+R

= M21 +M (R) + R

= MAI + MA-l(R)+... +M(R) + R

= M\,I+Rx.

Hence, R\ is a set of coset representatives of IIMAI. 0

For D E Z+, let n = 2D+I 1, p(x) = 1 + x +... + xD, A = Z[x]/(p), and
7 = 2 , where = x + (p), as before. Kitto, Vince, and Wilson [13] prove the next
result (we omit the use of overbars for elements of A/NAA; Nr is the endomorphism
of A defined in Section 4.3).

Lemma 23 Each of the following subsets of A is a set of coset representatives of the
quotient ring A/NTA:

1. R={0,1,...,n-1};

2. S = {bDD +. + bj. + bo: bj E {O,1}, forj = O,1,...,D, not all bI = }.

Note that elements in the second set in Lemma 23 can be written as strings (of
length D+1) of binary digits, or bits. Thus, we can represent s = bDCD+.-. +b1+bo0 E
S by the string (bD ... b1b0). Conversely, if (bD .. b1bo) is the binary representation
of r E R, then r =- s mod NA, where s = bDD +... + bj + bo = (bD... bbio) E S.

Example. If D = 1, then A = Z[x]/(1 + x) is isomorphic to Z, R = {0,1, 2} and
S = {0, 1, }. Note that 2 = mod NAA.

Example. If D = 2, then A = Z[x]/(1 + x + x2), R = {0,1,2,3,4,5,6} and
S = {0,1, + 1, 2, 2 + }. The following equivalences (modulo NIA) hold between

the elements of R and the elements of S:

2 = ,

3 + +1,

4 = C2,
5 = 2+1,

6 +.

Note that the binary representation of 6 E R is (110), which represents the equiva-
lence class 1 2 + 1 + 0. 1. Similarly for the other elements of R.

We now consider sets of coset representatives of the quotients A, and A'. Let
U = (uo, u1i,..., UD-1) be an ordered set of D linearly independent vectors in RD and
let UD = ,1D_0 uj. Under the identifications ++ uj, for j = 0,1,..., D, the set

G D
G= bju : bj EO, 1}, for = 0, L., D, not all b = 1} (4.2)
j=O I
is a residue system for N : L -+ L, where L is the D-dimensional lattice generated
by U and N is the endomorphism of L specified by the matrix in equation 4.1 with
respect to the basis U. The set G is equivalent to the second set in Lemma 23 under
the identifications j ++-* uj. Recall that det N = 2D+1 1. By Lemma 22, the set

G, = Z Nk(G)
k=O

is a set of coset representatives of A\ = L/NAL.
In the frequency domain, let 1, = (ul, u,... i u_) be an ordered set of D

linearly independent vectors in RD that generate the D-dimensional lattice L' =
(NAL)* = (N*)A-L-, and let u' = E. 'o u such that
S1}
Ga = bju1: bj E {0,1}, for j = 01,..., D, not all b6 = 1 (4.3)
Ij=o I
is a residue system for N* : L\ -+* L'. Note that det N* = det N. By Lemma 22, the
set
A-I
G' = (N*)k(G')
k=O
is a set of coset representatives of A\ = L\1(N*)L'\.
Example. Let D = 2 and let L be the 2-dimensional lattice in the spatial domain
generated byU = (uo, ul) and let u2 = -(u0 +ul). We have G = {b2u2 +bu i+bouo:
bj {0,1}, not all bj = 1},
N=( -1 3

with respect to the ordered basis U, and det N = 7 = card(G). If A = 2, then L'
is the lattice (N2L)* = (N*)-2L* in the frequency domain. Suppose L'2 is generated
by U' = (uo,U ) and let u' = -(u' + u'). Then G' = {b2z + blu + bo0uo : bj E

{0,1}, not all bj = 1} is a set of coset representatives of A'2 = L'/(N*)2L'. If
uo = (1,0)t and ul = (-l/2,v/3/2)t, then L is a hexagonal lattice in the plane
(congruent to the lattice L2 in Section 4.2); the associated Voronoi cells are hexagons.
Figure 4.6 depicts a level 1 hexagonal aggregate (the set G) in the spatial domain.
The figure shows the seven lattice points and the corresponding hexagonal tiles. Note
that the ordered pair (uO, ul) defines a right-handed coordinate system in the plane.
Level two hexagonal aggregates in the spatial and frequency domains are shown in
Figures 4.7 and 4.8, respectively. Geometrically, the two domains differ in scale,
orientation, and chirality (or handedness).

Figure 4.6: Level 1 hexagonal aggregate in the spatial domain.

Figure 4.7: Level 2 hexagonal aggregate in the spatial domain.

0.6

Figure 4.8: Level 2 hexagonal aggregate in the frequency domain.

Example. If D = 3, then

( 2 0 1)
N= -1 2 0
0 -1 3
and detN = 15 = card(G). Choose u0 = (1,0,0)t,ui = (-1/3,2v=/3,O),u =
(-1/3,-v/2/3,2/v/6)t. These vectors generate a body-centered cubic lattice (con-
gruent to the lattice L; from Section 4.2). The associated Voronoi cells are truncated
octahedra. A level 1 truncated octahedral aggregate is shown in Figure 4.5.

4.5 The Cooley-Tukey Factorization for Permutohedral Aggregates
Fix D, A E Z+ and let n = 2D+1 1.
In the spatial domain, consider the short exact sequence

O ---+ B -4 A) -+ Q, ---- U,

where AA = L/NAL, BO is the subgroup N'L/NAL of AA, and Qa, = AA,/BO. For
1 < 7 A, let Gy = E"-= Nk(G), where G is a residue system for the endomorphism
N : L -+ L. If Ko = No(G9), then K n G,, = {0}, and G\ = K# + G, is a complete
set of coset representatives of A\. Hence, for every r E G\, there exist unique
ro E K# and r, E G, such that r = r6 + ra,,. We regard KO as the set of high-order
representatives in G\; similarly, we regard Ga as the set of low-order representatives
in G\.
In the frequency domain, let L'\ = (NAL)* and consider the short exact se-
quence
0-- B A',- Q 0,
where A' = L'/(N*)'L', B, is the subgroup (N*)OL'/(N*)'L of A', and Q' =
AI/B^. For 1 < 7-y < A, let G'v = E-'>=(N*)k(G!), where G' is a residue system for

the endomorphism N* : L' -* L'. If K', = (N*)O(G'), then K' n G = {0}, and
G' = K' + Gn is a complete set of coset representatives of A'. Hence, for every
s E G, there exist unique s. E K. and sp E G' such that s = S, + sp. We regard
K' as the set of high-order representatives in GC, and G6 as the set of low-order
representatives in G'.
If a E L(A\), then the discrete Fourier transform a' = F\(a) E L(a) can
be computed by the spatial domain Cooley-Tukey factorization (Theorem 13 and
equation 3.6) as

a' (s) o/1 = ;/ 1 a (r)(r, s,)) (r,-,),) (r,S),
raEGa rpEK /
for all s = s, +sp E K +G, = G' where = r + r, E K +G = G\. If
a' E L(A'), then the inverse DFT a = .T'(a') E L(AA) can be computed by the
frequency domain CT factorization (Corollary 14 and equation 3.7) as
aro) 1 -z(( 1 ) s
,() = *PEG a EK a'(K)(r,,S) s. (,, (, ),

for all r = r6 + ra, E K, + Ga = GA, where s = S, + so E K' + G = G'A.
The ideals Bp do not split in A\. For example, in dimension two, A,\ Z/7AZ,
which does not split since 7 is prime. Thus, the Good-Thomas algorithm (as described
by An et al. [1], for example) is not applicable to discrete functions defined on the
quotients A\. The Good-Thomas algorithm is related to the Chinese remainder
theorem and is used to compute the DFT of functions with rectangular domains, that
is, domains of the form Z,, x*..* x ZL. The algorithm exploits the linear separability
of the Fourier transform by first computing the one-dimensional transform over the
rows of an image and then computing the one-dimensional transform over the columns
of the resulting intermediate image (similarly for higher dimensions).

CHAPTER 5
THE GENERALIZED MATRIX PRODUCT

5.1 Introduction

The goal of this chapter is to introduce the generalized matrix product over

digital indexing sets. The generalized matrix product, formulated by Ritter [17],
generalizes standard matrix multiplication. By using this product, data reordering

stages in the FFT algorithm are eliminated. Data reordering is a feature of stan-

dard FFT algorithms. Zhu and Ritter [26] examined the computational advantages

of using the generalized matrix product in signal processing applications. Digital
indices address the elements of permutohedral aggregates and their geometric duals.
This form of addressing is known as generalized balanced ternary (GBT) and was
developed by Gibson and Lucas [8]. The generalized matrix product described by

Ritter and Zhu is not easily adapted to data arrays indexed by strings of digits, as in

the GBT scheme for permutohedral aggregates. While the generalized matrix prod-

uct and GBT addressing for permutohedral aggregates have been studied previously

and independently, the product defined in this chapter is tailored to processing data

arrays that are associated with permutohedral aggregates and their duals.
Let F be a ring and denote by F,,xm the set of I x m matrices with entries in F.
The definition of standard matrix multiplication, which is a map Fixm x F,,,mxr -+ Fxr,
requires the number of columns in the left factor to equal the number of rows in the
right factor. With p E Z+, the generalized matrix product of order p (or p-product)
is a map p : Fixmp x Fpqxr -+ Fq~mrn. Hence, the p-product A p B requires that

p divide the number of columns in the left factor and the number of rows in the
right factor. Note the dimensions of the resulting matrix. We will also make use

of the Hadamard product, a componentwise operation on matrices having the same

dimensions, defined as follows. For A, B E Flxm, the Hadamard product C = (cij) of
A = (a1j) and B = (bj), denoted by A o B, is given by cj = aijbj, for all indices i
and j. Hence, the Hadamard product is a map o : Fixm,,, x Fixm -+ Fixm.

Section 5.2 defines F-indices, indexing sets consisting of finite strings of digits
that address the elements of permutohedral aggregates and their geometric duals
(that is, coset representatives of A\ and A', respectively). Addressing with F-indices
can be implemented in terms of strings of bits, or binary digits, and yields very

efficient data processing algorithms (see Gibson and Lucas [8], and Kitto, Vince, and
Wilson [13]). Section 5.3 defines the generalized matrix product, or p-product, over
r-indices, and Section 5.4 covers digital addressing for permutohedral aggregates and
their duals.

5.2 Digital Indexing Sets

The elements of D-dimensional permutohedral aggregates and their duals can
be addressed using finite strings of digits. We refer to the underlying digital indexing
sets as F-indices, defined as follows.

Definition 4 (Gamma Indices) Let D E Z+ and n = 2D+1 1.

1. Define the set

together with the order relationship

2. For A E Z+, define the set

FA = {gA-g1A-2 ... gigo0 : gk E Fi, for k = 0,1,..., A -1}

and extend the order relationship on F1 to lexicographic order on r\, where the
rightmost digit go is the low-order, that is, least significant, digit;

3. For A,a,,3 E Z+ such that A = a + 0, define the set

FrF = {gh :g E r,,,3 E Fr}.

We refer to the sets F,\ as F-indices.

Remark 4 The following properties of F-indices follow directly from Definition 4.
Fix D E Z+ and let n = 2D+1 1.

1. The set F\ has nA elements;

2. If A, a, # E Z+ such that A = a +/3, then

F F, = Fr = FrF,.

5.3 The Generalized Matrix Product for Permutohedral Aggregates

The generalized matrix product, or p-product, over F-indices has a particularly
simple form. Fix D E Z+ and let n = 2+1 1. Let F be any ring and a,,3 E Z+.
Denote by Fr. xr, the set of all no x nO, F-valued matrices with rows indexed by the
set of indices F,,, and columns indexed by the set of indices Fp. In this notation,
x E Fr.,, is a column vector of length n' indexed by Fr, and y E Flxr, is a row

vector of length n3 indexed by Fr. If transposition is denoted by superscript t, then
x E Fr.xI implies xt E Fixr0.

Definition 5 (Generalized Matrix Product Over Gamma Indices) Let F be
any ring, D E Z+, and n = 2D+1 1. Let A E Fr xrp, where a,3 E Z+, and

B E Fr,xr,, where -7,5 E Z+-. If A E Z+ such that A < {/,7}, then the p-product of

A and B, with p = nA, is the n'+"- x nO-A+ matrix

C = A p B E Fr.r,_- xrg r,

defined by

C(gt,rm) = A(g,rs)B(st,m),
sErF
for all gt E rr.,-,rm E 1-A s.

Remark 5 For full matrices A and B as in Definition 5, computation of the n\-

product C = A E,, B requires n++"1+6-'- multiplications. The case A = = 7

corresponds to standard matrix multiplication.

Example. If D = 1, then n = 3 and Fi = {0,1,2}. If x E Rixra and y E Rrxi,

then
x = (x(000), x(001), x(002), x(010),..., x(221),x(222))

and

y = (y(00),y(O1),y(02),y(10),..., y(21),y(22))t.

If A = 2 and z = x ,\ y = x 32 y, then z E Rixr, and

z(r) = x(rs)y(s),
sEF'2
for r= 0,1,2.

Example. If D = 2, then n = 7 and F[ = {0,1,...,6}. Suppose A E Cr2, xr and
B E Crxr1. If C= A&7B (A= 1), then CE Crsxr2 and

C(gt,rm) = 1 A(g,rs)B(st,m),
sEri

for all (gt,rm) F2F3 x F2F1 = F5 x F3.

5.4 Digital Addressing for Permutohedral Aggregates

Our goal in this section is to express the elements of a permutohedral aggregate
(and of the corresponding geometric dual) entirely in terms of digits, or F-indices. We
will need to distinguish between digital addresses for elements of permutohedral ag-
gregates in the spatial domain and digital addresses for elements of the corresponding
duals in the frequency domain. To this end, let Fr = Fr and Fr, = Fr (A E Z+).
Fix D, A E Z+. Let L be the D-dimensional lattice generated by U =
(uo, u1,..., UD-1i) and let N be the endomorphism of L specified by the D x D
matrix in equation 4.1 with respect to the basis U. Let U,\ = (u', u',..., un-) be an
ordered set of D linearly independent vectors that generate the lattice L' = (N*)-AL.
Let UD = D_ o uj and u'= -_o u. Recall the sets of coset representatives
G, G', GA, and G', described in Section 4.4. Note that any element of F1, that is,
any digit, can be expressed in the form bD2D + ..-. + b12 + bo, where the bj E {0, 1},
for all j, and not all bj = 1. Furthermore, the bj are unique. Assign the digit

bD2D +. +b12+b0 E F1 to the vector bDUD + +b6u +bouo E G. This assignment
of digits establishes a one-to-one correspondence between Fx and G. Now, assign the
string of digits g\.. *.. gigo E FA to the vector r = N\- (r,\-) +.- + N(ri) +ro C G,\,
where gk E F1 is the digit corresponding to the vector rk E G, for k = 0, 1,..., A 1.
In this manner, we establish a one-to-one correspondence between the digital in-
dexing set F\ and G\ (a permutohedral aggregate in the spatial domain). Write

9k = rk whenever gk E F1 is the digit corresponding to rk E G and refer to the string
9gA\- gi9o = N'-A(gA-I) + + N(gi) + go as the digital address of r.
In the frequency domain, assign the digit bDo2 + ... + b12 + bo E to bDuoD +
S+ blux + bouo E GC'. This assignment of digits defines a one-to-one correspondence

between Fr and G'. Assign the string of digits h,\-\l h1ho E r, to the vector
s = NA-l(s\i) + + N(si) + so E G', where hk E F' is the digit assigned to
sk E G', for k = 0,1,..., A 1. This assignment of digital strings defines a one-
to-one correspondence between the digital indexing set Fr and G' (a permutohedral
aggregate in the frequency domain). Write hk = Sk whenever hk E Fr is the digit
assigned to Sk E G'; the string hA-, ... h1ho = EiI(N*) k (hk) is the digital address
ofs.
Recall from Section 4.4 that the finite address of a lattice point, if it exists,
is unique. Hence, the digital addressing scheme we have described assigns unique
addresses to the elements of permutohedral aggregates in the spatial domain; similarly
for their counterparts in the frequency domain. Also, note that the elements of Li-
or the corresponding vectors in G-can be expressed as bit patterns of the form
(bD ... bibo). Similarly for the elements of r' and the corresponding vectors in G'.

Example. Let D = 2 and A = 3. The digital address of the vector N2(u2 + Uo) +
N(u2)+uo = N2(1.u2+0.u+l1.uo)+N(1.u2+0.ui+O.uo)+(O.u2+0.ui+l1.-uo) =
N2(101) + N(100) + (001) = N2(5) + N(4) + 1 E G3 is 541 E r3. Note, however, that
541 E rF is the digital address of the vector (N')2(ul + Uo) + N*(u') + U'o E G'3 in
the frequency domain.

Example. Figure 5.4 shows an indexed level 3 hexagonal aggregate in the spatial
domain.

Figure 5.1: Indexed level 3 hexagonal aggregate in the spatial domain.

CHAPTER 6
DIGITAL ALGORITHMS FOR PERMUTOHEDRAL AGGREGATES

6.1 Introduction

This chapter brings together the concepts developed in earlier chapters to

give statements of the discrete Fourier transform, the Cooley-Tukey factorization,

and the fast Fourier transform-together with the respective inverse transforms-for
complex-valued functions defined on permutohedral aggregates and their geometric

duals. The algorithms in this chapter are expressed in terms of the generalized matrix

product over digital indexing sets.
The two theorems in Section 6.2 define the discrete Fourier transform (DFT)

and the inverse DFT for permutohedral aggregates. The theorems and lemmas in

Section 6.3 are digital Cooley-Tukey (CT) factorization algorithms for permutohedral

aggregates. Section 6.4 gives statements of the fast Fourier transform (FFT) and the

inverse FFT algorithms for permutohedral aggregates addressed by digital indexing

sets.

Throughout, assume D, A E Z+ fixed and let n = 2D+1 1. Recall that
n) is the cardinality of a D-dimensional, A-level permutohedral aggregate. We will

use the digital indexing sets FA and FP to address spatial and frequency domain

data, respectively, as discussed in Section 5.4. As stated earlier, very efficient data

processing results from the use of digital indexing sets to address permutohedral

aggregates and their geometric duals.

Define the map
G, x G' -- U(nA)

by
(r,s) = exp[-27ri(r, s)],
where G\ is a set of coset representatives of A\ = L/NAL, G' is a set of coset
representatives of A' = L'\/(N*)A'L', and (r,s) is the inner product of r,s E RD
(see Proposition 17); (.,.) is a bilinear form on L x L'. For a E L(A\), define
a : GA -4 C by a(r) = a(f), for all r E GA; for a' E L(A'), define a': G'\ -+ C by
a'(s) = a'(), for all s 6 G'.. If g E rA is the digital address of r E GA, and h E r,
is the digital address of s E G', then define (g, h) = (r, s). Hence, in the context of
digital addressing, (., .) is a map r\ x Fr -4 U(nA).

6.2 The Discrete Fourier Transform and Inverse DFT

Let a,/# E Z+ such that A = a + P3. For the direct transform, define e0 E
Cixr rr (a row vector of length n2") by

e,4(vaq) = (N/(q'),v') E U(nA), q" E r,V r :,.

Let e = el. Note that e\(h\gA) = (g\,hA) and e(hg) = (N-1(g),h). For the inverse
DFT, define e' E Cixr.r; (a row vector of length n2a) by

e'(qvc') = (NA(q'),v')-1 = (NIA(qa),v) U(nA), q rC r'.

Let e' = e'. Note that e'(g\hA) = (gA, hA) and e'(gh) = (N1-(g), h).

Theorem 24 (Digital DFT for Permutohedral Aggregates) If D, A E Z+, n=
2D+1 1, a E L(AA), and a' = .F\(a), then

1
a' = ;-h2 (e\ Gx a)t.

Proof. The discrete Fourier transform of a E L(A)) is determined by

a'(s) E a(r)(r,s), s E G. (6.1)
rEG\
Let g" E r,\ be the digital address of r E G\ and hA E 1F the digital address of
s E G'. Then equation 6.1 is equivalent to

aG'(hA) = n/' a(g\)(g\, h\)
g9'Er
n/ 9A E rA n
1~I E e.\(h'g')a(gA), hA E rj\*

In terms of the nA-product, the last summation can be written

a' = (e\ a)t (6.2)

to yield a column vector of length n. 0

The nA-product in equation 6.2 requires no more than n2A complex multipli-
cations. A similar proof yields the corresponding result for computation of the digital
inverse DFT, which recovers a E L(A\) from a' E L(A'A).

Theorem 25 (Digital Inverse DFT for Permutohedral Aggregates) If D, A E
Z+,n = 2D+1 1, a' E L(A'J), and a = ."(a'), then

a = 1 (e a,)t.

6.3 Cooley-Tukey Factorization Algorithms

Let a, / E Z+ such that A = a + fl. For the Cooley-Tukey (CT) factorization
of the direct transform, define the n" x no matrix Da,0 E Cr. xr, by

DC,(qO',wf) = (qa,wf) E U(n'), q E rw E rf.

For factorization of the inverse transform, define the n6 x n' matrix D' E Cr, xr-
by
D',wq") = (qr,ww) E U(n'), wO 6 ,'( ,q' E F,.

Lemma 26 Fix D,A E Z+ and let a, 3 E Z such that A = a +3.

1. If p E Fr and q E Fr,, then N0(p-) + qa = pqa E r FA;

2. ifv E r' and w1 E r', then (N*)(v') + w, = vau E r,.

Proof. Part (1). If p3 = p-i -pipo and qO = q,_ " qiqo, then

paq = P0-I . pipoq-I ... q qo0
a+13-1 --
= Z k(pk.)+ k(qk)
k=a k--=O

k=O k=0
'8-1 CI-1
= NE Nk (pk)+ E Nk(qk)
k=O k=O
= NA(9)+ q.

Part (2) is proved similarly. 0

Lemma 27 Fix D,A E Z+. If r E L and s E L'\, then (NA(r),s) = 1.

Proof. We have NA(r) E NAL. Since s E L' = (N'L)*,(NA(r),s) E Z, which
implies (N'(r),s) = exp[-27ri(N(r),s)] = 1. 0

Theorem 28 (Digital CT Factorization) If D, A, a, 0 E Z+ with A = a + f3,
n = 2n+l 1, and a E L(AA), then

eA e,\ a = e, en (D0,, o (eO pn a)).

Proof. The nA-product e% E,. a is equivalent to the summation

E e\(hAgA)a(gA), hA E r.
9" Er.
Ifg\ = 9qo' E rFr, = FA and h\ = v=wo E 'FrF = rF, then

e,(hVgA) = (gA, h\)
= (pyqaIvaw)

= (Na(pl) + q+ (N*)~ (va) + w) (by Lemma 26)

= (Na (P),(N*)o(0'))* (N(p,),w'") (qQ,(qN*)fi(v0)) (q-,w)

= (NA(p9),v). (N"(po),w). (No(q'),v0) (qa, wo)
= 1 eo(wlo) e,(vq') D.c,,(q', w),

the last equality by Lemma 27. Hence,

Z eA(hAg)a(g) = e(v"wpq)a(pq) (6.3)
9AEr, qEra p ErI,
= E e,(vaq') [D'.,(q"',wI)Co(q`,w/o)], (6.4)
qa'Era
where
Co(q,uw) = ) ep(wVp)a(p9 q'). (6.5)
pOEFo
Summation 6.4 is equivalent to e. E,,. (D ,, o Co) and equation 6.5 is equivalent to
C# = e# (,o a. Hence,

eA Gn a = @n (Da,, o Cp)

= o e,1 (D,,, o (eo (O a)),

a row vector of length n. 0

The inside no-product in the CT factorization algorithm ea (no (Da,# o ,(eED
a)) requires no more than n2" complex multiplications; the Hadamard product takes

no more than nA complex multiplications. Hence, the total number of complex mul-
tiplications required by the CT factorization algorithm is bounded by n2" + n20 + n'.

Remark 6 Theorems 24 and 28 imply

a' = 1 (e,, ,,. (D(,0 o (ep ,, a)))t.

Arguments similar to those in the proof of Theorem 28 establish the digital
inverse CT factorization algorithm.

Theorem 29 (Digital Inverse CT Factorization) If D, A, a, /3 E Z+ with A =
a + /3, n = 2D+1 1, a' E L(A'), and a = .1(a'), then

el eLA a' = e' n (DI, o (el ena a')).

Remark 7 Theorems 25 and 29 imply
1
a=- (ee ^ (DN,. o (e' (Bn, a')))'.

Theorems 28 and 29 are digital versions of the CT factorization for permutohe-
dral aggregates. However, the next two lemmas, which use variations of the matrices
D,,,3 and D,a,, respectively, yield digital FFT algorithms (direct and inverse) for
permutohedral aggregates. Define T,,+i,-i E Cr+, xr'__ by

T.+,,- 1(pMq',w-1) = (N'(pM),w'01) E U(n'), pMq' E r ++,w' -1 e r,

and define Tp++,_1 Cr.+ xr.-i by

T+1,Q-1_(VMWqa-1) = (q-', (N*)P(vM))-l E U(n'), VMW c F +,q'q -1 E I,-1i.

In the signal processing literature, the entries in matrices DQ,,3,D#,,,T +i,p-, and
T#+, 1-_ are known as twiddle factors.

Lemma 30 Let D,A E Z+, n = 2+1 1, and a E L(AA). If 3 E Z+ with 2 < t < A,
then
eo E)p a = e D (To+1,-1 o (e-_i .-i a)),

where a = A -3.

Proof. The na-product e nep a is equivalent to the summation

E eo(w'p')a(pAq), q' E r,,w E r,.

If P = P-'PM E r-,r, = ro and w3 = WMW3-1 E rr_ = rL, then

ep(wlp) = (Ng(p),w)3)
= (N'('-pM),wMw'3^-1)
= (Na(N(p-1) + PM), (N*)3-(WM) + wf-1)
= (Ng+'(p-') + Na(pM), (N*)-I(wM) + wl-1)
= (Na+l(p-l), (N*)/-l(WM)) (N+'(p9-1), w-1)
(Na(pM),(N*)I-l(wM)). (N2(pM), w-1)
= (NA(p-1),WM). (N 1(p-1),w-1)
(N'-'(pM),wM) (N'(pM),W -1)
= 1- e-i(w-'p-') e(wMpM) o+r ,_-,(pm0, "),

the last equality by Lemma

E ew(w'p3)a(pAq')
p Erp

27. Hence,

= E E e(WMW -^'8-'pM)a(^P-lPMq')
PMErI pAe-1Ero-_
= E [Ti + ,w-)(PMq'w -1l)-l(pMqc,w -1)]
e(wpMprM),
"e(WMPM),

(6.6)

64

where
C#(pMqaw8l)= S e#j(w-lWpa-)a(p-PMq). (6.7)
PO-'Ero-i
Summation 6.6 is equivalent to e e, (Ta+l,p-i o Cp-i) and equation 6.7 is equivalent
to C,-I1 = ej-1 n,-1 a, a matrix of dimensions no x nO. Hence,

ep Enf9a = e (n (Tr+i,p.-i oC#_i)

= e e, (TM+4,,6-1 o (ef_-i (nP-i a)).

Similar arguments establish the corresponding inverse digital algorithm.

Lemma 31 LetD, A E Z+,n = 2D+1-1, anda' E L(A'). If a E Z+ with < a < A,
then
e no a' = e' E (T_+1OI o (e'1_, E-, a')),

where 3 = A a.

6.4 The Fast Fourier Transform and Inverse FFT Algorithms

The following algorithm is a radix-n, decimation-in-space, fast Fourier trans-
form for permutohedral aggregates addressed by digital indices.

Theorem 32 (Digital FFT for Permutohedral Aggregates) Let D, A E Z+, n =
2D+1 l,a E L(AA), and a' = F(a). If

C, = e (. a, (6.8)

and

C'o = e ,, (TA-2+j,.-1 o C0-i), forf/ = 2,3,..., A,

(6.9)

65

then

a=Lc

Proof. We claim Cp = e, &,, a, for 1 hypothesis. Arguing by mathematical induction, if C, = e- en-y a, for 1 < 7 < A 1,
then

C-f+l = e E) (TA-_,o 0 CU,)

= e n (T-y,-y o (e- ,, a))

= e.Y+1 .n-+i a, (by Lemma 30)

for 1 < 7-y < A 1. Hence, CD = e3 nO a, for 1 < /3 < A. In particular, CA = e,\ ED a,
which implies a' = n-C/2 C,\, by Theorem 24. 0

A similar argument proves the corresponding digital inverse FFT for permu-
tohedral aggregates, a radix-n, decimation-in-frequency algorithm.

Theorem 33 (Digital Inverse FFT for Permutohedral Aggregates) Let D, A E
Z+,n = 2D+1 1,a' E L(A'), and a = .7' (a'). If

C = e' ,n a',

and
'= e'. (T o C.), fora = 2,3,..., A,

then
] t
a = (C\) .

Computation of the DFT of a function a E L(A\) consisting of n- lattice
points or cells requires A stages, as specified by equations 6.8 and 6.9. Computation

of the n-product in each of these equations requires no more than nA+1 complex
multiplications. The Hadamard product in each of the A 1 stages specified by
equation 6.9 requires at most n' complex multiplications. Hence, the total number of
complex multiplications required to compute the DFT of a using the FFT algorithm is
bounded by An+l +(A 1)nx. The algorithm, therefore, has complexity O(M log M),

where M = n. Furthermore, since the FFT algorithm (Theorem 32) is written
entirely in terms of the generalized matrix product and Hadamard multiplication, no
explicit reordering of the data is required. Data reordering is a common feature of

standard FFT algorithms. Similar remarks apply to the computation of the digital
inverse FFT algorithm (Theorem 33).

CHAPTER 7
SUMMARY AND SUGGESTIONS FOR FURTHER RESEARCH

In this dissertation we showed how the additive structure of permutohedral
aggregates can be exploited to derive fast algorithms for the computation of the
discrete Fourier transform of functions whose domains are aggregates of any given

level in any given dimension. In these algorithms, the data is addressed in terms of

an efficient digital addressing scheme known as generalized balanced ternary (GBT)

addressing. The algorithms described here are not based on the Chinese remainder
theorem, which is used to derive the fast Fourier transform (FFT) for functions de-
fined on rectangular arrays. Furthermore, we defined the generalized matrix product
for permutohedral aggregates and stated the FFT and the inverse FFT for permuto-
hedral aggregates entirely in terms of digital addressing and the generalized matrix

product. The notation described in the text renders very succint algorithms that are

amenable to computer implementation.
Note, however, that we have not directly exploited the multiplicative structure

of permutohedral aggregates. A number of algorithms exist that compute the discrete

Fourier transform for domains that possess a multiplicative structure (see Tolimieri,
An, and Lu [20]). Hence, we propose the following problems for further research:

1. characterization of ideals and quotients of permutohedral aggregates in terms

2. derivation of the discrete Fourier transform of periodic and decimated functions
defined on permutohedral aggregates;

68

3. derivation of DFT algorithms that exploit the multiplicative structure of per-

mutohedral aggregates;

4. extension of the definition of the generalized matrix product for permutohedral

aggregates to more general addressing schemes;

5. description of self-similar tiles derived from permutohedral aggregates and deriva-

tion of a discrete wavelet transform and a multiresolution analysis for such tiles.

REFERENCES

[1] M. An, I. Gertner, M. Rofheart, and R. Tolimieri, Discrete fast Fourier transform
algorithms: a tutorial survey, Advances in Electronics and Electron Physics, Vol.
80, Academic Press, Inc., New York, NY, 1991, pp. 1-67.
[2] L. Auslander and R. Tolimieri, Is computing with the Fourier transform pure
or applied mathematics?, Bulletin (New Series) of the American Mathematical
Society, Vol. 1, No. 6, November 1979, pp. 847-897.
[3] C. Bandt, Self-similar sets 5. Integer matrices and fractal tilings of R", Pro-
ceedings of the American Mathematical Society, Vol. 112, No. 2, June 1991, pp.
549-562.
[4] J. H. Conway and N. J. A. Sloane, Sphere packing, lattices, and groups,
Springer-Verlag, New York, NY, 1988.
[5] J. W. Cooley and J. W. Tukey, An algorithm for the machine calculation of
complex Fourier series, Math. Comp., Vol. 19, No. 2, 1965, pp. 297-301.
[6] C. L. DeVito, Functional analysis and linear operator theory, Addison-Wesley
Publishing Company, New York, NY, 1990.
[7] D. E. Dudgeon and R. M. Mersereau, Multidimensional digital signal processing,
Prentice Hall, Englewood Cliffs, NJ, 1984.
[8] L. Gibson and D. Lucas, Spatial data processing using generalized balanced
ternary, Proceedings of the IEEE Conference on Computer Vision and Pat-
tern Recognition, June 1982, pp. 566-571.
[9] R. C. Gonzalez and R. E. Woods, Digital image processing, Addison-Wesley,
[10] T. W. Hungerford, Algebra, Springer-Verlag, Berlin, 1980.
[11] M. T. Heideman, D. H. Johnson, and C. S. Burrus, Gauss and the history of the
fast Fourier transform, IEEE ASSP Magazine, October 1984.
[12] J. Kappraf, Connections: the geometric bridge between art and science, McGraw-
Hill, Inc., New York, NY, 1991.
[13] W. Z. Kitto, A. Vince, and D. C. Wilson, An isomorphism theorem between the
p-adic integers and a ring associated with a tiling of N-space by permutohedra,
Discrete Applied Mathematics 52, 1994, pp. 39-51.

[14] D. Lucas, A multiplication in N-space, Proceedings of the American Mathemat-
ical Society 74, 1979, pp. 1-8.
[15] R. M. Mersereau and T. S. Speake, The processing of periodically sampled mul-
tidimensional signals, IEEE Transactions on Acoustics, Speech, and Signal Pro-
cessing, ASSP-31, No. 1, 1983, pp. 188-194.
[16] D. P. Petersen and D. Middleton, Sampling and reconstruction of wave-number
limited functions in n-dimensional euclidean spaces, Information and Control,
Vol. 5, 1962, pp. 279-323.
[17] G. X. Ritter, Heterogeneous matrix products, in Image algebra and morpholog-
ical image processing II, Vol. 1568 of Proceedings of SPIE, San Diego, CA, July
1991, pp. 92-100.
[18] G. X. Ritter and J. N. Wilson, Handbook of computer vision algorithms in image
algebra, CRC Press, Boca Raton, 1996.
[19] G. Strang, Wavelet transforms versus Fourier transforms, Bulletin of the Amer-
ican Mathematical Society, Vol. 28, No. 2, April 1993, pp. 288-305.
[20] R. Tolimieri, M. An., and C. Lu, Algorithms for discrete Fourier transform and
convolution, Springer-Verlag, New York, 1989.
[21] R. Tolimieri, M. An., and C. Lu, Mathematics of multidimensional Fourier
transform algorithms, 2nd ed., Springer-Verlag, New York, 1997.
[22] A. Vince, Replicating tesselations, SIAM Journal of Discrete Mathematics, Vol.
6, No. 3, August 1993, pp. 501-521.
[23] A. B. Watson and A. J. Ahumada, Jr., A hexagonal orthogonal-oriented pyramid
as a model of image representation in the visual cortex, IEEE Transactions on
Biomedical Engineering, Vol. 36, No. 1, 1989, pp. 97-106.
[24] H. Zhu, The generalized matrix product and its applications in signal processing,
Ph.D. Thesis, 1993.
[25] H. Zhu and G. X. Ritter, The generalized matrix product and the wavelet trans-
form, Journal of Mathematical Imaging and Vision, Vol. 3, No. 1, 1993, pp.
95-104.
[26] H. Zhu and G. X. Ritter, The p-product and its applications in signal processing,
SIAM Journal of Matrix Analysis and Applications, Vol. 16, No. 2, 1995, pp.
579-601.

BIOGRAPHICAL SKETCH

Jaime L. Zapata was born in Bogota, Colombia. He received the Bachelor
of Science degree in engineering with High Honors from the University of Florida in
August 1985, and the Master of Science degree in biomedical engineering from the

University of Miami in August 1986. From 1986 to 1990, he worked as a software
engineer with Coulter Electronics, a major biomedical engineering firm in Hialeah,

Florida. His research interests include image algebra, image processing, wavelets and
multiresolution analysis, and bioinformatics and biocomputing.

I certify that I have read this study and that in my opinion it cnforms to
acceptable standards of scholarly presentation and is fully adequate/ i scope and
quality, as a dissertation for the degree of Doctor of Philosophy.

Gerwrd X. Ritter, Chairman
P essor of Mathematics

I certify that I have read this study and that in my opinion it conforms to
acceptable standards of scholarly presentation and is fully adequate, in scope and
quality, as a dissertation for the degree of Doctor of Philosophy.
694s-e. L&P
David C. Wilson
Professor of Mathematics

I certify that I have read this study and that in my opinion it conforms to
acceptable standards of scholarly presentation and is fully adequate, in scope and
quality, as a dissertation for the degree of Doctor of Philosophy.

Murali Rao
Professor of Mathematics

I certify that I have read this study and that in my opini n it conforms to
acceptable standards of scholarly presentation and is fully adequ te, in scope and
quality, as a dissertation for the degree of actor of Philosophy. ]

'Ame sE. Keesling f
/4 or of Mathematics

I certify that I have read this study and that in my opinion it conforms to
acceptable standards of scholarly presentation and is fully adequate, in scope and
quality, as a dissertation for the degree of Doctor of Philosophy.

J ph N. Wilson
Assistant Professor of Computer and
Information Science and Engineering

This dissertation was submitted to the Graduate Faculty of the Department of
Mathematics in the College of Liberal Arts and Sciences and to the Graduate School
and was accepted as partial fulfillment of the requirements for the degree of Doctor
of Philosophy.

May 2000

UNIVERSITY OF FLORIDA

3 1262 08555 361W
I .. . |

I. / "'