Circulant Matrices and Their Application to Vibration Analysis

This paper provides a tutorial and summary of the theory of circulant matrices and their application to the modeling and analysis of the free and forced vibration of mechanical structures with cyclic symmetry. Our presentation of the basic theory is distilled from the classic book of Davis (1979, Circulant Matrices , 2nd ed., Wiley, New York) with results, proofs, and examples geared speciﬁcally to vibration applications. Our aim is to collect the most relevant results of the existing theory in a single paper, couch the mathematics in a form that is accessible to the vibrations analyst, and provide examples to highlight key concepts. A nonexhaustive survey of the relevant literature is also included, which can be used for further examples and to point the reader to important extensions, applications, and generalizations of the theory. [DOI: 10.1115/1.4027722]


Introduction
"The theory of matrices exhibits much that is visually attractive.Thus, diagonal matrices, symmetric matrices, (0, 1) matrices, and the like are attractive independently of their applications.In the same category are the circulants."Philip J. Davis The modeling and analysis of structural vibration is a mature field, especially when vibration amplitudes are small and linear models apply.For such models, the powerful tools of modal analysis and superposition allow one to decompose the system and its response into a set of uncoupled single degree-of-freedom (DOF) systems, each of which captures the motion of the overall system in a given normal mode.For geometrically simple continua, the natural frequencies, vibration modes, and response to known excitations can be analytically determined.When discrete models are developed, matrix methods are readily applied for both natural frequency and response analyses.For general system models with no special properties, which occur in the majority of applications, largescale computational models must be developed, typically using finite element methods.However, certain classes of systems possess special properties, such as symmetries, which aid in the analysis by enabling significant reduction of the finite element models (Fig. 1(a)).This is particularly true for systems with cyclic symmetry.
The main goal of this tutorial is to consider the vibrations of structural systems with cyclic symmetry, also known as rotationally periodic systems.A useful geometric view of these systems is that of a circular disk (a pie) split into N equal sectors (i.e., equally sized pieces of the pie), each of which contains an identical mechanical structure with identical coupling to forward-nearestneighbors and to ground (Fig. 1(b)).
The theory of circulants also applies to more general forms of coupling with non-nearest-neighbors, for example through a base substructure, as long as the rotational symmetry is preserved.These structural arrangements arise naturally in certain types of rotating machines.Turbomachinery examples include bladed disks, such as fans, compressors, turbines, and impeller stages of aircraft, helicopter engines, power plants, as well as propellers, pumps, and the like.Rotational periodicity also arises in some stationary structures such as satellite antennae.In some systems, such as planetary gears and rotors with pendulum vibration absorbers, the overall system is not necessarily cyclic, but subcomponents of it may be.When perfect cyclic symmetry of a model is assumed, special modal properties exist that significantly facilitate its vibration analysis.
The nature of rotationally periodic systems imposes a cyclic structure on their mass and stiffness matrices, which are block circulant for systems with many DOFs per sector (Fig. 1(a)) and circulant for the special case of a single DOF per sector (Fig. 1(b)).By denoting the stiffness of the internal elements of each sector by K 0 and the coupling stiffness between sectors as -K 1 , the stiffness matrices of rotationally periodic structures with nearestneighbor coupling have the general form where K 0 and K 1 are themselves matrices for the complex model shown in Fig. 1(a) and scalars for the simplest prototypical model shown in Fig. 1(b).A key property of K is that the elements of each row are obtained from the previous row by cyclically permuting its entries.That is, for j ¼ 2; 3; …; N, row j is obtained from row j -1 by shifting the elements of row j -1 to the right by one position and wrapping the right-end element of row j -1 into the first position.This is precisely the form of a circulant matrix, which is formally defined in Sec.2.2.The mass matrix of a rotationally periodic structure with nearest-neighbor coupling is block diagonal and also shares this cyclic property.The size of the elements of K is equal to the number of DOFs per sector, and is denoted by M. Thus, a system with N sectors and M DOFs per sector has a total of NM DOFs.The most important utility of the theory of circulants in analyzing rotationally periodic systems is that they enable a NM-DOF system to be decomposed to a set of NM-DOF uncoupled systems using the appropriate coordinate transformation.Admittedly, the same can be accomplished using brute-force methods to uncouple the entire system using modal analysis, but such an approach overlooks fundamental properties that are crucial to understanding the free and forced response of these systems and requires significantly more computational power.This is the central motivation for understanding and utilizing circulants to analyze cyclic systems.
The vibration modes of rotationally periodic systems consist of multiple pairs of repeated natural frequencies (eigenvalues) that lead to pairs of degenerate normal modes (eigenvectors).The number and nature of such pairs depend on whether N is even or odd.Each mode pair is characterized as a pair of standing waves (SWs) with different spatial phases, or a pair of traveling waves, labeled as a forward traveling wave (FTW) and backward traveling wave (BTW) when following the terminology used in applications to rotating machinery.The choice of formulation is based on convenience for a given application, which depends on the nature of the system excitation.For example, the excitation frequency is proportional to the engine speed for many cyclic rotating systems, which leads to the so-called engine order (e.o.) excitation, and often the spatial nature of the excitation (in the rotating frame of reference) is in the form of a traveling wave.When such excitation is applied to systems with cyclic symmetry, the response also has special properties that can be easily uncovered by making use of the system traveling wave vibration modes.
The strength of the intersector coupling is an important parameter in rotationally periodic systems.When the intersector coupling is strong, the frequencies of the mode pairs are well separated.In contrast, weak intersector coupling yields closely spaced frequencies, high modal density, and large sensitivity to cyclic-symmetry-breaking imperfections.A wave representation of the response [2] shows that the strength of the coupling determines frequency passband widths, wherein unattenuated propagation of waves takes place.Weak intersector coupling leads to narrow passbands, and the passbands widen as the coupling strength increases.Another important parameter for cyclically symmetric structures is the total number of sectors.The modal density is larger for large N, which corresponds to more natural frequencies within each frequency passband.In all cases, the modes are spatially distributed, or extended, for models of cyclic systems.That is, the pattern of displacements in a modal response is uniformly spread around the circumference of the structure.
The special properties of systems with cyclic symmetry extend to nonlinear systems, where they are expressed quite naturally in terms of symmetry groups: the cyclic group, in particular [65][66][67][68].The group theoretic formulation can also be applied to the linear vibration problems considered in this paper [69,70], but the approach presented here is more approachable to readers with a standard engineering background in linear algebra.
The extension to systems with small imperfections that perturb the cyclic symmetry has led to important results related to mode localization, which arises in systems with high modal density caused by weak intersector coupling or a large number of sectors.In particular parameter regimes, the mode shapes are highly sensitive to small, symmetry-breaking imperfections among the nominally identical sectors, and the spatial nature of the vibration modes can become highly localized.For these cases, the vibration energy is focused in a small number of sectors, and sometimes even a single sector.This behavior, which stems from the seminal work of Anderson on lattices [71], was originally recognized to be relevant to structural vibrations by Hodges and Woodhouse [72,73] and Pierre and Dowell [74], and has been extensively studied from both fundamental [75] and applied [76][77][78] points of view.The phenomenon of mode localization is also observed in the forced response and has practical implications for the fatigue life of bladed disks in turbomachinery [79,80].It is interesting to note that localization also arises in nonlinear systems with perfect symmetry, where the dependence of the system natural frequencies on the amplitudes of vibration naturally leads to the possibility of mistuning of frequencies between sectors if their amplitudes are different [81][82][83][84][85][86][87][88].
Another topic central to vibration analysis that relies on the theory of circulants is the discrete Fourier transform (DFT) [89,90].The DFT was known to Gauss [91], and is the most common tool used to process vibration signals from experimental measurements and numerical simulations.The DFT and inverse DFT (IDFT) provide a computationally convenient means of determining the frequency content of a given signal.Because the mathematics of circulants is at the heart of the computation of the DFT, we include a brief introduction to the relationship between the DFT and IDFT, and its connection to the theory of circulants.
The goal of this paper is to provide a detailed theory of circulant matrices as it applies to the analysis of free and forced structural vibrations.Much of the material was developed as part of the Ph.D. research of the lead author [21,22,[92][93][94][95]. References to other relevant work are included throughout this paper, but we do not claim to provide an exhaustive survey of the relevant literature.The remainder of the paper is organized as follows.Section 2 gives a quite exhaustive and self-contained treatment of the theory of circulants, which is distilled from the seminal work by Davis [96].We adopt a presentation style similar to that of Ottarsson [97], one that should be familiar to an analyst in the vibrations engineering community.This section is meant to act simultaneously as a detailed reference and tutorial, including proofs of the main results and simple illustrative examples.Section 3 provides three examples that make use of the theory, including ordinary circulants and the more general block circulant matrices.Particular attention is given to cyclic systems under traveling wave engine order excitation because this type of system forcing appears naturally in many relevant applications of rotating machinery.The apprised reader, or the reader who wishes to learn by example, can skip directly to Sec. 3, depending on their background, and revisit Sec. 2 as warranted.The paper closes with a brief summary in Sec. 4.

The Theory of Circulants
This section details the theory and mathematics of circulant matrices that are relevant to vibration analysis of mechanical structures with cyclic symmetry.The basic theory is distilled from the seminal work by Davis [96] and is presented using mathematics and notation that should be familiar to the vibrations engineer.Selected topics from linear algebra are reviewed in Sec.2.1 to introduce relevant notion and support the theoretical development of circulant matrices in Secs.2.2-2.8.This material is included for completeness; the apprised reader can skip directly to Secs.2.2 and 2.3, where circulant and block circulant matrices (also referred to as circulants and block circulants) are defined.Representations of circulants are discussed in Sec.2.4.Diagonalization of circulants and block circulants is discussed at length in Sec.2.5, which begins with a treatment of the Nth roots of unity in Sec.2.5.1 and the Fourier matrix in Sec.2.5.2.It is subsequently shown how to diagonalize the cyclic forward shift matrix in Sec.2.5.3 a circulant in Sec.2.5.4,and a block circulant in Sec.2.5.5.Some generalizations of the theory are discussed in Sec.2.6, including the diagonalization of block circulants with circulant blocks.Relevant mathematics of the DFT and IDFT are summarized in Sec.2.7.Finally, the circulant eigenvalue problem (cEVP) is discussed in Sec.2.8, including the eigenvalues and eigenvectors of circulants and block circulants, their symmetry characteristics, and connection to the DFT process.
2.1 Mathematical Preliminaries.Definitions and relevant properties of special operators and matrices are discussed in Secs.2.1.1 and 2.1.2,respectively, including the direct (Kronecker) product, and Hermitian, unitary, cyclic forward shift, and flip matrices.This is followed in Sec.2.1.3with a treatment of matrix diagonalizability.
2.1.1Special Operators.Let C denote the set of complex numbers and Z þ be the set of positive integers.
DEFINITION 1 (Direct Sum).For each i ¼ 1; 2; …; N and p i 2 Z þ , let A i 2 C piÂpi .Then the direct sum of A i is denoted by and results in the block diagonal square matrix of order p 1 þ p 2 þ Á Á Á þ p N , where each zero matrix 0 has the appropriate dimension.᭝ It is convenient to define the operator diagðÁÞ that takes as its argument the ordered set of matrices A 1 ; A 2 ; …; A N and results in the block diagonal matrix given in Definition 1, that is, Then the direct product (or Kronecker product) of a and b T is the square matrix

!
Then the direct product of A and B is given by Some important properties of the direct product are as follows: (1) The direct product is a bilinear operator.If A and B are square matrices and a is a scalar, then (2) The direct product distributes over addition.If A, B, and C are square matrices with the same dimension, then (3) The direct product is associative.If A, B, and C are square matrices, then (4) The product of two direct products yields another direct product.If A, B, C, and D are square matrices such that AC and BD exist, then ðA BÞðC DÞ ¼ ðACÞ ðBDÞ (5) The inverse of a direct product yields the direct product of two matrix inverses.If A and B are invertible matrices, then where ðÁÞ À1 denotes the matrix inverse.
(6) The transpose or conjugate transpose of a direct product yields the direct product of two transposes or conjugate transposes.If A and B are square matrices, then where ðÁÞ H ¼ ðÁÞ T is the conjugate transpose and ðÁÞ denotes complex conjugation.where detðÁÞ and trðÁÞ denote the matrix determinant and trace.
2.1.2Special Matrices.The definitions and relevant properties of selected special matrices are summarized.Hermitian and unitary matrices are defined first (see Table 1), followed by a brief treatment of two important permutation matrices: the cyclic forward shift matrix and the flip matrix.The details of circulant matrices and the Fourier matrix, which are employed extensively throughout this work, are deferred to Secs.2.2, 2.3, and 2.5.2.

DEFINITION 3 (Hermitian Matrix
᭝ The elements of a Hermitian matrix H satisfy h ik ¼ h ki for all i; k ¼ 1; 2; …; N. Thus, the diagonal elements h ii of a Hermitian matrix must be real, while the off-diagonal elements may be complex.If H ¼ H T then H is said to be symmetric.

DEFINITION 4 (Unitary Matrix
᭝ Real unitary matrices are orthogonal matrices.If a matrix U is unitary, then so too is U H .To see this, consider ðU H Þ H ðU H Þ ¼ UU H ¼ I, from which it follows that Finally, if U is unitary and nonsingular, then U H ¼ U À1 .A general permutation matrix is formed from the identity matrix by reordering its columns or rows.Here, we introduce two such matrices: the cyclic forward shift matrix and the flip matrix. DEFINITION 5 (Cyclic Forward Shift Matrix).The N Â N cyclic forward shift matrix is given by   cyclically shifts the entries of a by one entry to the right.That is, the ith entry of a is shifted to entry i þ 1, except for entry N ¼ 3, which is placed in position 1 of a. DEFINITION 6 (Flip Matrix).The N Â N flip matrix is given by NÂN which is populated with ones in the (1, 1) position and along the subantidiagonal, and zeros otherwise.᭝ COROLLARY 1.Let j N be the flip matrix.Then where I N is the N Â N identity matrix.ٗ 2.1.3Matrix Diagonalizability.Matrix diagonalization is the process of taking a square matrix and transforming it into a diagonal matrix that shares the same fundamental properties of the underlying matrix, such as its characteristic polynomial, trace, and determinant.This section defines matrix diagonalizability in terms of similarity, provides necessary conditions for a matrix to be diagonalizable, and summarizes relevant properties of diagonalizable matrices.Diagonalization of circulant matrices is deferred to Sec. 2.5.
DEFINITION 7 (Similarity Transformation).Let Q be an arbitrary nonsingular matrix.Then B ¼ Q À1 AQ is a similarity transformation and B is said to be similar to A.
It therefore suffices to say that A and B are similar.A summary of selected additional linear transformations is provided in Table 2.If B is orthogonally (resp.unitarily) similar to A, then we say that A and B are orthogonally (resp.unitarily) similar matrices.THEOREM 1.If A and B are similar matrices, then they have the same characteristic equation and hence the same eigenvalues.ٗ Theorem 1 guarantees that the eigenvalues of a matrix are preserved under a similarity transformation.A proof can be found in any standard textbook on linear algebra [98,99].Because unitary Q, the eigenvalues are also preserved under orthogonal and unitary transformations.
THEOREM 2. Let the matrices A and B be similar.Then if denotes a finite polynomial in t with arbitrary coefficients c k ðk ¼ 1; 2; …; NÞ, the matrix polynomials p(A) and p(B) are similar.ٗ Proof.Let Q be an arbitrary nonsingular matrix.Then which completes the proof.If p(t) ¼ t k with k > 0 in Theorem 2, then we have the following result.
᭝ Thus, a matrix is diagonalizable if it is similar to a diagonal matrix.If A is diagonalizable by Q, we say that Q diagonalizes A and that Q is the diagonalizing matrix.
THEOREM 3.An N Â N matrix A is diagonalizable if it has N linearly independent eigenvectors.ٗ Proof.Suppose A has N linearly independent eigenvectors and denote them by q 1 ; q 2 ; …; q N .Let k i be the eigenvalue of A corresponding to q i for each i ¼ 1; …; N. Then if Q is the matrix that has as its ith column the vector q i , it follows that

Circulant Matrices
DEFINITION 9 (Circulant Matrix).A N Â N circulant matrix (or circulant, or ordinary circulant) is generated from the N-vector fc 1 ; c 2 ; …; c N g by cyclically permuting its entries, and is of the form ᭝ Thus, a circulant matrix is defined completely by the generating elements in its first row, which are cyclically shifted to the right by one position per row to form the subsequent rows.The set of all such matrices of order N is denoted by C N .A matrix contained in C N is said to be a circulant of type N.
It is convenient to define the circulant operator circðÁÞ that takes as its argument the generating elements c 1 ; c 2 ; …; c N and results in the array given in Definition 9, that is, An N Â N circulant is also characterized in terms of its (i, k) entry by ðCÞ ik ¼ c kÀiþ1ðmodNÞ for i; k ¼ 1; 2; …; N. If a matrix is both circulant and symmetric, its generating elements are c 1 ; …; cN 2 ; cNþ2 2 ; cN 2 ; …; c 3 ; c 2 ; N even c 1 ; …; cNÀ1 2 ; cNþ1 2 ; cNþ1 2 ; cNÀ1 2 ; …; c 3 ; c 2 ; N odd ( These matrices are equal if and only if the equalities are satisfied.Then C can be written as which implies that Ar Next we introduce block circulant matrices, which are natural generalizations of ordinary circulants.

Block Circulant Matrices.
A block circulant matrix is obtained from a circulant matrix by replacing each entry c k in Definition 9 by the M Â M matrix C i for i ¼ 1; 2; …; N.
DEFINITION 11 (Block Circulant Matrix).Let C i be a M Â M matrix for each i ¼ 1; 2; …; N. Then a NM Â NM block circulant matrix (or block circulant) is generated from the ordered set fC 1 ; C 2 ; …; C N g, and is of the form ᭝ A block circulant is therefore defined completely by its generating matrices.The matrix array given by Definition 11 is said to be a block circulant of type (M, N).The set of all such matrices is denoted by BC M;N .A matrix C 2 BC M;N is not necessarily a circulant, as the following example demonstrates.
Example 6.Let  is a block circulant of type (2, 4), but it is not a circulant.
Next we give a necessary and sufficient condition for a matrix to be block circulant.
THEOREM 5. Let r N be the cyclic forward shift matrix of dimension N and I M be the identity matrix of dimension M. Then a NM Â NM matrix C is a block circulant of type (M, N) if and only if ٗ The proof of Theorem 5 follows similarly to that of Theorem 4 by replacing the scalar elements c ik with M Â M matrices C ik for i; k ¼ 1; 2; …; N. The reader can verify that the matrix C in Example 6 satisfies the condition in Theorem 5, but not that in Theorem 4.
A block circulant, block symmetric matrix of type (M, N) has generating matrices of the same form as Eq. ( 10), and is obtained by replacing each entry c k by the M Â M matrix C k for k ¼ 1; 2; …; N. The set of all such matrices is denoted by BCBS M;N .The matrix C in Example 6 is recognized to be a block symmetric, block circulant matrix of type (2,4), that is, it is contained in BCBS 2;4 , which is a subset of BC 2;4 .
where each successive power cyclically permutes the generating elements.This enables the representation of a general circulant in terms of a finite matrix polynomial involving the cyclic forward shift matrix and its powers.COROLLARY 3. Let C 2 C N be a circulant matrix of type N with generating elements c 1 ; c 2 ; …; c N .Then C can be represented by the matrix sum where r N is the N Â N cyclic forward shift matrix.ٗ Example 7. The matrix A ¼ circða; b; cÞ from Example 5 can be represented by the matrix sum where I 3 and r 3 are the 3 Â 3 identity and cyclic forward shift matrices.Corollary 3 is exploited in Sec.2.5.4 to diagonalize a general circulant matrix, and can be generalized to represent a general block circulant matrix in terms of the cyclic forward shift matrix and its powers.
COROLLARY 4. Let C 2 BC M;N be a block circulant matrix of type (M, N) with generating matrices C 1 ; C 2 ; …; C N .Then C can be represented by the matrix sum where r N is the cyclic forward shift matrix.ٗ Corollaries 3 and 4 motivate the following result, which captures the representation of circulant and block circulant matrices in terms of the cyclic forward shift matrix, and facilitates their diagonalization in Secs.2.5.4 and 2.5.5.
DEFINITION 13.Let t and s be arbitrary square matrices.Then the function qðt; sÞ ¼ X N k¼1 t kÀ1 s is a finite sum of direct products.᭝ COROLLARY 5. Let r N be the N Â N cyclic forward shift matrix.Then a circulant matrix with generating elements c 1 ; c 2 ; …; c N and a block circulant matrix with generating elements C 1 ; C 2 ; …; C N can be represented by the matrix sums where the function qðÁÞ is given by Definition 13. ٗ What is meant by the notation qðr N ; c k Þ, for example, is to substitute t with r N and s with c k in Definition 13 and then perform the summation observing any indices k introduced by the substitution.
2.5 Diagonalization of Circulants.Any circulant or block circulant matrix can be represented in terms of the cyclic forward shift matrix according to Corollary 5.The diagonalization of a general circulant begins, therefore, by finding a matrix that diagonalizes r N .Together with some basic results from linear algebra (these are summarized in Sec.2.1), this leads naturally to the diagonalization of an arbitrary circulant.Regarding a suitable diagonalizing matrix, there are a number of candidates [10][11][12]16,17,100], but all feature powers of the Nth roots of unity or their real/imaginary parts.In this work, we employ an array composed of the distinct Nth roots of unity (Sec.2.5.1) and their integer powers in the form of the complex Fourier matrix (Sec.2.5.2).A unitary transformation involving the Fourier matrix is used to diagonalize the cyclic forward shift matrix (Sec.2.5.3),general circulant matrices (Sec.2.5.4), and general block circulant matrices (Sec.2.5.5).

Nth Roots of Unity.
A root of unity is any complex number that results in 1 when raised to some integer N 2 Z þ [101].More generally, the Nth roots of a complex number z o ¼ r o e jho are given by a nonzero number z ¼ re jh such that z N ¼ z o or r N e jNh ¼ r o e jho (12) where j ¼ ffiffiffiffiffiffi ffi À1 p .Equation ( 12) holds if and only if r N ¼ r o and and the Nth roots are Equation (14) shows that the roots all lie on a circle of radius ffiffiffiffi r o N p centered at the origin in the complex plane, and that they are equally distributed every 2p=N radians.Thus, all of the distinct roots correspond to k ¼ 0; 1; 2; …; N À 1. DEFINITION 14 (Distinct Nth Roots of Unity).The distinct Nth roots of unity follow from Eq. ( 14) by setting r o ¼ 1 and h o ¼ 0 and are denoted by Example plots of the distinct Nth roots of unity are shown in Fig. 2, where w ðkÞ N are arranged on the unit circle in the complex plane (centered at the origin) for N ¼ 1; 2; …; 9. Note that w ð0Þ if N is even.The remaining roots appear in complex conjugate pairs.Thus, the distinct Nth roots of unity are symmetric about the real axis in the complex plane.

2.5.2
The Fourier Matrix.This section introduces the complex Fourier matrix and its relevant properties, including the symmetric structure of the N-vectors that compose its columns.A key result is that the Fourier matrix is unitary, which is systematically developed and proved.
DEFINITION 17 (Fourier Matrix).The N Â N Fourier matrix is defined as  Clearly the Fourier matrix is symmetric, but generally it is not Hermitian.It can be written element wise as where is the angle subtended from the positive real axis in the complex plane to the ith power of w N , which is also the ði þ 1Þth of the N roots of unity according to Definition 14 and the numbering scheme in Fig. 2. COROLLARY 7. The matrix E H N is obtained from E N by changing the signs of the powers of each element.ٗ Proof.The Fourier matrix is unaffected by transposition because it is symmetric.Thus, the (i, k) element of E H N is where the identity N can be obtained from E N by changing the sign of the powers of the Nth roots of unity.
It is shown in Sec.2.8 that all circulant matrices contained in C N share the same linearly independent eigenvectors, the elements of which compose the N columns (or rows) of E N .DEFINITION 18.Let w N be the primitive Nth root of unity.Then for i ¼ 1; 2; …; N the columns of the Fourier matrix E N are defined by where u i is given by Eq. ( 16).᭝ Example 10.The third column of the Fourier matrix E 4 is given by The columns e i of the Fourier matrix E N ¼ ðe 1 ; …; e N Þ exhibit a symmetric structure with respect to the index i ¼ (N þ 2)/2 for even N.In Example 9, for instance, the vectors e 1 and e ðNþ2Þ=2 ¼ e 3 are real and distinct, and the vectors e 2 and e 4 appear in complex conjugate pairs.This same structure generally holds for any E N with even N.To see this, consider where the integers 6q correspond to vector pairs relative to the index i ¼ ðN þ 2Þ=2 and the identity N is employed.The case of q ¼ 0 corresponds to i ¼ ðN þ 2Þ=2 and yields the real vector which has the same value for each element with alternating signs from element to element.For q 6 ¼ 0 the terms w 6q N are complex conjugates according to the proof of Corollary 7 such that e ðNþ2Þ=2 ð Þ þ q and e ðNþ2Þ=2 ð Þ À q are complex conjugate pairs.There are (N À 2)/2 such pairs corresponding to q ¼ f1; 2; …; ðN À 2Þ=2g in Eq. (17).Finally, the case of i ¼ 1 always yields the real vector for even and odd N. A similar formulation for odd N shows that e 1 is real and distinct and the remaining (N -1)/2 vectors appear in complex conjugate pairs.This is shown by example in Sec.3.3 in the context of vibration modes for a cyclic structure with a single DOF per sector.
A key feature of the Fourier matrix is that it is unitary.This is essentially a statement of orthogonality of each column of E N and is captured by the following lemmas.
LEMMA 1 (Finite Geometric Series Identity).Let N 2 Z þ and q 2 C. Then X sþNÀ1 r¼s q r ¼ q s ð1 À q N Þ 1 À q for any s 2 Z and q 6 ¼ 1. ٗ Proof.Consider the finite geometric series Subtraction of the second equation from the first results in ð1 À qÞ X sþNÀ1 r¼s q r ¼ q s ð1 À q N Þ from which the proof is established by division of the term (1q) because q 6 ¼ 1 by restriction.Lemma 1 is used to establish the following result, which is required to show that the Fourier matrix is unitary.The orthogonality condition is fundamental to the diagonalization of circulants in Secs.2.5.3-2.5.5, and

¼ e j 2p
N ðiÀkÞ and note that q N ¼ 1.If iÀk ¼ mN, then q ¼ e j2pm ¼ 1 for any integer m, and it follows that For the case of i À k 6 ¼ mN it follows from Lemma 1 that which completes the proof.Example 11.Consider the orthogonality condition given by Lemma 2 for s ¼ 0.
which is numerically equal to N, as expected.If instead we set i À k ¼ 5 but N ¼ 4, then Lemma 2 allows for representations of the N Â N identity, flip, and cyclic forward shift matrices in terms of certain conditions on their indices relative to N. COROLLARY 8.For i; k ¼ 1; …; N and any integer m, the (i, k) elements of the N Â N identity, flip, and cyclic forward shift matrices can be represented by the summations where d ik is the Kronecker delta.ٗ The reader can verify Corollary 8 for the special case of N ¼ 3 by inspection of the arrays in Fig. 3.
We are now ready to state the key result required to diagonalize a general circulant matrix.Corollaries 7 and 8 are used to show that the Fourier matrix is unitary.
THEOREM 6 (Unitary Fourier Matrix).The Fourier matrix E N is unitary.
The following corollaries follow from Theorem 6. COROLLARY 9.If E N is the unitary Fourier matrix, then which implies that E N is unitary.Example 13.Consider the unitary matrix E 4 from Example 12. Corollary 9 guarantees that E 4 is also unitary.Thus, 1 Àj À1 j where I N is the N Â N identity matrix.Thus, e H i e k ¼ d ik .Corollary shows that the columns (and rows) of the Fourier matrix are mutually orthogonal.The result can also be proved by expanding the product e H i e k according to Eq. ( 15) and invoking Lemma 2, as is done in the proof of Theorem 6.The same result also follows by expanding . .
corresponds to d 33 ¼ 1 in Corollary 10.However, the product where I N and I NM are identity matrices of dimension N and NM, respectively.COROLLARY 12.If e i denotes the ith column of the Fourier matrix E N , I M is the identity matrix of dimension M, and d ik is the Kronecker delta, then the NM Â M matrices e i I M are such that which completes the proof.
Corollary 12 can also be obtained directly from Corollary 11 by writing . .and expanding these matrices similarly to Eq. ( 20), and setting the result equal to Next we derive a relationship between the Fourier and flip matrices.
THEOREM 7. Let E N and j N denote the N Â N Fourier and flip matrices, respectively.Then For any integer m and i; k ¼ 1; 2; …; N, the (i, k) entry of E N E N is given by follows from complex conjugation and transposition of j N ¼ E N E N , and by invoking the properties j H A number of properties follow directly from Theorem 7. COROLLARY 13.Let E N and j N be the N Â N Fourier and flip matrices.Then where I N is the N Â N identity matrix.ٗ Property (a) of Corollary 13 says that the flip and Fourier matrices commute or, since E N is unitary, that j N is invariant under a unitary transformation with respect to E N .Thus, j N is not diagonalizable by E N .Properties (b) and (c) give alternative definitions of the flip and Fourier matrices.Moreover, because the power of a diagonal matrix is obtained by raising each diagonal element to the power in question, if follows that the eigenvalues of j N are 61 and those of E N are 61 and 6j, each with the appropriate multiplicities.

Diagonalization of the Cyclic Forward Shift Matrix.
In light of Corollary 5, diagonalization of a general circulant or block circulant matrix begins by diagonalizing the cyclic forward shift matrix.THEOREM 8. Let E N be the N Â N Fourier matrix and r N be the N Â N cyclic forward shift matrix.Then ÀðrÀ1ÞðkÀ1Þ N (Eq.15) The desired result follows by multiplying from the left by E H N , multiplying from the right by E N , and invoking Theorem 6.
Theorem 8 implies that r N is unitarily similar to a diagonal matrix whose diagonal elements are the distinct Nth root of unity (i.e., Definition 16).Because the eigenvalues of a matrix are preserved under such a transformation (this is guaranteed by Theorem 1), it follows that where aðÁÞ denotes the matrix spectrum.The eigenvectors of the circulant matrix r N are the linearly independent columns of where X N given by Definition 16. ٗ COROLLARY 15.Let e i be the ith column of the Fourier matrix E N and r N be cyclic forward shift matrix.Then for any n 2 Z þ and i; k ¼ 1; 2; … N, where d ik is the Kronecker delta.
ٗ Corollary 15 shows that the columns of the Fourier matrix are orthogonal with respect to the cyclic forward shift matrix and its powers.It is shown in Sec.2.5.4 that they are in fact orthogonal with respect to any circulant matrix.

Diagonalization of a Circulant.
Corollary 16.Let the matrix X N be populated with the distinct Nth roots of unity according to Definition 16 and s be an arbitrary square matrix.Then where qðÁÞ is given by Definition 13. ٗ Proof.Consider the representation which is diagonal when s is a scalar and block diagonal when s is a matrix.THEOREM 9 (Diagonalization of a Circulant).
is a diagonal matrix.For i ¼ 1; 2; …; N, the diagonal elements are where w N is the primitive Nth root of unity and the function qðÁÞ is given by Definition 13. ٗ Proof.Consider the representation where the last step follows directly from the proof of Theorem 2 and the polynomial term Theorem 9 shows that the Fourier matrix E N diagonalizes any N Â N circulant matrix.As discussed in Sec.2.8, the columns e 1 ; e 2 ; …; e N of E N are the eigenvectors of C. The scalars k i in Theorem 9 are the eigenvalues of C. Unlike the eigenvectors, they depend on the elements (i.e., generating elements) of C. In fact, Eq. ( 21) has the same form as the DFT of a discrete time series, which is clear by comparing it to Definition 20 of Sec.2.7.In this case, c k ðk ¼ 1; 2; …; NÞ and k i ði ¼ 1; 2; …; NÞ are analogous to a discrete signal and its DFT, and are related by which has the same form as the matrix-vector representation of the DFT defined by Eq. ( 27) of Sec.2.7.Thus, the eigenvalues k i can be calculated directly using Eq. ( 21) or Eq. ( 22), which are equivalent.
If the circulant matrix in Theorem 9 is also symmetric, the eigenvalues are real-valued and certain ones are repeated, as shown in Corollary 17.
COROLLARY 17.If C 2 SC N is a symmetric circulant, Eq. ( 21) reduces to where the generating elements are given by Eq. ( 10).ٗ Proof.If N is even, the generating elements of C are given by c 1 ; c 2 ; …; cN 2 ; cNþ2 2 ; cN 2 ; …; c 3 ; c 2 and Eq. ( 21) reduces to where the identity is employed.If N is odd, the generating elements of C are given by c 1 ; c 2 ; …; cNÀ1 2 ; cNþ1 2 ; cNþ1 2 ; cNÀ1 2 ; …; c 3 ; c 2 and Eq. ( 21) reduces similarly to the case for even N, which is left as an exercise for the reader.
The cosinusoidal nature of k i in Corollary 17 implies that k 1 is distinct for a symmetric circulant matrix C 2 SC N , but the remaining elements k i ¼ k Nþ2Ài appear in repeated pairs.However, k Nþ2=2 ð Þ is also distinct if N is even.Example 15.Let C ¼ circð4; À1; 0; À1Þ 2 SC 4 be a symmetric circulant matrix.Then the product   is a diagonal matrix.The diagonal elements can be computed directly using Eq. ( 21) for N ¼ 4. For example, which is recognized to be the third diagonal element of the matrix product Because C is a symmetric circulant matrix, Corollary 17 can also be used.Observing that N ¼ 4 is even, it follows that    a certain symmetry about the so-called "Nyquist" component, which is analogous to the DFT of a real-valued sequence.This is because Eq. ( 22) represents the DFT of the generating elements c 1 ; c 2 ; …; c N .COROLLARY 18.Let e i be the ith column of the Fourier matrix E N and C be a N Â N circulant matrix.Then for i; k ¼ 1; 2; …; N, is a NM Â NM block diagonal matrix.For i ¼ 1; 2; …; N, the M Â M diagonal blocks are where w N is the primitive Nth root of unity and the function qðÁÞ is given by Definition 13. ٗ Proof.Consider the representation where    23).Then if each C i is symmetric, it follows that K i is symmetric for i ¼ 1; 2; …; N. ٗ Proof.If each C k is symmetric for k ¼ 1; 2; …; N, then so too are the matrices C k w ðkÀ1ÞðiÀ1Þ N for each i because w ðkÀ1ÞðiÀ1Þ N is a scalar.Moreover, the sum and difference of two symmetric matrices is again symmetric.If follows that is symmetric for i ¼ 1; 2; …; N. COROLLARY 20.Let C 2 BC M;N be a NM Â NM block circulant matrix.Then for i; k ¼ 1; 2; …; N, where K i is defined by Eq. ( 23).ٗ Example 19.Consider C ¼ circðA; B; 0; BÞ 2 BC 2;4 from Examples 6 and 18.Then ! is the third 2 Â 2 block of the block diagonal matrix obtained in Example 18.
2.6 Generalizations.Let ðÁÞ Ã denote an arbitrary operation that takes a square matrix as its argument and returns another square matrix with the same dimension.For example, the operation could denote a matrix inverse such that ðÁÞ Ã ¼ ðÁÞ À1 .Let ðÁÞ # be another arbitrary matrix operation with the same restrictions (i.e., returns another square matrix with the same dimension).Then if C 2 BC M;N has generating matrices C 1 ; C 2 ; …; C N , it follows that for any matrices A 2 C NÂN and B 2 C MÂM .The importance of this result is that C can be decomposed into a summation of direct products of two separate equivalence transformations, one that operates on the cyclic forward shift matrix and the other on the generating matrices of C.This decomposition justifies the diagonalizing matrix used in Sec.2.5, motivates some generalizations of Theorem 10, and aids in proving orthogonality relationships for the cyclic eigenvalue problems described in Sec.2.8.
In light of Corollary 14, it is clear that the choice of A ¼ E N and ðÁÞ Ã ¼ ðÁÞ H accomplishes block diagonalization of a matrix C 2 BC M;N .Then if B ¼ I M , the appropriate diagonalizing matrix to block decouple C without operating on its generating matrices is E N I M (see Theorem 10).However, if B and ðÁÞ # are kept general, we have the following result.
THEOREM 11.Let C 2 BC M;N have M Â M generating matrices C 1 ; C 2 ; …; C N and E N be the N Â N Fourier matrix.Then for an arbitrary matrix B 2 C MÂM and operator ðÁÞ # , is a block diagonal matrix, where is the ith M Â M diagonal block for i ¼ 1; 2; …; N. ٗ COROLLARY 21.Let C 2 BC M;N be a NM Â NM block circulant matrix, B 2 C MÂM be an arbitrary square matrix, and ðÁÞ # denote an arbitrary operation that takes a square matrix as its argument and returns another square matrix with the same dimension.Then for i; k ¼ 1; 2; …; N, where W i is defined by Eq. ( 24).ٗ Theorem 11 is useful if there exists an equivalence transformation B # C k B that simplifies each of the generating matrices.is a NM Â NM diagonal matrix, where is the pth diagonal element of the ith M Â M block for i ¼ 1; 2; …; N and p ¼ 1; 2; …; M. ٗ Corollary 22 shows that the NM-dimensional eigenvectors q ðpÞ i of C are the columns of E N E M .This is in contrast to rotationally periodic structures, where q is partitioned into N M-vectors corresponding to each sector and decomposed into a set of N reduced-order eigenvalue problems described in Sec.

Relationship to the Discrete Fourier Transform.
Before turning to the circulant eigenvalue problem, we consider a somewhat tangential but relevant subject on the DFT and its inverse, which are central to the analysis of experimental data in a wide range of fields, including mechanical vibrations.In this detour we show that computation of the DFT is, in fact, a multiplication of an N-vector of discrete signal samples by the Fourier matrix and a constant c f .The IDFT is similarly defined using the Hermitian of the Fourier matrix and a constant c i , where c f c i ¼ 1=N.More importantly, it is shown that the DFT computation is exactly analogous to the determination of the eigenvalues of a circulant matrix given its generating elements.We begin by defining the DFT sinusoids, which provide a convenient means of representing the DFT of a discretized signal, and then develop the relationships of interest for the DFT and the IDFT.We present only the basic results as they relate to the theory and mathematics of circulants.The reader can find a vast literature on related topics [89][90][91][102][103][104][105][106][107][108][109][110].
DEFINITION 19 (DFT Sinusoids).Let w k N denote the distinct Nth roots of unity, where w N ¼ e j 2p N is the primitive root.Then the DFT sinusoids are The inner product of S i (r) and S k (r) is given by which shows that the DFT sinusoids are orthogonal.The Fourier matrix in Definition 17 can be written element wise as for each i; k ¼ 1; 2; …; N. It follows that is a representation of the Fourier matrix in terms of the DFT sinusoids, where W N is the DFT sinusoid matrix.The DFT is formally defined next and subsequently reformulated in terms of a matrix multiplication involving W N .DEFINITION 20 (DFT).Let x(i) denote a finite sequence with indices i ¼ 1; 2; …; N. Then for k ¼ 1; 2; …; N, the DFT of x(i) is the sequence where w N is the primitive Nth root of unity and S i ðkÞ ¼ w ik N is a DFT sinusoid.
᭝ The DFT preserves the units of x(i).That is, if x(i) has engineering units EU, then the units of X(k) are also EU.This is clear from Definition 20, where the exponential function is dimensionless.
If the discrete samples x(i) and corresponding sequence of DFTs X(k) are assembled into the configuration vectors x N ¼ ðxð1Þ; xð2Þ; …; xðNÞÞ T X N ¼ ðXð1Þ; Xð2Þ; …; XðNÞÞ T ) then the DFT of x N can be represented in the matrix-vector form where E N follows from Eq. ( 26).Thus, the DFT computation is simply a multiplication of the configuration vector x N with the DFT sinusoid matrix W N and constant c f ¼ 1. Equation (27) has exactly the same form as Eq.(22), where the generating elements of a circulant matrix are analogous to the sequence of signals x(i) and the resulting eigenvalues are analogous to the sequence of DFTs X(k).
Example 21.Consider the configuration vector of samples x 4 ¼ ð0; 1; 2; 3Þ T .The DFT of x 4 follows from Eq. ( 27) with N ¼ 4 and is given by which is the third element of X 4 , as expected.Results for k ¼ 1, 2, 4 follow similarly.
If X N is known, the N-vector x N is recovered from Eq. ( 27) by multiplying from the left by E H N and invoking Theorem 6.Then and it follows that is a matrix-vector representation of the IDFT of X N .That is, the IDFT computation is simply a matrix multiplication of X N with the Hermitian of W N and constant c i ¼ 1=N such that c f c i ¼ 1=N.
In light of Corollary 7, Eq. ( 28) can be written as Expanding each row yields . . .
which provides an alternative representation of the IDFT.DEFINITION 21 (IDFT).Let the sequence X(k) be defined according to Definition 20.Then for each i ¼ 1; 2; …; N, the IDFT of X(k) is the sequence

040803
where w N is the primitive Nth root of unity and S i ðkÞ ¼ w ik N is a DFT sinusoid.
᭝ The IDFT preserves the units of X(k).If X(k) has engineering units EU, then the units of x(i) are also EU.This is clear from Definition 21, where the exponential function is dimensionless, as is the number N.
There are other suitable definitions of the DFT/IDFT pair.For example, the signs of the exponents are irrelevant as long as they are opposite in the DFT and IDFT.To see this, suppose that the sign of the exponential in Definition 20 is negative instead of positive.Then each entry S k (r) in W N is replaced with S k ðrÞ to produce W N (see Corollary 7), and the DFT is instead given by The corresponding IDFT takes the form which follows in the same way as Eq. ( 28) by replacing E N with E N and invoking Corollary 9. Similarly, the multiplicative constants that define the DFT/IDFT pair are arbitrary as long as the product of the constants is equal to 1/N.This is important for applications involving a transformation of data to the frequency domain for analysis and then back to the time domain for results.However, the multiplicative constant is irrelevant if the analysis goal is only to identify periodicities in a data set (e.g., frequencies that correspond to amplification of a structural response).Thus, Definitions 20 and 21 form a representation of the DFT/ IDFT pair, where Eqs. ( 27) and ( 28) are the corresponding matrix-vector forms.The DFT and IDFT pair can be written in the general form (31a) where c i and c f are such that c f c i ¼ 1=N but are otherwise arbitrary.Common choices are fc f ; c i g ¼ f1= ffiffiffi ffi N p ; 1= ffiffiffi ffi N p g or fc f ; c i g ¼ f1; 1=Ng, where the multiplicative constants in the matrix-vector formulation (i.e., Eqs. ( 27) and (28) or Eqs. ( 29) and ( 30)) are adjusted accordingly.
If the generally complex sequences x(i) and X(k) are dissected into their real and imaginary parts, the DFT/IDFT pair is represented by where the subscripts R and I denote the real and imaginary parts and each sequence pair ðx R ðiÞ; x I ðiÞÞ and ðX R ðkÞ; X I ðkÞÞ is real for i; k ¼ 1; 2; …; N.An alternative representation of the DFT and IDFT is obtained by substituting Eq. (32) in Definition 20, writing the complex exponential in terms of sines and cosines, and collecting the real and imaginary terms.The result is which together form the DFT of x(i) according to Eq. (32b).A similar relationship for the IDFT is obtained by substituting Eq. ( 32) into Definition 21 and collecting the real and imaginary terms.Then together form the IDFT of X(k) according to Eq. (32a).The DFT/ IDFT pair defined by Eqs. ( 32)-( 34) are equivalent to Definitions 20 and 21.

2.8
The Circulant Eigenvalue Problem.The eigenvalues and eigenvectors of circulants and block circulants with circulant blocks have thus far been inferred from Theorem 9 (circulant matrix) and Corollary 22 (block circulant matrix with circulant blocks), where these matrices are fully diagonalized.Determination of the eigenvalues and eigenvectors of a block circulant matrix with generally noncirculant and nonsymmetric blocks is systematically described here, which reinforces the results already obtained for the special cases of C 2 C N and C 2 BC M;N with circulant generating matrices.The standard cEVP is discussed in Sec.2.8.1, where an eigensolution is obtained for a general matrix C 2 BC M;N .Section 2.8.2 generalizes these results to handle (M,K) systems with system matrices contained in BC M;N , which arise naturally in vibration studies of rotationally periodic structures.The structure of the eigenvalues and eigenvectors is discussed in Sec.2.8.3, which builds upon the relationship to the DFT described in Sec.2.7.Section 2.8.4 introduces fundamental orthogonality conditions for the special case of block circulants with symmetric generating matrices.
where 0 NM is a NM Â 1 vector of zeros and I NM is the identity matrix of dimension NM.The problem can be simplified considerably by exploiting the block circulant nature of C. To this end, partition q ¼ q 1 ; q 2 ; …; q N ð Þ T into M Â 1 vectors q i ði ¼ 1; …; NÞ and introduce the change of coordinates where I M is the M Â M identity matrix (same dimension as the generating matrices of C) and u ¼ u 1 ; u 2 ; …; u N ð Þ T is composed of N Mvectors u i .Substituting Eq. ( 36) into Eq.( 35) and multiplying from the left by the unitary matrix which follows from Theorem 10, and from in light of Corollary 11.Thus, the single eigenvalue problem defined by Eq. ( 35) is decomposed into the N reduced-order standard EVPs where 0 M is a M Â 1 vector of zeros and each K i is defined by Eq. ( 23) in terms of the generating matrices of C. Because the transformation defined by Eq. ( 36) is unitary, it preserves the eigenvalues of C. Thus, the NM eigenvalues of C are the eigenvalues of the N M Â M matrices K i , which follow from the characteristic polynomials denotes the pth eigenvalue of the ith matrix K i for p ¼ 1; 2; …; M, then the attendant M Â 1 eigenvector u ðpÞ i is obtained from Eq. ( 37).The corresponding NM Â 1 eigenvector of C follows from Eq. ( 36) and is given by The eigenvectors e 1 and e 3 are real and correspond to the distinct eigenvalues k 1 ¼ 2 and k 3 ¼ 6.The eigenvectors e 2 and e 4 are complex conjugates according to Eq. ( 17) and correspond to the repeated eigenvalue k The eigenvalues of C are obtained from the reduced-order matrices K i ði ¼ 1; 2; 3; 4Þ.For example, the eigenvalues of K 3 follow from the characteristic polynomial and are given by aðK 3 Þ ¼ f3; 5g.Similarly, aðK 1 Þ ¼ fÀ1; 1g and are distinct sets of eigenvalues and aðK 2 Þ ¼ aðK 4 Þ are repeated sets.Section 2.8.3 discusses the symmetry characteristics and multiplicities of these groups of eigenvalues.The eigenvectors of C follow from the eigenvectors of each K i according to Eq. (39).For example, for the eigenvector u ð1Þ 3 ¼ ð1; 1Þ T of K 3 with eigenvalue k ð1Þ 3 ¼ 3, the corresponding eigenvector of C is The entire set of NM ¼ 8 reduced-order eigenvectors and eigenvalues is given by u The eigenvectors are made orthonormal according to Eqs. ( 40) and ( 41) by setting each a The determinant of any matrix is the product of its eigenvalues, which yields the following result.
COROLLARY 24.Let C 2 BC M;N be a block circulant matrix.Then the determinant of C is given by where k i is the ith eigenvalue of C. ٗ The determinant of C 2 C N follows from Corollary 24 by setting M ¼ 1.In light of Eq. ( 21), for an ordinary circulant matrix with generating elements c 1 ; c 2 ; …; c N .For the special case of C 2 SC N , the ith eigenvalue k i in Corollary 24 can be replaced by the result given in Corollary 17.Similarly, k i can be replaced with Eq. ( 25) if C 2 BC M;N has circulant generating matrices.Example 26.Consider C ¼ circðA; B; 0; BÞ 2 BC 2;4 from Example 25, where it is shown that the eigenvalues are given by the set aðCÞ ¼ fÀ1; 1; 1; 3; 3; 5; 1; 3g.The determinant of C follows from Corollary (24) and is given by which is simply the product of the eigenvalues of C.

Generalized Circulant Eigenvalue
Problem.The generalized cEVP involves determination of the scalar eigenvalues k and NM Â 1 eigenvectors q of a cyclic system (M,K) such that where are block circulant matrices contained in BC M;N .This type of problem arises in the study of cyclic vibratory mechanical systems composed of N sectors with M DOFs per sector, where M and K are NM Â NM block circulant mass and stiffness matrices.Several examples are discussed in Sec. 3, where the eigenvectors q are the system mode shapes and the eigenvalues k correspond to the system natural frequencies.The generalized cEVP defined by Eq. ( 44) is handled in the same way as the standard cEVP of Sec.2.8.1.To this end, partition q ¼ q 1 ; q 2 ; …; q N ð Þ T into M Â 1 vectors q i ði ¼ 1; 2; …; NÞ, transform to a new set of coordinates u ¼ u 1 ; u 2 ; …; u N ð Þ T by substituting Eq. ( 36) into Eq.( 44), and left-multiply the result by the unitary matrix E H N I M .It follows that where the decomposed follow from Theorem 10.Equation ( 45) represents a set of N reduced-order generalized EVPs which are analogous to the reduced-order standard EVPs defined by Eq. (37).The eigenvalues are determined from the characteristic polynomials detð e and are denoted by k ¼ k ðpÞ i for p ¼ 1; 2; …; M. Each k ðpÞ i also satisfies Eq. ( 44) because the transformation to new coordinates via Eq.( 36) is unitary, and hence preserves the system eigenvalues.The reduced-order eigenvectors u i ¼ u ðpÞ i are obtained from Eq. ( 47) for each k ðpÞ i .The eigenvectors of the full system (M,K) are given by q ðpÞ i ¼ e i u ðpÞ i , which is the same as Eq. ( 39).The relationships defined by Eqs. ( 40) and ( 41) also hold for the generalized cEVP to make the eigenvectors orthonormal.
If M and K are ordinary circulants (i.e., M ¼ 1), then the system eigenvectors reduce to i ¼ e i 1 ¼ e i , which is the same as Eq. ( 42) and shows that e 1 ; e 2 ; …; e N are the orthonormal eigenvectors of all generalized eigensystems defined by M; K 2 C N .If M 1 ; M 2 ; …; M N and K 1 ; K 2 ; …; K N are the generating elements of M and K, respectively, the corresponding eigenvalues are which follows from Eq. ( 45) by replacing each e M k with M k and e K k with K k .For the special case of M ¼ I N with generating elements 1; 0; …; 0 (i.e., for the standard cEVP), Eq. ( 49) reduces to the form shown in Theorem 9, as expected.
It is clear from this formulation that the same basic steps are followed to solve the standard and generalized cEVPs.For the standard cEVP, we say that k ðpÞ i and q ðpÞ i are the eigenvalues and eigenvectors of the matrix C. For the generalized cEVP, the eigenvalues and eigenvectors correspond to the system (M,K), not the individual matrices M and K.
Example 27.Consider the generalized cEVP defined by Eq. ( 44).Let the generating elements of M 2 C 2 be M 1 ¼ 2 and M 2 ¼ À1 and the generating elements of K 2 C 2 be K 1 ¼ 5 and Then the eigenvalues of the system (M, K) follow from Eq. ( 49) with N ¼ 2 and are given by The corresponding orthonormal eigenvectors are given by which follow from Definition 18 according to Eq. ( 42).

Eigenvalue and Eigenvector
Structure.For real-valued generating elements, such as those that arise in models of physical systems with cyclic symmetry, the eigenvalues of a circulant matrix C 2 C N are endowed with certain symmetry characteristics, and the same is true for the eigenvalues of systems M; K 2 C N .However, we do not require the circulants to be symmetric, as it is assumed in Corollary 17.We begin by systematically describing the eigenvalue structure of an ordinary circulant with real generating elements.The results are generalized by inspection to handle block circulants and (M,K) systems composed of circulant or block circulant matrices.
The eigenvalues k i of a circulant matrix C 2 C N with generating elements c 1 ; c 2 ; …; c N are given by Eq. ( 21), or the equivalent matrix-vector form in Eq. (22).It is convenient to re-index each c k and k i such that yðpÞj p¼0;1;2;…;NÀ1 ¼ c k j k¼1;2;3;…;N YðrÞj r¼0;1;2;…;NÀ1 ¼ k i j i¼1;2;3;…;N follow from the formulation of Eq. ( 33) in Sec.2.7 for a realvalued signal where, recall, the eigenvalue expression defined by Eq. ( 21) has exactly the same form as the DFT in Definition 20.
That is, Eq. ( 51) also represents the real and imaginary parts of the DFT of a real-valued signal, where the sequence of generating elements y(p) is analogous to a real signal and the eigenvalues Y(r) are analogous to its DFT.It is shown that the symmetry of the DFT about the so-called Nyquist component also exists for the eigenvalues of a circulant matrix with real generating elements.
As is customary in signal processing, we restrict the formulation to even N.The case of odd N also yields symmetric eigenvalues, but with multiplicity of the Nyquist component.This is handled by example in Sec. 3 (for instance, see Fig. 9).For even N, the zeroth eigenvalue Y(0) and "Nyquist" eigenvalue Y N=2 ð Þare always real because are such that the imaginary parts vanish. 2The remaining eigenvalues appear in complex conjugate pairs, as the following corollaries show.COROLLARY 25.Let y(p) be the real-valued generating elements of a circulant matrix for p ¼ 0; 1; 2…; N À 1 and Y(r) denote its eigenvalues according to Eq.where the property sinðÀhÞ ¼ À sinðhÞ is employed.It follows that which completes the right-hand side of the proof.To prove the left-hand side, consider and a similar expansion shows that Y I ðN À rÞ ¼ ÀY I ðrÞ.It follows that ٗ The magnitudes jYðrÞj and arguments ffYðrÞ of the eigenvalues Y(r) are listed in Table 3 for the special case of even N ¼ 8, where YðrÞ ¼ YðN À rÞ follows from complex conjugation of the lefthand equality given by Corollary 25.The zeroth eigenvalue Y(0) ¼ Y 0 is real and generally distinct, as is Y N=2 ð Þ¼Y N=2 for even N, but the remaining eigenvalues generally appear in complex conjugate pairs.It follows that the eigenvalue magnitudes (frequencies) are symmetric about the Nyquist component Y N=2 , as are the arguments (phase angles) but with the opposite sign for indices r > N=2.
The A similar formulation shows that the eigenvalues of real-valued (M,K) systems exhibit the same symmetry characteristics because the numerator and denominator of Eq. ( 49) have exactly the same form as Eq. ( 21).If the scalars M k and K k in Eq. ( 49) are reindexed and restricted to be real-valued, then Corollaries 25 and 26 are generalized accordingly.
For signal processing applications, the symmetry characteristics summarized in Table 3 have practical significance because the subset of ðN þ 2Þ=2 frequency-domain components 0 r N=2 is endowed with the same basic "information" contained in all N time-domain signal samples for 0 p N À 1.That is, when transformed to the frequency domain by the DFT process, any real-valued signal has a zero-frequency or dc component (r ¼ 0), distinct magnitude and phase information for 0 r N=2, and repeated magnitude and phase information for r > N=2.For circulant matrices that describe physical systems with cyclic symmetry, Y 0 and Y N=2 correspond to standing wave vibration modes and the remaining eigenvalues are associated with traveling wave modes.This is discussed in Sec. 3.
Equation (42) shows that e 1 ; e 2 ; …; e N are the eigenvectors of any circulant matrix C, and the same is true for (M,K) systems composed of circulant matrices.Each eigenvector e i is associated with an eigenvalue Y(r) (i.e., k i ) according to the indices defined by Eq. ( 50).For the special case of real generating elements, the eigenvalues summarized in Table 3 have exactly the same symmetry characteristics as the vectors e 1 ; e 2 ; …; e N , which are discussed in Sec.2.5.2.For example, if N is even, the real eigenvalues Y N/2 (i.e., k (N+2)/2 ) and Y(0) (i.e., k 1 ) correspond to the real eigenvectors e Nþ2=2 ð Þ and e 1 defined by Eqs. ( 18) and (19), respectively.The remaining eigenvalues appear in complex conjugate pairs and are associated with the complex conjugate eigenvectors according to Eq. (17).
If the scalar sequences y(p) and Y(r) for a circulant matrix are replaced by a sequence of matrices y(p) and Y(r), it is clear that the formulation given above also holds for Eq. ( 23), which defines the eigenvalues for block circulant matrices.In this case, the groups of eigenvalues associated with each K i are endowed with the symmetry properties given in Table 3.This is confirmed by Example 25 where, using the indexing scheme of that section, aðK 1 Þ ¼ fÀ1; 1g and aðK 3 Þ ¼ f3; 5g are distinct sets of eigenvalues and aðK 2 Þ ¼ aðK 4 Þ ¼ f1; 3g are repeated.

Eigenvector Orthogonality.
Here, we consider eigenvector orthogonality relationships for block circulant matrices C and systems (M,K) contained in BC M;N for the special case of symmetric generating matrices.However, we do not restrict either C or (M,K) to be symmetric or block symmetric.We require only that the generating matrices C 1 ; C 2 ; …; C N of C are symmetric for the standard cEVP, which guarantees that each K i is symmetric according to Corollary 19.Similarly, we require that the generating matrices M 1 ; M 2 ; …; M N and K 1 ; K 2 ; …; K N of (M,K) are symmetric for the generalized cEVP, which implies that ð e K i ; e M i Þ are symmetric.Symmetric generating matrices commonly arise in models of rotating flexible structures, including the ones considered in Sec. 3, where the sector models and intersector coupling are described by symmetric matrices.Then the usual orthogonality relationships hold for the reduced-order eigenvectors defined in Secs.orthogonal with respect to K i for the standard cEVP.Proofs can be found in any standard textbook on linear algebra [98,99].Together with generic orthogonality conditions on the basic circulant structure of a matrix C 2 BC M;N , this gives rise to an orthogonality condition on the eigenvector qðpÞ i ¼ e i ũðpÞ i with respect to C. Orthogonality of qðpÞ i with respect to the system M; K 2 BC M;N is handled similarly.
The results of this section, and the requirement of symmetric generating matrices, are meant to highlight how orthogonality of an eigenvector qðpÞ i ¼ e i ũðpÞ i is essentially dissected into the orthogonality of e i with respect to the circulant structure of C 2 BC M;N and orthogonality of ũðpÞ i with respect to the generating matrices (e.g., K i ), where the latter requires symmetric C 1 ; C 2 ; …; C N .It should be noted that none of the results in Sec. 2, aside from this section, require symmetric generating matrices.In particular, Theorems 9 and 10, upon which block reduction of the cEVPs in Secs.
corresponding to the symmetric matrix ! from Example 25, where the reduced-order eigenvectors are in orthonormal form.The corresponding reduced-order modal matrix is denoted by e U 3 ¼ ðũ It follows from Corollary 27 that the diagonal matrix Orthogonality of an eigenvector qðpÞ i ¼ e i ũðpÞ i with respect to C 2 BC M;N is essentially decomposed into orthogonality of e i with respect to the circulant structure of C and orthogonality of each ũðpÞ i with respect to the symmetric matrices K i .These individual orthogonality conditions are captured by Corollaries 18 and 27, which lead to the following fundamental result.COROLLARY 28.Let C 2 BC M;N be a block circulant matrix with symmetric generating matrices, e i be the ith column of the N Â N Fourier matrix E N , and ũðpÞ i be the pth reduced-order orthonormal eigenvector corresponding to K i defined by Eq. (23).Then is numerically equal to k 2 according to Eq. ( 54).However, the product Example 30.Consider C ¼ circðA; B; 0; BÞ 2 BC 2;4 from Examples 6, 18, 19, 20, and 25, where it is shown that ffiffi ffi 2 p ðÀ1; 1; 1; À1; À1; 1; 1; À1Þ T is an orthonormal eigenvector of C corresponding to the eigenvalue k ð2Þ 3 ¼ 5.Because the generating matrices A, B, 0, B are symmetric (see Example 6), Corollary 28 guarantees that each qðpÞ i ¼ e i ũðpÞ i is mutually orthogonal with respect to C. For example, it follows from Eq. ( 53) that 3 .However, the matrix product ffiffi ffi 2 p ðÀ1; 1; 1; À1; À1; 1; 1; À1Þ T ffiffi ffi 2 p ðÀ1; 1; 1; À1; À1; 1; 1; À1Þ T The orthogonality conditions used in Example 30 for the special case of an ordinary circulant does not require that the circulant structure is symmetric (i.e., C need not be contained in SC N ).The reader can verify that Eq. ( 54) also holds for nonsymmetric matrices C 2 C N by inspection of Example 16.
The fundamental orthogonality relationship given by Corollary 28 also holds for M; K 2 BC M;N systems with symmetric generating matrices, as the following corollary shows.
Corollary 29.Let M; K 2 BC M;N be block circulant with symmetric generating matrices, e i be the ith column of the N Â N Fourier matrix E N , and ũðpÞ i be the pth M Â 1 reduced-order orthonormal eigenvector corresponding to the ith system ð e M i ; e K i Þ defined by Eq. (46).Then ; for i; k ¼ 1; 2; …; N and p; q ¼ 1; 2; …; M, where k ðpÞ i is the eigenvalue associated with ũðpÞ i and d ik is the Kronecker delta.ٗ In practice, of course, the M Â 1 reduced-order eigenvectors ũðpÞ i are not known a priori.Instead, Theorem 10 is used in Sec. 3 to decouple the NM-DOF system equations into a set of N reduced-order M-DOF standard vibratory problems, from which the system eigenvalues (natural frequencies) and eigenvectors (normal modes) are extracted.

Example Applications
In this section we apply the results developed in Sec. 2 to vibration models of systems with cyclic symmetry.For each model considered, we begin by formulating the equations of motion (EOM) and then use the theory of circulants to diagonalize or block diagonalize the governing equations.This is achieved by a coordinate transformation that exploits the special relationship between circulant matrices and the Fourier matrix.The process also shows how external forces are projected on the resulting block diagonal EOM.The special case of traveling wave engine order excitation is also presented in some detail because it appears in many relevant applications of rotating machinery.
Three examples are presented.The first example (Sec.3.1) considers the structure of the EOM for a general cyclic system with N sectors, M DOFs per sector, and arbitrary excitation.It is shown how to block diagonalize the system equations via a modal transformation involving the Fourier matrix, which results in NM-DOF reduced-order vibratory systems.If engine order excitation is assumed (Sec.3.2), it is shown that the steady-state forced response of the NM-DOF system can be obtained from a single M-DOF harmonically forced, reduced-order system in modal space.The mathematical and physical details of engine order excitation are discussed, including its temporal and spatial duality.The second example (Sec.3.3) considers a cyclic system with one DOF per sector under engine order excitation.This system is fully diagonalized by the Fourier matrix.The example is presented in detail, showing the nature of the natural modes and frequencies, and the resonant response to excitation of various engine orders.The third example (Sec.3.4) has two DOFs per sector and demonstrates the block diagonalization process for a perfectly cyclic system with specified sector models, as opposed to the general sector models in the first example.In each sector, one DOF is due to flexure, and thus has a constant frequency, while the other DOF is a centrifugally driven pendulum whose frequency is proportional to the rotor speed.The coupling between these DOFs leads to some interesting behavior in both the free and forced vibration of the system, which is discussed in Refs.[21,22] and [92][93][94][95].More importantly, this example shows the process of handling multiple DOFs per sector, which easily extends from two to M DOFs per sector using the theory presented in Sec. 2.

General Cyclic System
3.1.1Equations of Motion.Consider the general cyclic vibratory system shown schematically in Fig. 4, which consists of N sectors with coupling (elastic and damping) to adjacent neighbors.The topology diagram only indicates nearest-neighbor coupling, but more general coupling is admissible as long as the cyclic symmetry is preserved.If there are M DOFs per sector, each M Â 1 vector q i describes the dynamics of the ith sector for i 2 f1; 2; …; Ng N .Then the linear EOM takes the form where q ¼ ðq 1 ; q 2 ; …; q N Þ T is a NM Â 1 configuration vector, the system matrices are block circulant with M Â M blocks, and overdots denote differentiation with respect to time.If the M Â 1 vector f i denotes the component of forcing on the ith sector, then b f ¼ ðf 1 ; f 2 ; …; f N Þ T is a NM Â 1 general forcing vector.The NM Â NM system mass, damping, and stiffness matrices are of the form where the generating matrices M i , C i , and K i depend on the M Â M sector mass, damping, and stiffness matrices and the intersector coupling (stiffness and damping).Equation ( 55) is a general model for any linear, lumped-parameter, conservative, nongyroscipic, cyclic vibratory system with N sectors and M DOFs per sector.For example, a linearized lumped-parameter model of a bladed disk assembly under engine order excitation is captured by Eq. ( 55), where N is the number of blades, M is the number of DOFs per blade, f has the special properties described in Sec.3.2, and the system matrices depend on the structural details of each blade (i.e., sector) and its connectivity to adjacent blades and rotating hub.

Modal Transformation.
Of course, one can apply standard techniques [99] to investigate the free and forced response of the model given by Eq. ( 55).However, this requires solving an NM Â NM eigenvalue problem to determine the modal properties, or inversion of an NM Â NM impedance matrix to determine the response to harmonic excitation.Such an approach may be prohibitive or computationally expensive for practical models with Fig. 4 Topology diagram of a general cyclic system many sectors and many DOFs per sector, nor does it highlight or exploit the underlying features of the cyclic system.It is precisely these special features that are brought to light by the properties of the circulants or block circulants that describe the system.Solving for the system response is significantly facilitated by a modal transformation that exploits the cyclic symmetry among the N sectors.Specifically, it is straightforward to block decouple the EOM into a set of N systems, each with M DOFs.To this end, we introduce the change of coordinates where u ¼ ðu 1 ; u 2 ; …; u N Þ T is a NM Â 1 vector of modal coordinates.Each u i is M Â 1 and describes the sector dynamics in modal space, where the EOM are block decoupled (as described below).In this formulation the physical coordinates q are expressed in terms of the modal coordinates u i and the Fourier elements of length N, which accounts for the overall cyclic nature of the EOM.Substituting Eq. ( 57) into Eq.( 55) and multiplying from the left by The block diagonal structure on the left-hand side of Eq. ( 58) follows from Theorem 10, where the M Â M modal mass, damping, and stiffness matrices associated with the pth mode follow from Eq. ( 23) and are given by The forcing terms on the right-hand side of Eq. ( 58) follow from the decomposition and make no assumptions on the nature of the applied forcing.Thus, the solution to the NM-DOF matrix EOM given by Eq. ( 58) reduces to solving NM-DOF uncoupled systems for the modal solutions u p .These N reduced-order EOMs of order M offer a substantial computational savings compared to the full NM-DOF system defined by Eq. ( 55).The reduced equations can be solved directly or with standard modal analysis by solving a set of N generalized eigenvalue problems, each of order M, from which one can form the global eigenvectors of the system.Specifically, solving the N eigenvalue problems associated with Eq. ( 61) yields a set of normalized M Â M modal matrices U p , which define a change to reduced modal coordinates via u p ¼ U p s p .The s p are modes that define the behavior of the internal sector dynamics and account for the manner in which these are coupled to each other in a given Fourier mode of the overall system.The EOM for the s i form a set of M uncoupled equations, and these are related to the physical coordinates of the original system by q ¼ ðE N I M Þu where u ¼ ðU 1 s 1 ; U 2 s 2 ; …; U N s N Þ T .This process is general, and it allows one to decouple the full EOM in two steps, one which accounts for the global cyclic nature of the system, and the other which accounts for the details of the sector model.Equation ( 61) can be simplified even further if the system is subjected to the so-called engine order excitation.The mathematics and physics of this type of excitation are developed in the next section, which closes with a treatment of the general cyclic system governed by Eq. ( 61) under engine order excitation.

Engine Order Excitation.
Traveling wave engine order excitation arises in rotating machinery and is a primary source of forced vibration response in bladed disk assemblies [11,111].A mathematical model of this common form of excitation is developed in Sec.3.2.1 and its traveling wave characteristics are described in Sec.3.2.2.While this material is known, the discussion that follows is unique because it provides physical insights and a systematic explanation of the temporal and spatial duality of engine order excitation.The general cyclic system of Sec.3.1 is reconsidered in Sec.3.2.3under engine order excitation, where it is shown that the steady-state forced response of the NM-DOF system reduces to that of a single sector in modal space.

Mathematical Model.
Ideally, the steady axial gas pressure in a jet engine might vary with radius but is otherwise uniform in the circumferential direction, thus resulting in an identical force field on each blade in a particular fan, compressor, or turbine within the engine.In practice, however, flow entering an engine inlet invariably meets static obstructions, such as struts, stator vanes, etc., in addition to rotating bladed disk assemblies in its path to the exhaust.Even in steady operation, therefore, the flow slightly upstream of these bladed assemblies is nonuniform in pressure, temperature, and so on.This results in a static pressure (effective force) field on the blades that vary circumferentially, a notional example of which is shown in Fig. 5.
Consider, for example, an engine in steady operation with n evenly spaced struts slightly upstream (or downstream) of a bladed assembly.As explained in Ref. [3], these obstructions produce a circumferential variation upon the mean axial gas pressure that is essentially proportional to cos nh, where h is an angular position.Thus, a blade rotating through this static pressure field experiences a force proportional to cos nXt, where X is the constant angular speed of the bladed disk assembly and t is time.An adjacent blade experiences the same force, but at a constant fraction of time later.This type of excitation is defined as engine order (e.o.) excitation and n is said to be the order of the excitation.
To be more precise, the axial gas pressure of a steady flow through a jet engine may be described by the function pðhÞ ¼ pðh þ 2pÞ, where h is an angular coordinate measured 040803-28 / Vol.66, JULY 2014 Transactions of the ASME relative to a fixed origin on the machine.That is, the pressure field is rotationally periodic and can therefore be expanded in a Fourier series with terms of the form p o cos nh.If the angular position of the ith blade relative to the same origin is defined by where N is the total number of blades and N ¼ f1; 2; …; Ng is the set of blade, or sector numbers, it follows that the total effective force exerted on blade i due to the nth harmonic of the pressure field pðhÞ is captured by Upon complexifying, is a model for the nth predominant component of the excitation.Eq. ( 63) has period T ¼ 2p=nX, strength F, and is said to have angular speed X.The so-called interblade phase angle is defined by where n 2 Z þ and u i is the angle subtended from blade 1 to blade i and is defined by Eq. ( 16).Equation ( 63) is defined as nth engine order, or traveling wave excitation.The traveling wave characteristics of this type of excitation are considered next.

Traveling Wave
Characteristics.The function defined by Eq. ( 63) is continuous in time and discretized in space via the index i.This gives rise to two interpretations of engine order excitation relative to the rotating hub, one discrete and the other continuous.These can be visualized in Fig. 6, which shows a dissection of the excitation amplitudes along time and sector axes.In the first and usual sense, Eq. ( 63) is a discrete temporal variation of the dynamic loading applied to individual blades.That is, under an engine order n excitation, each sector is harmonically forced with strength F and frequency nX, but with a fixed phase difference relative to its nearest neighbors.Physically, one can think of this as placing N different observers at the discrete sectors and having the ith observer record the excitation strength applied to sector i as a function of time.Their recorded time traces would resemble those shown in Fig. 6(a).In the second and more general sense, Eq. ( 63) can be viewed as a continuous spatial variation of the excitation strength relative to the rotating hub (along the sector axis) that evolves with increasing time, i.e., it is a propagating waveform, or traveling wave.If a single observer was placed on the rotating hub and recorded the strength of this traveling wave as a function of i (taken here to be continuous), it would resemble the curve shown in Fig. 6(b).In this context, the instantaneous loading applied to individual blades is obtained by essentially "sampling" the continuous traveling wave at each sector i 2 N and, as time evolves, these sampled points define N time profiles of the force amplitudes, which is equivalent to the discrete temporal interpretation described above.However, the latter interpretation illuminates some important traveling wave characteristics of the engine order excitation that are otherwise difficult to explain, and in what follows these are systematically described.
To explain the traveling wave mathematically, it is convenient to define the function which is a harmonic waveform with wavelength 2p=u k .Then for i 2 N and noting that / i ¼ u nþ1 ði À 1Þ, Eq. ( 63) can be written in real form as which is a harmonic function with a wavelength of 2p=u nþ1 ¼ N=n (u nþ1 is the wave number) and angular frequency nX.Equation (66) shows that engine order excitation is a TW in the negative i-direction (descending blade number) with speed C ¼ nX=u nþ1 ¼ NX=2p, measured in sectors per second.An example plot of this continuous BTW is shown in Fig. 6(b) and, as described above, the applied loads can be obtained from this figure by continuously "sampling" the waveform at the discrete sector numbers as time evolves.Then the engine order excitation applied to the individual blades consists of a wave composed of these N discrete points, examples of which are shown in Figs.7(a)-7(d).Interestingly, this gives rise to discrete SW or even FTW applied dynamic loads (depending on the value of n relative to N), even though Eq. ( 66) is strictly a BTW relative to the rotating hub.These additional possibilities arise due to aliasing of the "sampled points" just as it occurs in elementary signal processing theory [112,113].Before characterizing the traveling and standing waveforms, it is shown that one need only consider engine orders n 2 N .The traveling wave nature of the discrete applied loads (i.e., SW, BTW, or FTW) depends only on the value of n relative to N. To see this, let and assume n ¼ n þ mN for some integer m.Then That is, if n ¼ n corresponds to a SW, BTW, or FTW, then so does n þ mN for any m 2 Z þ .In this sense, the traveling wave nature of the applied dynamic loads is seen to alias relative to N.These features are characterized for engine orders where it is understood that the results apply to any n > N simply by taking n modulo N, as appropriate.The subsets N O;E BTW , N O;E FTW , and N O;E SW are defined in Table 4 and discussed below.For the special case of n ¼ N the rotating blades become entrained with the excitation because / ðNÞ i ¼ 2pnði À 1Þ with Fig. 5 The axial gas pressure pðhÞ: ideal and (notional) actual conditions i; n 2 Z þ , and hence each is forced with the same strength and phase.As illustrated in Fig. 7(d), this is a SW excitation where each blade is harmonically forced according to F i ðtÞ ¼ F cos nXt.Entrainment also occurs when n ¼ N/2 if N is even, in which case / ðN=2Þ i ¼ pnði À 1Þ, where (i À 1) is odd (resp.even) for even (resp.odd) sector numbers i 2 N .Accordingly, all blades with odd sector numbers are driven by F i ðtÞ ¼ F cos nXt, as are the blades with even sector numbers, but with a 180-deg phase shift.
As shown in Fig. 7 4.These sets can be visualized in Figs.7(i) and 7(ii) for odd and even N, respectively.
The manner in which cyclic systems respond to this type of excitation is considered in the examples that follow.In Sec.3.2.3,we prove the most important result related to the forced response of cyclic systems, namely, that in the case of perfect symmetry each engine order excites only a single mode.This is clear mathematically and will be explored with more physical insight in the examples presented in Secs.3.3 and 3.4.

Forced Response of a General Cyclic System Under
Engine Order Excitation.The general cyclic system governed by Eq. ( 55) is reconsidered here under engine order excitation.Using the notation of Sec.3.1, a model for the nth engine order excitation is where f is a constant M Â M vector of sector force amplitudes, t is time, u i is the angle subtended from sector 1 to sector i and is defined by Eq. ( 16), n is the order of the excitation, nX is the angular frequency of the excitation, and X is the angular speed of the system relative to the excitation.Under this type of excitation, the system forcing vector becomes  f 0 ¼ ðe jnu 1 ; e jnu 2 ; …; e jnu N Þ T ¼ ðe j/ 1 ; e j/ 2 ; …; e j/ N Þ T (69) is a vector of constant intersector phase angles / i ¼ nu i .Thus, the pth modal forcing term on the right-hand side of Eq. ( 61) reduces to which is a direct product of the scalar product e H p f 0 with the nth order harmonic excitation fe jnXt imparted to each sector.This dissection of the modal forcing term highlights an orthogonality condition that is shared by all cyclic systems under engine order excitation.In particular, e H p f 0 ¼ which follows from Lemma 2 and shows that the force vector f 0 is mutually orthogonal to all but one of the modal vectors e p .This result is obtained more directly from Corollary 10 by observing that which is the same orthogonality condition given by Eq. ( 71).Thus, for p 2 N , the only excited mode is for an engine order n 2 Z þ excitation.Then Eq. ( 70) becomes and the forced response of the NM-DOF matrix EOM given by Eq. ( 58) reduces to solving a single, M-DOF system in modal space.Assuming harmonic motion, the steady-state modal response is given by where Table 4 Sets of engine orders nðmod NÞ ‰ N corresponding to BTW, FTW, SW dynamics loads applied to the blades for odd and even N.These can be visualized in Figs.7(i) and 7(ii).
is the ðn þ 1Þth modal impedance matrix.All other steady-state modal responses are zero because only mode p ¼ n þ 1 is excited.
In light of the decomposition shown in Eq. ( 39), the forced response in physical coordinates follows from Eq. ( 57) and is given by q ss ðtÞ ¼ e nþ1 u ss nþ1 ðtÞ (77) where e nþ1 is the ðn þ 1Þth column of the Fourier matrix and u ss nþ1 ðtÞ is given by Eq. ( 76).Expanding Eq. ( 77) into its sector components yields q ss 1 ðtÞ q ss 2 ðtÞ . . .q ss i ðtÞ . . .
where / i ¼ 2pðn=NÞði À 1Þ.The response of each sector is identical but simply shifted in time by a constant phase relative to its nearest neighbors.

Cyclic System With
One DOF Per Sector.This example considers the simplest prototypical model for vibrations of a bladed disk with only one DOF per sector, nearest-neighbor coupling, and perfect symmetry.Despite the simplicity of the model, we show that the resonant response can be quite complicated when the system is subjected to engine order excitation with multiple orders.We begin by formulating the EOM, and then consider a direct (traditional) approach to deriving the response to traveling wave excitation.The forced response is then derived using the modal analysis based on the Fourier matrix, which decouples the EOM.The nature of the natural frequencies and modes is considered next, which sets the stage for examining the resonance behavior of the system when subjected to engine order excitation.This provides a quite general view of the forced response of these models.

Equations of Motion
. The undamped cyclic system to be considered is shown in Fig. 8 in dimensionless form.It consists of a cyclic chain of N identical and identically coupled single-DOF oscillators with unit mass, the dynamics of which are captured by the dimensionless transverse displacements q i for i 2 N .The oscillators are uniformly attached around the circumference of a stationary rigid hub via linear elastic elements with unit stiffness and unit effective length.Adjacent masses are elastically coupled via linear springs, each with nondimensional stiffness .It is assumed that the elastic elements are unstressed when the oscillators are in a purely radial configuration, that is, when q i ¼ 0 for each i 2 N .An individual oscillator, together with the forward-nearest-neighbor elastic coupling, forms one fundamental sector and there are N such sectors in the overall system.The oscillator chain has cyclic boundary conditions such that q 0 ¼ q N and q Nþ1 ¼ q 1 .Finally, the system is subjected to engine order excitation (Sec.3.2) according to f i ðtÞ ¼ fe j/ i e jnrs ; i 2 N (80 where f is the strength of the excitation, the interphase blade angle / i is defined by Eq. ( 64), n 2 Z þ is the excitation order, r is the angular speed, and s is time (all dimensionless).The linear dynamics of the ith sector are obtained using Newton's laws and are governed by where overdots denote differentiation with respect to dimensionless time s.In Eq. ( 81), the q i61 terms arise from the left-nearestneighbor (i À 1) and right-nearest-neighbor (i þ 1) elastic coupling.By stacking the N coordinates q i into the configuration vector q ¼ ðq 1 ; q 2 ; …; q N Þ T , the governing EOM for the overall N-DOF system takes the form where f 11 ¼ ðfe j/ 1 ; fe j/ 2 ; …; fe j/ N Þ T is the system forcing vector, which accounts for the constant phase difference in the dynamic loading from one sector to the next.The N Â N matrix captures the nondimensional stiffness of each sector relative to the hub (additive unity along its diagonal) and the intersector coupling (t 2 along the super-diagonal and subdiagonal).The elements Àt 2 appearing in the (1, N) and (N, 1) positions of K 11 are due to the cyclic boundary conditions q 0 ¼ q N and q Nþ1 ¼ q 1 .In the absence of these cyclic coupling terms, the system represents a finite chain of N oscillators.Thus, in addition to being symmetric, Eq. ( 83) is also a circulant and can be written as where 1 þ 2t 2 ; Àt 2 ; 0; …; 0; Àt 2 are the N generating elements.
In the absence of coupling (that is, if t ¼ 0) K 11 is diagonal and Eq. ( 82) represents a decoupled set of N harmonically forced, single-DOF oscillators.
The forced response of Eq. ( 82) is considered next with emphasis on a modal analysis whereby the fully coupled system (that is, one in which t 6 ¼ 0) is reduced to a set of N single-DOF oscillators, only one of which is harmonically excited.The approach taken here, and a generalization in which each sector has multiple DOFs, is applied to the linear system in Sec.3.4 to block decouple the system matrices as it is done in Sec.3.1.

Forced
Response.The steady-state forced response of Eq. ( 82) can be obtained using standard techniques [99] and, for nonresonant forcing, is given by where I is the N Â N identity matrix.However, this requires inversion of the impedance matrix K 11 À n 2 r 2 I, which is computationally expensive for a large number of sectors, and it offers little insight into the basic vibration characteristics.In what follows, a transformation based on the cyclic symmetry of the system is exploited to fully decouple the single N-DOF system to a set of N single-DOF oscillators from which the steady-state response is easily obtained.The procedure is similar to the usual modal analysis from elementary vibration theory.However, a key difference is that the transformation matrix (and hence the system mode shapes) is known a priori and, because the transformation is unitary (thus preserving the system eigenvalues), the natural frequencies are obtained after the transformation is carried out.Moreover, due to orthogonality conditions between the normal modes and forcing vector, the steady-state response of the overall system reduces to finding the forced response of a single harmonically forced, single-DOF oscillator in modal space, which offers a clear advantage over the direct computation of Eq. ( 85).
As described in Sec.3.2.2,engine order excitation can be regarded as traveling wave dynamic loading.It is therefore reasonable to expect steady-state solutions of the same type.We begin with a simple way to show the existence of such a response, and then systematically describe it using the modal analysis techniques described in Secs.3.1.2and 3.2.3.
Existence of a Traveling Wave Response.It is natural to search for steady-state solutions of the form q ss i ðsÞ ¼ Ae j/ i e jnrs ; i 2 N which has the same traveling wave characteristics as the engine order excitation described in Sec.3.2.2.Equation ( 86) assumes that each sector responds with the same amplitude A, but with a constant phase difference relative to its nearest neighbors, and together all N such solutions form a traveling wave response among the sectors.By mapping this trial solution into Eq.( 81) and dividing through by the common term e j/ i e jnrs , it follows that where the identity / i61 À / i ¼ 6u nþ1 is employed.Solving for the amplitude A yields from which it follows that is one of the N natural frequencies of the coupled system corresponding to mode p ¼ n mod N þ 1. Equation (88) shows that mode n þ 1 is excited, but the reason is not clear from this approach.A modal analysis that considers the fully coupled system is required to systematically describe the response characteristics of the cyclic system under engine order excitation.
Modal Analysis.Theorem 9 guarantees that circulant matrices, such as the stiffness matrix defined by Eq. ( 84), are diagonalizable via a unitary transformation involving the Fourier matrix, and in what follows this property is exploited to fully decouple the matrix EOM defined by Eq. (82).To this end, the change of coordinates qðsÞ ¼ EuðsÞ or q i ðsÞ ¼ e T i uðsÞ; i 2 N is introduced, where E is the N Â N complex Fourier matrix (Definition 17), e i is its ith column (Definition 18), and u ¼ ðu 1 ; u 2 ; …; u N Þ T is a vector of modal, or cyclic coordinates.Substituting Eq. ( 89) in Eq. ( 82) and multiplying from the left by E H yields where  60) by replacing the identity matrix I M with unity and the vector b f with f 11 .Equation ( 89) is a unitary (similarity) transformation and hence the system natural frequencies are preserved, which is guaranteed by Theorem 1.For each p 2 N , the dimensionless natural frequencies follow from Eq. ( 21) and are given by where w N is the primitive Nth root of unity and the identity w 91) is a decoupled set of N single-DOF harmonically forced modal oscillators of the form Thus, the single N-DOF system given by Eq. ( 82) is transformed to a system of N decoupled single-DOF systems defined by Eq. ( 93).The steady-state, nonresonant modal response of the pth decoupled system follows from Eq. ( 93) using standard techniques [99].Assuming harmonic motion, the solution is where the identity ðp À 1Þu i ¼ u p ði À 1Þ is employed.Equation (95) shows that there are N possible resonances, depending on the details of the modal forcing terms e H p f 11 .However, only a single mode survives under an engine order excitation of order n, which is clear from the orthogonality condition described in Sec.3.2.3.Noting that f 11 ¼ f f 0 ¼ ffiffiffi ffi N p f e nþ1 , it follows from Eq. ( 72) that which shows that the force vector f 11 is mutually orthogonal to all but one of the modal vectors e p .That is, for p 2 N , the only excited mode is for an engine order n 2 Z þ excitation.Thus, Eq. ( 95) reduces to nþ1 À ðnrÞ 2 e j/ i e jnrs ; i 2 N where the identity ði À 1Þu nþ1 ¼ / i is employed and from Eq. ( 92).Equation ( 98) is the same result as that obtained from Eq. ( 88).Indeed, the process described here is significantly more laborious than the direct approach, but many general features can be gleaned from the analysis.The eigenfrequency characteristics (Sec.3.3.3),normal modes of vibration (Sec.3.3.4),and resonance structure (Sec.3.3.5)are systematically described based on the modal decomposition results of this section.

