Cyclic-symmetry reduction ========================= A rotor built from :math:`N` geometrically-identical sectors spanning :math:`360°` admits a *cyclic-symmetry* modal reduction: the full-rotor eigenproblem decouples into :math:`\lfloor N/2 \rfloor + 1` smaller eigenproblems, one per **harmonic index**. This page derives the reduction from the rotational-symmetry property of :math:`K` and :math:`M`, shows how femorph-solver maps it onto the real, base-sector linear-algebra primitives the :mod:`femorph_solver.solvers.cyclic` backend works with, and lists the references the kernel cites. The driving API (sector mesh → cyclic faces → per-harmonic modal results) lives in :class:`~femorph_solver.CyclicModel`; see :doc:`/user-guide/solving/cyclic` for the user-facing walk-through. Setting ------- A rotor with :math:`N` identical sectors has assembled matrices :math:`K, M \in \mathbb{R}^{Nn \times Nn}`, where :math:`n` is the per-sector free-DOF count. Index the sectors :math:`j = 0, 1, \dots, N-1` and partition each matrix accordingly: .. math:: K = \begin{bmatrix} K_{00} & K_{01} & 0 & \cdots & K_{0,N-1} \\ K_{10} & K_{11} & K_{12} & \cdots & 0 \\ 0 & K_{21} & K_{22} & \cdots & 0 \\ \vdots & & & \ddots & \vdots \\ K_{N-1,0} & 0 & \cdots & K_{N-1,N-2} & K_{N-1,N-1} \end{bmatrix}, with :math:`K_{ij} = 0` whenever :math:`|i - j| > 1 \pmod{N}` — i.e. each sector only couples to its two cyclic neighbours through the shared cyclic face. Cyclic symmetry of the geometry implies .. math:: K_{j,j} = K_{0,0}, \qquad K_{j,j+1} = K_{0,1}, \qquad K_{j,j-1} = K_{0,1}^\top \quad \forall j. — every diagonal block is the same self-stiffness :math:`A \equiv K_{00}`, and every off-diagonal block is the same inter-sector coupling :math:`B \equiv K_{01}`. The same structure holds for :math:`M`; call its blocks :math:`A^M, B^M`. Block-circulant decomposition ----------------------------- Matrices of the form .. math:: K = \begin{bmatrix} A & B & 0 & \cdots & 0 & B^\top \\ B^\top & A & B & \cdots & 0 & 0 \\ \vdots & & & \ddots & & \vdots \\ B & 0 & \cdots & 0 & B^\top & A \end{bmatrix} are **block circulant**. The discrete Fourier transform on the sector index block-diagonalises them — *exactly*, not approximately. Define the per-sector phase shift :math:`\theta_k = 2\pi k / N` for :math:`k \in \{0, 1, \dots, N-1\}` and the harmonic basis vectors .. math:: \mathbf{e}_k = \frac{1}{\sqrt{N}} \bigl[ 1,\; e^{i \theta_k},\; e^{i 2 \theta_k},\; \dots,\; e^{i (N-1) \theta_k} \bigr]^\top. The block-Fourier projection .. math:: \tilde{K}_k \;=\; A + e^{i \theta_k} B + e^{-i \theta_k} B^\top gives :math:`N` reduced :math:`n \times n` complex-Hermitian problems, each independent of the others: .. math:: \tilde{K}_k\, \tilde{\boldsymbol{\varphi}}_k \;=\; \omega^2_k\, \tilde{M}_k\, \tilde{\boldsymbol{\varphi}}_k, \qquad k = 0, 1, \dots, N-1, with :math:`\tilde{M}_k = A^M + e^{i \theta_k} B^M + e^{-i \theta_k} {B^M}^\top` defined identically. The full-rotor mode is recovered by re-expansion: .. math:: \boldsymbol{\varphi}_k^{(j)} \;=\; e^{i j \theta_k}\, \tilde{\boldsymbol{\varphi}}_k, \qquad j = 0, 1, \dots, N-1. — the same base-sector mode shape repeated around the rotor with a per-sector phase advance :math:`\theta_k`. Physical harmonic indices ------------------------- The :math:`N` Fourier indices :math:`\{0, 1, \dots, N-1\}` contain redundancy: the spectrum at :math:`k` is the complex conjugate of the spectrum at :math:`N - k`. For a real rotor, this collapses to the **physical harmonic indices** .. math:: k \in \{0, 1, \dots, \lfloor N/2 \rfloor\}. * :math:`k = 0` — every sector moves in phase; axisymmetric (umbrella) modes; mode shape is real. * :math:`k = N/2` (only when :math:`N` is even) — adjacent sectors move :math:`180°` out of phase; mode shape is real. * :math:`0 < k < N/2` — travelling-wave modes that come in forward / backward pairs of identical frequency, often called the **standing-wave doublets**. The two members of a doublet are the real and imaginary parts of a single complex eigenvector. The "nodal diameter" terminology used in turbomachinery (forward / backward :math:`k`-nodal-diameter modes) is exactly this harmonic index. The number of distinct frequencies a cyclic modal solve returns is therefore not :math:`N \cdot n_\text{modes}` but :math:`(\lfloor N/2 \rfloor + 1) \cdot n_\text{modes}`, counting each :math:`0 < k < N/2` mode once with multiplicity two. Real :math:`2n` augmentation ---------------------------- The reduced problem :math:`\tilde{K}_k\, \tilde{\boldsymbol{\varphi}}_k = \omega_k^2\, \tilde{M}_k\, \tilde{\boldsymbol{\varphi}}_k` is complex-Hermitian. Most mature sparse eigensolvers (ARPACK, LOBPCG, PRIMME) ship a real-symmetric path optimised for the SPD case; the complex- Hermitian path is less mature and slower per iteration. Grimes, Lewis & Simon (1994) pointed out that any complex- Hermitian generalised eigenproblem can be lifted to a real-symmetric :math:`2n \times 2n` problem with the same spectrum. Write :math:`\tilde{K}_k = K_R + i K_I`, :math:`\tilde{M}_k = M_R + i M_I` (with :math:`K_R, M_R` symmetric and :math:`K_I, M_I` antisymmetric since the original :math:`A, B` are real) and form .. math:: K_k^{\,2n} = \begin{bmatrix} K_R & -K_I \\ K_I & K_R \end{bmatrix}, \qquad M_k^{\,2n} = \begin{bmatrix} M_R & -M_I \\ M_I & M_R \end{bmatrix}. Both lifted matrices are real symmetric (and SPD, when the original sector matrices are SPD on the constraint quotient). Their generalised eigenvalues are exactly the eigenvalues of the complex-Hermitian problem, each with double multiplicity that mirrors the real / imaginary part decomposition. femorph-solver uses this real-:math:`2n` lift on every reduced harmonic-index problem so the same shift-invert ARPACK / SPD- direct path the standard modal solve uses applies unchanged. Practical mapping in femorph-solver ----------------------------------- The user-facing path is: 1. Author the **base-sector** mesh as a regular :class:`~femorph_solver.Model` carrying one wedge of the rotor. No need to assemble the full rotor. 2. Identify the **cyclic faces** — the two θ-bounding surfaces of the wedge. Either supply them explicitly via the ``low_face`` / ``high_face`` keyword arguments, or let :class:`~femorph_solver.CyclicModel` auto-detect them by matching every node on the low face to its rotated image on the high face within a :math:`pair\_tol` tolerance. 3. Call :meth:`CyclicModel.solve_modal ` with the harmonic indices of interest (default: every physical index :math:`0 \le k \le N/2`). One per-:math:`k` reduced eigenproblem runs internally, on the :math:`2n`-lifted real-symmetric system. 4. Each :class:`~femorph_solver.solvers.cyclic.CyclicModalResult` carries the harmonic index :math:`k`, the per-mode frequencies, and the base-sector mode shapes. Full-rotor mode-shape recovery is handled by the result's expansion helpers. The compute saving over assembling the full rotor scales as :math:`O((\lfloor N/2 \rfloor + 1) \cdot n^p / N \cdot n^p) = O((\lfloor N/2 \rfloor + 1) / N^p)` for an :math:`O(n^p)` solver — close to the obvious "factor of :math:`N`" speedup, plus an extra factor :math:`p-1` from the better factorisation cost on the smaller system. For a typical :math:`N = 36`-sector turbine rotor this is ~50× faster than the full-rotor solve. Validity assumptions -------------------- Cyclic-symmetry reduction is exact when: * The sector geometry repeats *identically* under the per-step rotation — no mistuning, no per-sector mass / stiffness variation. Real turbomachinery is mistuned at the 0.1-1 % level, which the cyclic reduction *cannot* capture; mistuned analysis requires a separate technique (component mode synthesis or full-rotor solve). * The cyclic faces of adjacent sectors mate one-to-one. Mesh generators sometimes produce slightly-misaligned faces that the femorph-solver auto-pairing can match within tolerance, but mismatch beyond that fails the pairing check explicitly. * The boundary conditions on the non-cyclic faces of every sector are identical. A bracket bolted to the rotor at one sector position breaks cyclic symmetry. Violations of these assumptions don't produce a wrong answer silently — the pairing check catches the geometry mismatch, the modal-solve sanity check catches BC asymmetry, and the mistuning case is intercepted by the user reading these constraints up front. References ---------- * **D. L. Thomas**, "Dynamics of rotationally periodic structures," *International Journal for Numerical Methods in Engineering* 14 (1979), 81–102 — the foundational paper for cyclic-symmetry FE reduction in modal analysis. * **A. Wildheim**, "Excitation of rotationally periodic structures," *Journal of Applied Mechanics* 46 (1979), 878–882 — companion paper covering the harmonic-excitation case. * **R. Grimes, J. Lewis, H. Simon**, "A shifted block Lanczos algorithm for solving sparse symmetric generalized eigenproblems," *SIAM J. Matrix Anal. Appl.* 15 (1994), 228–272 — origin of the real-:math:`2n` augmentation that lets a real-symmetric Lanczos drive the complex-Hermitian per-harmonic problem. * **R. R. Craig & A. J. Kurdila**, *Fundamentals of Structural Dynamics*, 2nd ed., 2006, §16.3 — accessible textbook treatment of cyclic / dihedral group representations applied to FE eigenproblems. * **A. V. Srinivasan**, "Flutter and resonant vibration characteristics of engine blades," *Journal of Engineering for Gas Turbines and Power* 119 (1997), 742–775 — application to bladed-disc / rotor problems with mistuning. See also -------- * :doc:`eigenproblem` — the underlying generalised eigenproblem theory. * :doc:`/user-guide/solving/cyclic` — driving API. * :doc:`/getting-started/known-limitations` — mistuning and per-sector variation are explicitly out of scope. * :class:`femorph_solver.CyclicModel` — the user-facing cyclic wrapper.