.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/analyses/cyclic/example_cyclic_modes_disk.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_gallery_analyses_cyclic_example_cyclic_modes_disk.py: Cyclic-symmetry modal sweep — annular disk ========================================== A 12-sector annular disk has its global eigenspectrum factored exactly by **harmonic index** :math:`k \in \{0, 1, \ldots, N/2\}` (Thomas 1979, Wildheim 1979). Each :math:`k` produces a reduced complex-Hermitian eigenproblem on a single base sector; femorph-solver solves the real-symmetric 2n-augmentation form (Grimes, Lewis & Simon 1994) so :func:`scipy.sparse.linalg.eigsh` can return mode pairs at every harmonic. This example builds a 12-sector disk, runs :meth:`Model.cyclic_modal_solve` over every harmonic, and prints the lowest non-rigid frequency for each — illustrating how the harmonic sweep recovers the full free-vibration spectrum at :math:`1 / N` the DOF count of an explicit full-rotor model. References ---------- * Thomas, D. L. (1979) "Dynamics of rotationally periodic structures," *Journal of Sound and Vibration* 66 (4), 585–597. * Wildheim, S. J. (1979) "Excitation of rotationally periodic structures," *Journal of Applied Mechanics* 46, 891–893. * Grimes, R. G., Lewis, J. G. and Simon, H. D. (1994) "A shifted block Lanczos algorithm for solving sparse symmetric generalized eigenproblems," *SIAM J. Matrix Analysis and Applications* 15, 228–272 (real 2n augmentation for the complex-Hermitian sub-problem). * Bathe, K.-J. (2014) *Finite Element Procedures*, 2nd ed., §10.3.4 (cyclic-symmetry modal analysis). .. GENERATED FROM PYTHON SOURCE LINES 35-44 .. code-block:: Python from __future__ import annotations import numpy as np import pyvista as pv import femorph_solver from femorph_solver import ELEMENTS .. GENERATED FROM PYTHON SOURCE LINES 45-47 Geometry — one base sector of a 12-sector annular disk ------------------------------------------------------ .. GENERATED FROM PYTHON SOURCE LINES 47-60 .. code-block:: Python N_SECTORS = 12 ALPHA = 2.0 * np.pi / N_SECTORS R_INNER = 0.40 R_OUTER = 0.55 T_DISK = 0.02 N_T, N_R, N_Z = 6, 4, 1 E = 2.0e11 NU = 0.3 RHO = 7850.0 .. GENERATED FROM PYTHON SOURCE LINES 61-63 Build the base sector mesh -------------------------- .. GENERATED FROM PYTHON SOURCE LINES 63-75 .. code-block:: Python theta = np.linspace(0.0, ALPHA, N_T + 1) r = np.linspace(R_INNER, R_OUTER, N_R + 1) z = np.linspace(0.0, T_DISK, N_Z + 1) R, T, Z = np.meshgrid(r, theta, z, indexing="ij") X = R * np.cos(T) Y = R * np.sin(T) grid = pv.StructuredGrid(X, Y, Z).cast_to_unstructured_grid() print(f"sector mesh: {grid.n_points} nodes, {grid.n_cells} HEX8 cells") print(f"full-rotor equivalent: {N_SECTORS * grid.n_points} nodes") .. rst-class:: sphx-glr-script-out .. code-block:: none sector mesh: 70 nodes, 24 HEX8 cells full-rotor equivalent: 840 nodes .. GENERATED FROM PYTHON SOURCE LINES 76-78 Wrap as a Model with steel material ----------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 78-85 .. code-block:: Python m = femorph_solver.Model.from_grid(grid) m.assign( ELEMENTS.HEX8, material={"EX": E, "PRXY": NU, "DENS": RHO}, ) .. GENERATED FROM PYTHON SOURCE LINES 86-92 Identify low / high cyclic-face nodes by angular position --------------------------------------------------------- The base sector spans :math:`\theta \in [0, \alpha]`. Nodes at :math:`\theta = 0` form the **low** face; nodes at :math:`\theta = \alpha` form the **high** face. .. GENERATED FROM PYTHON SOURCE LINES 92-114 .. code-block:: Python pts = np.asarray(grid.points) node_nums = np.asarray(grid.point_data["ansys_node_num"], dtype=np.int64) angle = np.arctan2(pts[:, 1], pts[:, 0]) tol = 1e-6 low_mask = np.abs(angle) < tol high_mask = np.abs(angle - ALPHA) < tol low_nodes = node_nums[low_mask] high_nodes = node_nums[high_mask] # Pair the two faces by (radius, z) coordinates so the cyclic # constraint links matching points across the sector boundary. low_pts = pts[low_mask] high_pts = pts[high_mask] low_radii_z = np.column_stack([np.linalg.norm(low_pts[:, :2], axis=1), low_pts[:, 2]]) high_radii_z = np.column_stack([np.linalg.norm(high_pts[:, :2], axis=1), high_pts[:, 2]]) order_low = np.lexsort(low_radii_z.T[::-1]) order_high = np.lexsort(high_radii_z.T[::-1]) low_nodes_paired = low_nodes[order_low] high_nodes_paired = high_nodes[order_high] print(f"cyclic faces: {len(low_nodes_paired)} paired nodes per face") .. rst-class:: sphx-glr-script-out .. code-block:: none cyclic faces: 10 paired nodes per face .. GENERATED FROM PYTHON SOURCE LINES 115-122 Run the harmonic sweep ---------------------- Solve the cyclic-symmetry modal problem on the base sector at every harmonic :math:`k = 0, 1, \ldots, N/2`. Each :class:`CyclicModalResult` carries the lowest n_modes eigenvalues for that harmonic. .. GENERATED FROM PYTHON SOURCE LINES 122-131 .. code-block:: Python n_modes = 8 results = m.cyclic_modal_solve( low_face=low_nodes_paired, high_face=high_nodes_paired, n_sectors=N_SECTORS, n_modes=n_modes, ) .. GENERATED FROM PYTHON SOURCE LINES 132-134 Tabulate the lowest non-rigid frequency at each harmonic -------------------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 134-147 .. code-block:: Python print() print(f"{'k':>3} {'f_min [Hz]':>12} {'mode shape':>16}") print(f"{'-':>3} {'-' * 12:>12} {'-' * 16:>16}") for r_k in results: freqs = np.asarray(r_k.frequency, dtype=np.float64) # Discard rigid-body modes (f < 1 Hz) and report the lowest # bending-like frequency instead. elastic = freqs[freqs > 1.0] f_min = float(elastic[0]) if elastic.size else float("nan") mode_kind = "rigid (k=0)" if r_k.harmonic_index == 0 else f"k={r_k.harmonic_index} pair" print(f"{r_k.harmonic_index:>3} {f_min:12.2f} {mode_kind:>16}") .. rst-class:: sphx-glr-script-out .. code-block:: none k f_min [Hz] mode shape - ------------ ---------------- 0 328.21 rigid (k=0) 1 438.02 k=1 pair 2 76.12 k=2 pair 3 218.20 k=3 pair 4 421.23 k=4 pair 5 684.49 k=5 pair 6 1008.36 k=6 pair .. GENERATED FROM PYTHON SOURCE LINES 148-154 Render the lowest non-rigid mode at the most interesting harmonic ----------------------------------------------------------------- k = 2 typically yields the lowest non-rigid mode on an annular disk — the classical "diametrical-2" pattern with two stationary radial nodes on the disc's rim. .. GENERATED FROM PYTHON SOURCE LINES 154-179 .. code-block:: Python target_k = 2 chosen = next((r for r in results if r.harmonic_index == target_k), results[1]) phi = np.asarray(chosen.mode_shapes)[:, 0].reshape(-1, 3) # k > 0 mode shapes are complex; take the real part for the # stationary snapshot. The travelling-wave companion lives in # the imaginary part — see example_cyclic_modes_bladed_rotor.py. phi_real = np.real(phi) warp = grid.copy() warp.points = grid.points + (T_DISK * 4.0 / np.max(np.abs(phi_real))) * phi_real warp["UZ"] = phi_real[:, 2] f_chosen = float(np.asarray(chosen.frequency)[0]) plotter = pv.Plotter(off_screen=True, window_size=(720, 480)) plotter.add_mesh(grid, style="wireframe", color="grey", opacity=0.4) plotter.add_mesh( warp, scalars="UZ", cmap="coolwarm", show_edges=False, scalar_bar_args={"title": f"k = {chosen.harmonic_index} mode 1 — f = {f_chosen:.1f} Hz"}, ) plotter.view_isometric() plotter.camera.zoom(1.05) plotter.show() .. image-sg:: /gallery/analyses/cyclic/images/sphx_glr_example_cyclic_modes_disk_001.png :alt: example cyclic modes disk :srcset: /gallery/analyses/cyclic/images/sphx_glr_example_cyclic_modes_disk_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.357 seconds) .. _sphx_glr_download_gallery_analyses_cyclic_example_cyclic_modes_disk.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_cyclic_modes_disk.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_cyclic_modes_disk.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: example_cyclic_modes_disk.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_