Eigenfrequency Characteristics.
The dimensionless natural frequencies follow from Eq. ( 92) and are given by which clearly exhibit the effect of cyclic coupling.For the special case of t ¼ 0, the sectors are dynamically isolated and each has the same natural frequency x p ¼ 1.There are repeated natural frequencies for nonzero coupling ðt 6 ¼ 0Þ, a degeneracy that is due to the circulant structure of K.This is captured by the cyclic term which is obtained by projecting the powers of the Nth roots of unity onto the real axis (see Fig. 2).Multiplicity of the eigenfrequencies can also be visualized in Fig. 9, which shows the dimensionless natural frequencies in terms of the number of nodal diameters (n.d.) in their attendant mode shapes versus: the mode number p for the special case of N ¼ 10 sectors; the wave type (i.e., BTW, FTW, or SW); the number of n.d.; and the sector number i. 3 Results are shown for weak coupling (WC), strong coupling (SC), odd N (Fig. 9(a)), and even N (Fig. 9(b)).These cyclic features are described in terms of mode numbers where each subset is defined in Table 5.A description of the BTW, FTW, and SW designations of these sets is deferred to Sec 3.3.4.The natural frequency corresponding to mode p ¼ 1 2 P O;E SW (zero harmonic of Eq. ( 100)) is distinct, but the remaining natural frequencies appear in repeated pairs, except for the case of even N, in which case the p ¼ ðN þ 2Þ=2 2 P E SW frequency (N/2 harmonic) is also distinct.There are (N À 1)/2 such pairs if N is odd, and these correspond to mode numbers in P O BTW and P O FTW , respectively.For even N there are (N À 2)/2 repeated natural frequencies corresponding to mode numbers in P E BTW and P E FTW .Finally, if k 2 P O;E BTW then the mode number of the corresponding repeated eigenfrequency is N þ 2 À k 2 P O;E FTW .The normal modes

