SDT-rotor         Contents     Functions              PDF Index

2.3  Fourier analysis of structures

• cyclic structure basics
• Fourier transform for shaft computations
• Solutions in periodic media

For more details, you can refer to[5] that is available on-line.

2.3.1  cyclic structure basics

For cyclic system with N sectors (angle α=2π/N), a point load placed at angle nα is associated with the harmonic load

(2.9)

this can be used to simply compute the k diameter response. Using the symmetry of the spectrum, one can recover the full spatial response by inverse fourier transform. For N odd

(2.10)

for N even

(2.11)

The displacements y(n) of each sector are however expressed in local coordinates, when applying conditions a transformation θ^{α n} to global coordinates is thus needed

(2.12)

For a point on the axis, the in plane response is given by two components u(0),v(0), the response at other points of the series is given by (2.10) but the point being coincident on has {u(0),v(0)}=(θ^α)ⁿ{u(n),v(n)}. It follows that Y(k)=0 for all k≠ 1.

The intersector continuity condition says that the nodes that are common to two sectors have the same motion. By convention the left nodes in SDT rotor are those with the lowest clockwise polar angle, the left nodes of sector 0 have thus equal motion than the right nodes or sector N−1. Thus c_L−c_R(Y(k)e^{−jα k})=0. Which leads to the constraint

(2.13)

The Fourier transform being a linear relation, one can actually rewrite the relation as a {y}=[E]{Y}. For N even

(2.14)

Note that E^TE is a diagonal matrix with [N 2N ⋯ 2N N] on the diagonal.

Posing I_NQ is the identity matrix whose size is the number of physical DOF of the sector and ⊗ is the Kronecker product, the fourier DOFs Q_k (sorted as in (2.14)) and physical DOFs on the whole disk are related by

(2.15)

When performing harmonic computations, one typically uses the fact that the model matrices are identical for each sector. As a result, the transformed matrices are block-diagonal, which is the basis for the ability to compute full responses based on independent computation of each Fourier solution Q_k.

When adding a matrix to sector n, the relation to to the Fourier DOFs is given by a row E_n, the product E_n^TE_n is not block diagonal as a result one has coupling between the Fourier harmonics.

2.3.2  Fourier transform for shaft computations

From mono-harmonic modes stored in SDT, you can obtain the physical response using for complex value shape pairs

(2.16)

and for real valued vectors associated with 0 and N/2 diameters.

(2.17)

2.3.3  Solutions in periodic media

This really does not apply to rotors but is implemented in fe_cyclic, and thus documented here. One considers a model whose properties are spatially periodic. For a physical response, known at regularly spaced positions nΔ x, one can compute the its fourier transform

(2.18)

U(κ) is a complex shape defined on the mesh of the repeated cell. One actually uses two cells to represent the real and imaginary parts of U. The wave number κ varies in the [0,2π/Δ x] interval (or any interval of the same length, since U(κ) is periodic in the wavelength domain).

Given the fourier transform U(κ), on can recover the physical motion by computing the inverse Fourier transform

(2.19)

For a mono-harmonic response (fixed wave number κ), the spatial transform is given by u(nΔ x) = (U(κ)e^{jn(κ Δ x)}), using the continuity condition linked to the fact that u_left(nΔ x)=u_right((n−1)Δ x), one thus has [c_L]{U(κ)}=[c_R]{U(κ)}e^{−j(κ Δ x)}. Hence in the real/imaginary format, one has the constraint equation

(2.20)

and one solves for U(κ) knowning the transform of the applied loads. For a discrete load on the first slice n=0, the fourier coefficients of the load is F(κ)=f(0).

One can usefully note that the wave length L=π/κ covers the full interval of positive lengths, that U(2π/δ x−κ)=Ū(κ), and that the half spectrum corresponds to a wavelength of Δ x. Using the property of symmetry allows the use of computations for wavelengths larger than Δ x only. In the fe_cyclic eig calls, one specifies the period as a number of steps pΔ x. The phase angle is thus κΔ x=2π/p and the spectrum is symmetric for p′=2/p.

mo1=femesh('testhexa8b divide 2 1 1');
mo1=fe_cyclic('build -1 1 0 0',mo1);

% symmetry along right edge
r2=fe_cyclic('build -1 0 1 0',mo1);
r2=fe_case(r2,'getdata','Symmetry');
mo1=fe_case(mo1,'FixDof','y_symmetry',r2.IntNodes(:,1)+.02);

mo1=stack_set(mo1,'info','EigOpt',[2 2 1e3]);
range=1./linspace(.01,.99,21);
def=fe_cyclic(['teig -all' sprintf(' %g',range)],mo1);
figure(1);plot(1./def.data(:,2),def.data(:,1),'x')
feplot(mo1,def)

[r1,i1]=unique(def.data(:,2),'first');
Ek=fe_caseg('enerm -bygroup',mo1,fe_def('subdef',def,i1));
if norm(Ek.Y*4-1)>sqrt(eps); error('Inconsistent energy');end
if 1==2 % manual verification
 [m,k,mdof]=fe_mknl(mo1,'assemble NoT');
 def=feutilb('placeindof',mdof,def);
 feutilb('dtkt',real(def.def),m)+feutilb('dtkt',imag(def.def),m)
end