3
A mode shape nodal diameter refers to a line of zero sector responses across which adjacent sectors respond out of phase.For example, in Fig. 11 82) can be decoupled via a unitary transformation involving the Fourier matrix E ¼ ðe 1 ; e 2 ; …; e N Þ.As a consequence, e p is the pth normal mode of vibration corresponding to the natural frequency x p .In what follows these mode shapes are characterized by investigating the free response of the system, and it is shown that they are of the SW, BTW, or FTW variety.
The free response of the system in its pth mode of vibration can be described by q ðpÞ ðsÞ ¼ a p e p e j xps where a p is a modal amplitude and the natural frequency x p is defined by Eq. ( 99).There is generally a phase angle as well, which is omitted because its presence does not affect the arguments that follow.Noting that element i of e p can be written as w ðpÀ1ÞðiÀ1Þ N ¼ e ju p ðiÀ1Þ , the free response of sector i can be written in real form as where C p ¼ x p =u p and the function U p ðvÞ is defined by Eq. ( 65).Equation ( 101) is a function of continuous time s and it is discretized according to the sector number i.In this way, it is endowed with the same discrete temporal and continuous spatial duality that is described in Sec.3.2.2 in the context of traveling wave engine order excitation.That is, it can be regarded as the time-response of individual (discrete) sectors, or a continuous spatial variation of displacements among the sectors that evolves with increasing time (i.e., a traveling wave).The propagating waveform is strictly a BTW in the negative i-direction (descending sector number) with wavelength 2p=u p ¼ N=ðp À 1Þ and speed C p , an illustration of which is shown in Fig. 10.However, depending on the value of p, this gives rise to SW, BTW, or FTW mode shapes, a property that follows analogously from the features described in Sec.3.2.2,where it is seen that Eq. ( 101) has the same form as Eq.(66).
For the special case of p ¼ 1 it is clear from Eq. ( 101) that each sector behaves identically with the same amplitude and the same phase because u 1 ¼ 0. An additional special case occurs when p ¼ (N þ 2)/2 if N is even.Then u ðNþ2Þ=2 ¼ p and each sector has the same amplitude but adjacent sectors oscillate with a 180-deg phase difference.In this case, the vibration modes p 2 P O;E SW correspond to SW mode shapes whose characteristics can be visualized in Figs.7(b) and 7(d) by replacing the amplitude F with a p .The remaining mode shapes correspond to repeated natural frequencies and are either BTWs or FTWs.In particular, the normal modes p 2 P O;E BTW (resp.p 2 P O;E FTW ) are backward (resp.forward) traveling waves and can be visualized in Fig. 7(a) (resp.Fig 7(c)).If mode k 2 P O;E BTW is a BTW corresponding to a natural frequency x k , then the attendant FTW mode is N þ 2 À k 2 P O;E FTW with the same natural frequency x Nþ2Àk ¼ x k .Figure 11 illustrates the normal modes of free vibration for a model with N ¼ 100 sectors.In this figure, the extent of the radial lines represents sector displacements.Those appearing outside the hub are to be interpreted as being positively displaced relative to their zero positions, and the opposite is true for lines inside the hub.Modes 1 and 51 are SWs, modes 2-50 are BTWs, and modes 52-100 are FTWs.Finally, the number of nodal diameters can be clearly identified in Fig. 11.For example, modes 4 and 98 have 3 n.d.The basic resonance structure shown in Fig. 12(a) for n 2 N essentially aliases relative to the total number of sectors, in the sense that the excited modes for n ¼ mN þ 1; …; ðm þ 1ÞN with m 2 Z þ are the same as those for n 2 N .This follows from the orthogonality condition given by Eq. ( 96) and is manifested in Eq. ( 97), which gives a relationship for the excited mode in terms of the engine order n and total number of sectors N. Because n > 0 by assumption (see Sec. 3.2) the first mode (p ¼ 1) is excited when n ¼ mN ¼ 10; 20; 30; …, the second mode (p ¼ 2) is excited when n ¼ 1 þ mN ¼ 1; 11; 21; …, and so on.Table 6 summarizes these conditions for a model with N ¼ 10 sectors and the corresponding resonance structure for n ¼ N À 1; …; 20N is shown in Fig. 13 While each engine order excites only a single mode, realistic excitation it composed of multiple harmonics (that is, orders), so that many modes can be excited.The nature of the natural frequencies and the order excitation lines leads to nontrivial resonance behavior even in the case of perfect symmetry.Of course, as noted elsewhere in this paper, imperfections that disturb the symmetry lead to even more complicated responses, in which every intersection between e.o. and natural frequency lines can lead to a resonance.These are especially important when the intersector coupling is small.x p sÞ with amplitude a p , wavelength 2p=u p 5 N=ðp 2 1Þ, and speed C p 5 x p =u p 3.4 Cyclic System With Two DOFs Per Sector.This example generalizes the one DOF per sector model of Sec.3.3 to a simple system with two DOFs per sector, which demonstrates the process of block diagonalizing the EOM when there are multiple DOFs per sector.The mathematics of the decoupling process described here applies equally as well to models with two or N DOFs per sector.Of course, the nature of the natural frequencies and mode shapes depend on the details of each sector model which, for the cyclic system considered here, is discussed in Refs.[92][93][94].There is much more to the topic of multiple DOFs per sector; the reader is referred to the works of Ottarsson [97,114] and Bladh [29,[115][116][117] for more details and more complex examples.

Equations of Motion.
The nondimensional bladed disk model shown in Fig. 14(a) consists of a rotationally periodic array of N identical, identically coupled sector models (Fig. 14(b)).The disk has radius d and rotates with a fixed speed r about an axis through C. Each blade is modeled by a simple pendulum with unit mass and length, the dynamics of which are captured by the normalized angles x i with i 2 N .The blades are attached to the rotating disk via linear torsional springs with unit stiffness, and adjacent blades are elastically coupled by linear springs with stiffness .It is assumed that the springs are unstretched when the blades are in a purely radial configuration, that is, when each x i ¼ 0. As shown in the inset of Fig. 14(b), each blade is fitted with a pendulum like, circular-path vibration absorber with radius c and mass l at an effective distance a along the blade length.The absorber dynamics are captured by the normalized pendulum angles y i , which are physically limited to jy i j 1 by stops that represent the rattling space limits imposed by the blade geometry.This feature is included for generality, but in all of what follows it is assumed that jy i j < 1, i.e., that impacts do not occur.Linear viscous damping is also included at the spring locations, but is not indicated in Fig. 14.Blade and interblade damping is captured by linear torsional and translational dampers with constants n b and n c , respectively, and the absorber damping is captured by a torsional damper with constant n a .Finally, the system is subjected to the traveling wave dynamic loading defined by Eq. ( 80), as shown in Fig. 14(b).
Sector Model.The EOM for each two-DOF sector are derived using Lagrange's method and linearized for small motions of the primary and absorber systems, that is, for small x i and y i .Then for each i 2 N , the dynamics of the ith sector are governed by [92,93] where Eq. (102a) describes the absorber dynamics and Eq.(102b) describes the blade dynamics.The indices i are taken modN such that x Nþ1 ¼ x 1 and x 0 ¼ x N , which are cyclic boundary conditions implying that the Nth blade is coupled to the first.In matrix-vector form, and for each i 2 N , Eq. ( 102) becomes where z i ¼ x i ; y i ð Þ T captures the sector dynamics, f ¼ f ; 0 ð Þ T is a sector forcing vector, and the elements of the sector mass, damping, and stiffness matrices are defined in Table 7.The matrices capture the interblade coupling and vanish if n c ¼ ¼ 0, in which case Eq. ( 103) describes the forced motion of N isolated blade/ absorber systems.
System Model.By stacking each z i into the configuration vector q ¼ z 1 ; z 2 ; …; z N ð Þ T , the governing matrix EOM for the overall 2N-DOF system takes the form where b M 2 BCBS 2;N is block diagonal with diagonal blocks M and b K 2 BCBS 2;N has generating matrices K þ 2K c ; ÀK c ; 0; …; 0; ÀK c .The matrix b C 2 BCBS 2;N is similarly defined by replacing K with C and K c with C c in b K.In terms of the circulant  (106) The 2N Â 1 system forcing vector is b f ¼ fe j/ 1 ; fe j/ 2 ; …; fe j/ N À Á T ¼ f 0 f (107) where the N Â 1 vector f 0 is defined by Eq. ( 69) and the interblade phase angle / i is given by Eq. ( 64).

Forced
Response.The forced response of the overall system defined by Eq. ( 105) can be handled directly using standard techniques [99].Its nonresonant solution in the steady-state follows in the usual way and is given by where b C is the system impedance matrix of dimension 2N Â 2N.However, Eq. ( 108) does not offer any insight into the system's modal characteristics and it requires computation of b Z À1 , which can be prohibitive for practical bladed disk models with many sectors and many DOFs per sector.We thus turn to a decoupling strategy that exploits the system symmetry and the theory developed in Sec. 2. The analysis follows similarly to that presented in Sec.3.1, except in this case the single coupled 2N-DOF system is transformed into a set of N block decoupled two-DOF systems.To this end, we introduce the change of coordinates where E is the N Â N complex Fourier matrix and e i is its ith column, is the Kronecker product, I is the 2 Â 2 identity matrix (the dimension of I corresponds to the number of DOFs per sector), and u ¼ u 1 ; u 2 ; …; u N ð Þ T is a vector of modal, or cyclic coordinates.Each u p is 2 Â 1 and describes the sector dynamics in modal space.Substituting Eq. ( 109) into Eq.( 105 where ðe H p IÞ b f is the pth 2 Â 1 block of ðE H IÞ b f.Equation ( 110) is analogous to the N M-DOF systems given by Eq. ( 61) for the general formulation in Sec.3.1, but in this case M ¼ 2 and engine order excitation is assumed from the onset.Figure 15 illustrates the transformation of the single 2N-DOF system given by Eq. (105) to a system of N block decoupled two-DOF forced oscillators defined by Eq. ( 110).
The 2 Â 2 mass, damping, and stiffness matrices associated with the pth mode follow from Theorem 10 and are given by where u p is defined by Eq. ( 16), the elements of M, C, and K are defined in Table 7 and the coupling matrices C c and K c are defined by Eq. ( 104).In light of Eqs.where 0 ¼ (0,0) T and the scalar product e H p f 0 vanishes except for p ¼ n þ 1.Because only mode p ¼ n þ 1 is excited, u nþ1 ðsÞ is the only nonzero modal response in the steady-state.
Assuming harmonic motion, and in light of Eq. ( 112), the pth steady-state modal response follows easily from Eq. ( 110) and is given by where w nðiÀ1Þ ¼ e j/ i is employed.From Eq. ( 115) it is clear that each blade/absorber combination behaves identically except for a constant phase shift from one sector to another, which is captured by the interblade phase angle / i .This approach offers a significant computational advantage over the direct solution to the full 2N-DOF system, as given by Eq. ( 108).

Conclusions
The goal of this paper is to provide the mathematical tools for handling circulant matrices as they apply to the free and forced vibration analysis of structures with cyclic symmetry.As demonstrated by past work in this area and the review provided here, the theory of circulants provides a useful description of the fundamental structure of the mode shapes and spectrum of systems with cyclic symmetry, including those of large scale.The theory also provides a convenient means for computing the vibration response of these systems, even when the idealized symmetry is broken by mistuning or by nonlinear effects.As with any mathematical tool, the overhead in learning it must provide appropriate benefit, whether in terms of fundamental understanding, insight, or ease of computation.We trust that the results presented here offer such benefits to readers interested in vibration analysis of cyclic systems.
It must be noted that no physical system has perfect symmetry, as assumed herein.This assumption must be examined in light of the system under consideration.One key to the suitability of a cyclically symmetric model is the intersector coupling.If the coupling is strong, so that the pairs of modes have well separated frequencies, then small imperfections will not alter the picture substantially, and one can consider each mode pair as robust against coupling to other modes.However, if the coupling is weak, so that the system frequencies are clustered near those of the isolated sector model, the possibility of localization is significantly increased.This topic has been investigated quite thoroughly, primarily in the context of the vibration of bladed disk assemblies with small blade mistuning; see, for example, [72][73][74][75][76][77][78][79][80][81][82][83][84][85][86][87][88].Also, nonlinear effects can couple linear modes under certain resonance conditions, even at small amplitudes [118].In cyclic systems this can occur for the pairs of modes with equal frequencies [119], and this possibility expands to groups of modes for the case of weak coupling, leading to extremely complicated behavior [82,84,88].In such cases, the tools from group theory can be applied to categorize the possible modes and forced response in terms of their symmetries [65].This topic, while interesting, is outside the scope of the present paper.

Fig. 1 (
Fig.1(a) Finite element model of a bladed disk assembly[1] and (b) general cyclic system with N identical sectors and nearest-neighbor coupling

5 . 2 . 5 ¼
with ones along the superdiagonal and in the (N, 1) position, and zeros otherwise.᭝ The cyclic forward shift matrix plays a key role in the representation and diagonalization of circulant matrices, which are discussed in Secs.2.4 and 2.Example Let a ¼ (a, b, c) be a three-vector.Then the operation ar 3 ¼ a b c ½ ðc; a; bÞ

DEFINITION 10 (
Generating Elements).Let the N Â N circulant matrix C be given by Definition 9. Then the elements of the N-vector fc 1 ; c 2 ; …; c N g are said to be the generating elements of C.

Example 3 .
The circulant array formed by the generating elements a, b, c, d can be written as circða; b; c; dÞ ¼ circulant matrix of type 4.

c 21 c
) which are necessarily repeated.Only (N þ 2)/2 generating elements are distinct if N is even and (N þ 1)/2 are distinct if N is odd.The set of all N Â N symmetric circulants is denoted by SC N .A matrix contained in SC N is said to be a symmetric circulant of type N. Example 4. The 5 Â 5 matrix circða; b; c; c; bÞ ¼ and circulant.Because N ¼ 5 is odd, it has (N þ 1)/2 ¼ 3 distinct elements.The matrix defined in Example 3 is not a symmetric circulant because its generating elements are distinct.Next, we give a necessary and sufficient condition for a square matrix to be circulant.THEOREM 4. Let r N be the cyclic forward shift matrix.Then a N Â N matrix C is circulant if and only if Cr N ¼ r N C. ٗ Proof.Let C be an N Â N matrix with arbitrary elements c ik for i; k ¼ 1; 2; …; N. Then Cr N ¼ c 1N c 11 c 12 Á Á Á c 1ðNÀ1Þ c 2N c 21 c 22 Á Á Á c 2ðNÀ1Þ 22 c 23 Á Á Á c 2N c 31 c 32 c 33 Á Á Á c 3N which is a N Â N circulant matrix with generating elements c 11 ; c 12 ; …; c 1N .Any matrix that commutes with the cyclic forward shift matrix is, therefore, a circulant.Theorem 4 also says that circulant matrices are invariant under similarity transformations involving the cyclic forward shift matrix.That is, C is similar to itself for a similarity transformation using r N .Example 5. Consider the 3 Â 3 matrix

DEFINITION 12 (
Generating Matrices).Let the NM Â NM block circulant C be given by Definition 11.Then the elements of the ordered set fC 1 ; C 2 ; Á Á Á ; C N g are said to be the generating matrices of C.
w ðkÞ N ¼ e j 2pk N for integers k ¼ 0; 1; 2; …; N À 1. ᭝ DEFINITION 15 (Primitive Nth Root of Unity).The primitive Nth root of unity is denoted by w N ¼ e j 2p N which corresponds to k ¼ 1 in Definition 14. ᭝ COROLLARY 6.The integer powers w k N of the primitive Nth root of unity are equivalent to the distinct Nth roots of unity w ðkÞ N for k ¼ 1; 2; …; N. ٗ Proof.Consider w k N ¼ ðe j 2p N Þ k ¼ e j 2pk N ¼ w ðkÞ N , which follows from Definition 14.Thus, the powers 1; w 1 N ; w 2 N ; …; w NÀ1 N are equivalent to the distinct Nth roots of unity 1; w

10 .
Let e i denote the ith column of the Fourier matrix E N and d ik be the Kronecker delta.Then because e 2 and e 3 are mutually orthogonal.COROLLARY 11.If E N is the N Â N Fourier matrix and I M is the identity matrix of dimension M, then the NM Â NM matrix E N I M is unitary.ٗ Proof.Consider the matrix product k ÞÞ follows from Corollary 16 and Definition 13.The desired result is obtained by multiplying from the left by E H N , multiplying from the right by E N , and invoking Theorem 6.

¼
diagð2; 4; 6; 4Þ matrix.The diagonal elements are the eigenvalues of C, which is stated explicitly for general C in the following corollary.They may be computed directly using Eq.(21), as it is done in Example 15, but Corollary 17 cannot be used because C is not symmetric.The eigenvalue magnitudes of the nonsymmetric matrix C in Example 16 are observed to exhibit the same multiplicity (and symmetry) as the eigenvalues in Example 15 for a symmetric circulant.It is shown in Sec.2.8.3 that the eigenvalues of any circulant matrix C 2 C N with real-valued generating elements exhibit

2. 5 . 5
Block Diagonalization of a Block Circulant.Theorem 9 is generalized to handle block circulants using the Fourier and identity matrices together with the Kronecker product.The choice of diagonalizing matrix E H N I M is discussed in Sec.2.6, where generalizations of Theorem 10 are considered.THEOREM 10 (Block Diagonalization of a Block Circulant).Let C 2 BC M;N and denote its M Â M generating matrices by C 1 ; C 2 ; …; C N .Then if E N is the N Â N Fourier matrix and I M is the identity matrix of dimension M, For example, if each C k is a circulant of type M, then the additional choice of B ¼ E M and ðÁÞ # ¼ ðÁÞ H fully diagonalizes a block circulant matrix C 2 BC N;M with circulant blocks.COROLLARY 22.Let C 2 BC M;N have generating matrices C 1 ; C 2 ; …; C N 2 C M and denote the generating elements of each C i by c ð1Þ i ; c ð2Þ i ; …; c ðMÞ i .Then r ¼ 0; 1; …; N À 1. ᭝ COROLLARY 23.The DFT sinusoids are orthogonal.ٗ Proof.Consider the DFT sinusoids S i ðrÞ ¼ e j 2p N ir S k ðrÞ ¼ e j 2p N kr ) ; i; k; r ¼ 0; 1; …; N À 1

2. 8 . 1
Standard Circulant Eigenvalue Problem.Consider the block circulant matrix C 2 BC M;N with arbitrary generating matrices C 1 ; C 2 ; …; C N .The standard cEVP involves determination of the scalar eigenvalues k and NM Â 1 eigenvectors q of C such that

) ( 50 )
which facilitates the results for real generating elements and clarifies their interpretation.If the eigenvalues are dissected according to YðrÞ ¼ Y R ðrÞ þ jY I ðrÞ, then the real and imaginary components (21).Then YðN À rÞ ¼ YðrÞ ¼ YðÀrÞ for r ¼ 0; 1; 2…; N À 1. ٗ Proof.The eigenvalues of y(p) are given by YðrÞ ¼ Y R ðrÞ þ jY I ðrÞ, where Y R ðrÞ ¼ ) also holds if N is odd, but the Nyquist component is repeated with multiplicity of two.
(b), this amounts to the same standing wave excitation as the n ¼ N case, except for a phase reversal in the excitation among adjacent blades.The engine orders corresponding to SW excitations for odd and even N are denoted by the sets N O;E SW & N , and all other values of n 2 N correspond to traveling waves.Engine orders n 2 N O;E BTW (resp.n 2 N O;E FTW ) correspond to BTW (resp.FTW) excitation, an example of which is shown in Fig. 7(a) (resp.Fig. 7(c)), where N O;E BTW and N O;E FTW are also defined in Table

Fig. 6
Fig. 6 Example illustration of the discrete temporal and continuous spatial variations of the traveling wave excitation defined by Eq. (63) in real form: (a) the discrete dynamic loads with amplitude F and period T 5 2p=nX applied to each sector; and (b) the continuous BTW excitation with wavelength N/n and speed C 5 NX=2p relative to the rotating hub

Fig. 7
Fig. 7 Engine orders n mod N corresponding to BTW, FTW and SW applied dynamic loading for (i) odd N and (ii) even N (see also Table 4); example plots of applied dynamic loading (represented by the dots) for a model with N 5 10 sectors and with (a) n 5 1 (BTW), (b) n 5 5 (SW), (c) n 5 9 (FTW), and (d) n 5 10 (SW).The BTW engine order excitation is represented by the solid lines.

Fig. 8
Fig. 8 Linear cyclic vibratory system with N sectors and one DOF per sector

Fig. 9
Fig. 9 Dimensionless natural frequencies x p in terms of the number of n.d.versus mode number p for WC and SC: (a) N 5 11 (odd) and (b) N 5 10 (even).Also indicated below each figure is, for general N, the number of n.d. at each value of p and also the mode numbers corresponding to SW, BTW, and FTW.

3. 3 . 5
Resonance Structure.In general, there may be a system resonance if the excitation frequency matches a natural frequency, that is, if nr ¼x p .These possible resonances are conveniently identified in a Campbell diagram, an example of which is shown in Fig.12(a) for engine orders n 2 N (the general case of n 2 Z þ is considered below), N ¼ 10, and ¼ 0:5.The natural frequencies are plotted in terms of the dimensionless rotor speed and several engine order lines nr are superimposed.Possible resonances correspond to intersections of the order lines and eigenfrequency loci.There are ðN þ 2Þ=2 such possibilities for each engine order if N is odd and ðN þ 1Þ=2 possible resonances if N is even.In light of Eq. (96), however, there is only a single resonance associated with each n under the traveling wave dynamic loading of Sec.3.2, which corresponds to mode p ¼ n mod N þ 1.The set of N resonances for a system excited by N engine orders ðn ¼ 1; 2; …; NÞ are indicated by the black dots in Fig. 12(a) and the corresponding frequency response curves jq ss i ðsÞj (for each n) are shown in Fig. 12(b) for a model with f ¼ 0.01.For example, a 3 e.o.excitation resonates mode 4 (p ¼ 4), which is a BTW with 3 n.d.Mode 8 (p ¼ 4) also has 3 n.d. and is excited by a 7 e.o.excitation.The TW and n.d.designations can be verified in Fig. 9.
(a).Each collection of resonance points n ¼ mN þ 1; …; ðm þ 1ÞN is qualitatively the same in structure.However, for m > 1 the resonances become increasingly clustered, which is shown in Fig.13(b) for n ¼ N; …; 2N.In terms of the sets defined in Tables4 and 5, an engine order n mod N 2 N O;E SW excites a SW mode p 2 P O;E SW .Similarly, an engine order n mod N 2 N O;E FTW (resp.n mod N 2 N O;E BTW ) excites a FTW (resp.BTW) mode p 2 P O;E FTW (resp.p 2 P O;E BTW ).

Fig. 10 A
Fig. 10 A backward traveling wave a p U p ði À 1 1 C p sÞ 5 a p cosðu p ði 2 1Þ 1x p sÞ with amplitude a p , wavelength 2p=u p 5 N=ðp 2 1Þ, and speed C p 5x p =u p

Fig. 11
Fig. 11 Normal modes of free vibration for a model with N 5 100 sectors.Mode 1 consists of a SW, in which each sector oscillates with the same amplitude and phase.Mode 51 also corresponds to a SW, but neighboring oscillators oscillate exactly 180 deg out of phase.Modes 2-50 (resp.52-100) consist of BTWs (resp.FTWs).

Fig. 15
Fig. 15 The topology of a bladed disk assembly fitted with absorbers in (a) physical space and (b) modal space.The modal transformation qðsÞ 5 ðE IÞuðsÞ reduces the cyclic array of N, two-DOF sector models ðB; AÞ, which together form a 2N-DOF coupled system, to a set of N, two-DOF block decoupled models ðB p ; A p Þ.

Table 1
Selected special matrices

Table 2
Selected types of linear transformations Equivalence P and Q are nonsingular for the case of N ¼ 4. DEFINITION 16.Let w N be the primitive Nth root of unity.Then

Downloaded From: http://appliedmechanicsreviews.asmedigitalcollection.asme.org/ on 02/19/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use Example
9. For the special case of N ¼ 4, the Fourier matrix is given by NÂNwhere w N is the primitive Nth root of unity and N 2 Z þ .᭝ Fig. 2 Example plots of the distinct Nth roots of unity 040803-8 / Vol.66, JULY 2014 Transactions of the ASME the relationship between the DFT and IDFT in Sec.2.7.LEMMA 2. Let w N be the primitive Nth root of unity with N 2 Z þ .Then , which is shown in Sec.2.8.In light of Corollary 2, we have the following results.COROLLARY 14.Let E N and r N be the N Â N Fourier and cyclic forward shift matrices.Then for any n which are given by Definition 18.In fact, all circulant matrices contained in C N share the same eigenvectors 040803-12 / Vol.66, JULY 2014 Transactions of the ASME Downloaded From: http://appliedmechanicsreviews.asmedigitalcollection.asme.org/ on 02/19/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-usee i ik where d ik is the Kronecker delta and k i is defined by Eq. (21) forC 2 C N or Corollary 17 if C 2 SC N .ٗThus, the columns of the Fourier matrix are mutually orthogonal (Corollary 10) and orthogonal with respect to any circulant matrix (Corollary 18), not just the cyclic forward shift matrix and its integer powers (Corollary 15).
which is recognized to be the third diagonal element k 3 of the matrix E H 4 CE 4 in Example 15.However, the scalar vanishes, as expected, because i 6 ¼ k in Corollary 18 such that d ik ¼ 0.
S 0 ð1Þ S 0 ð2Þ S 0 ð3Þ S 1 ð0Þ S 1 ð1Þ S 1 ð2Þ S 1 ð3Þ S 2 ð0Þ S 2 ð1Þ S 2 ð2Þ S 2 ð3Þ S 3 ð0Þ S 3 ð1Þ S 3 ð2Þ S 3 ð3Þ E 4 follows from Example 9 or by elementwise direct computation according to Definition 19.Alternatively, each element X(k) of X 4 can be computed using the summation given in Definition 20.If k ¼ 3, for example, S 2 ð1Þ S 2 ð2Þ S 2 ð3Þ S 3 ð0Þ S 3 ð1Þ S 3 ð2Þ S 3 ð3Þ all circulants contained in C N share the same N Â 1 eigenvectors e 1 ; e 2 ; …; e N , and each e i is orthonormal with respect to C. The eigenvalues of a matrix C 2 C N with generating elements c 1 ; c 2 ; …; c N are given by Eq. (21), or equivalently by Eq. (22).If C 2 SC N , the eigenvalues follow from Corollary 17.
The same basic results also hold for the nonsymmetric circulant C ¼ circð4; À1; 0; 1Þ 2 C 4 in Example 16, where it is shown that E H 4 CE 4 ¼ diagð4; 4 À 2j; 4; 4 þ 2jÞ.In this case, the eigenvalues are given by aðCÞ ¼ f4; 4 À 2j; 4; 4 þ 2jg and the corresponding eigenvectors are the same as those in Example 24 because all circulants contained in C N share the same linearly independent eigenvectors.
2.8.1 and 2.8.2.It is first shown that each ũðpÞ i is

Table 3
Magnitudes and arguments of Y(r) for even N 5 8 2.8.1 and 2.8.2 are based, are valid for arbitrary M Â M generating matrices.COROLLARY 27.Suppose each K i defined by Eq. (23) is symmetric.
T ; k Þ be the orthonormal modal matrix associated with K i and C 1 ; C 2 ; …; C N be the symmetric generating matrices of C. Corollary 19 guarantees that K i is symmetric because the generating matrices are symmetric.By setting B ¼ e U i and ðÁÞ # ¼ ðÁÞ T in Corollary 21, it follows that However, in light of the diagonal structure of each W i , only the diagonal elements survive in this expansion.That is, i ; …; e k ũðqÞ i ; …; e k ũðMÞ i Þ which produces an M Â M array with scalar elements Applied Mechanics Reviews JULY 2014, Vol.66 / 040803-25 Downloaded From: http://appliedmechanicsreviews.asmedigitalcollection.asme.org/ on 02/19/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-useiÞin the (p, q) position for p; q ¼ 1; 2; …; M.

Table 5
Sets of mode numbers p ‰ N corresponding to BTW, FTW, and SW normal modes of free vibration for odd and even N of Sec.3.3.4,mode 1 has 0 n.d., modes 2 and 100 have 1 n.d., modes 3 and 99 have 2 n.d., and so on.are described next, where it is shown that each can be categorized as a SW, BTW, or FTW.3.3.4Normal Modes of Vibration.It was shown that Eq. ( 040803-34 / Vol.66, JULY 2014 Transactions of the ASME Downloaded From: http://appliedmechanicsreviews.asmedigitalcollection.asme.org/ on 02/19/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use of vibration