.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/analyses/cyclic/example_rotor_sector.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_rotor_sector.py: Cyclic-symmetry modal on a rotor sector ======================================= A rotor with :math:`N` identical sectors spanning 360° has a spectrum that factors by *harmonic index* :math:`k \in \{0, 1, \dots, N/2\}`. Analysing one base sector plus the constraint .. math:: u_\text{high} = e^{i k \alpha}\, R(\alpha)\, u_\text{low}, \qquad \alpha = 2\pi / N, recovers the full-rotor spectrum at a fraction of the DOF count. This example builds a simple 4-sector rotor from a pyvista annular sector, then runs :func:`~femorph_solver.solvers.cyclic.solve_cyclic_modal` across all harmonic indices. .. GENERATED FROM PYTHON SOURCE LINES 21-31 .. code-block:: Python from __future__ import annotations import numpy as np import pyvista as pv from scipy.spatial import cKDTree import femorph_solver from femorph_solver.solvers.cyclic import solve_cyclic_modal .. GENERATED FROM PYTHON SOURCE LINES 32-34 Build a 90° annular sector -------------------------- .. GENERATED FROM PYTHON SOURCE LINES 34-50 .. code-block:: Python N_SECTORS = 4 ALPHA = 2.0 * np.pi / N_SECTORS N_T, N_R, N_Z = 12, 6, 3 theta = np.linspace(0.0, ALPHA, N_T + 1) r = np.linspace(0.8, 1.0, N_R + 1) z = np.linspace(0.0, 0.1, N_Z + 1) # (r, θ, z) order gives a CCW-from-above hex winding in the (r, θ) # plane — swapping r and θ would flip the winding inside-out. 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() .. GENERATED FROM PYTHON SOURCE LINES 51-53 Build the femorph-solver model and assemble -------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 53-62 .. code-block:: Python m = femorph_solver.Model.from_grid(grid) m.et(1, "SOLID185") m.mp("EX", 1, 2.1e11) m.mp("PRXY", 1, 0.30) m.mp("DENS", 1, 7850.0) K = m.stiffness_matrix() M = m.mass_matrix() .. rst-class:: sphx-glr-script-out .. code-block:: none /home/runner/_work/solver/solver/examples/analyses/cyclic/example_rotor_sector.py:54: DeprecationWarning: Model.et(...) is a MAPDL-dialect shortcut and has moved off the Model public surface. Use `APDL(model).et(et_id, name)` for line-by-line APDL deck porting, or the native `Model.assign("HEX8", material)` for new code. m.et(1, "SOLID185") /home/runner/_work/solver/solver/examples/analyses/cyclic/example_rotor_sector.py:55: DeprecationWarning: Model.mp(...) is a MAPDL-dialect shortcut and has moved off the Model public surface. Use `APDL(model).mp(prop, mat_id, value)` for line-by-line APDL deck porting, or the native `Model.assign(element, {prop: value, ...})` for new code. m.mp("EX", 1, 2.1e11) /home/runner/_work/solver/solver/examples/analyses/cyclic/example_rotor_sector.py:56: DeprecationWarning: Model.mp(...) is a MAPDL-dialect shortcut and has moved off the Model public surface. Use `APDL(model).mp(prop, mat_id, value)` for line-by-line APDL deck porting, or the native `Model.assign(element, {prop: value, ...})` for new code. m.mp("PRXY", 1, 0.30) /home/runner/_work/solver/solver/examples/analyses/cyclic/example_rotor_sector.py:57: DeprecationWarning: Model.mp(...) is a MAPDL-dialect shortcut and has moved off the Model public surface. Use `APDL(model).mp(prop, mat_id, value)` for line-by-line APDL deck porting, or the native `Model.assign(element, {prop: value, ...})` for new code. m.mp("DENS", 1, 7850.0) .. GENERATED FROM PYTHON SOURCE LINES 63-65 Identify cyclic-face node pairs ------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 65-87 .. code-block:: Python pts = np.asarray(grid.points) low = np.where(np.abs(np.arctan2(pts[:, 1], pts[:, 0])) < 1e-6)[0] high = np.where(np.abs(np.arctan2(pts[:, 1], pts[:, 0]) - ALPHA) < 1e-6)[0] c, s = np.cos(ALPHA), np.sin(ALPHA) R_sector = np.array([[c, -s, 0.0], [s, c, 0.0], [0.0, 0.0, 1.0]]) low_rot = pts[low] @ R_sector.T _, pair_idx = cKDTree(pts[high]).query(low_rot) high_paired = high[pair_idx] dof_map = m.dof_map() nn_col = dof_map[:, 0] dof_col = dof_map[:, 1] node_nums = np.asarray(grid.point_data["ansys_node_num"]) low_dofs = np.empty((len(low), 3), dtype=np.int64) high_dofs = np.empty((len(low), 3), dtype=np.int64) for i, (lo, hi) in enumerate(zip(low, high_paired)): nn_lo, nn_hi = int(node_nums[lo]), int(node_nums[hi]) for d in range(3): low_dofs[i, d] = int(np.where((nn_col == nn_lo) & (dof_col == d))[0][0]) high_dofs[i, d] = int(np.where((nn_col == nn_hi) & (dof_col == d))[0][0]) .. GENERATED FROM PYTHON SOURCE LINES 88-90 Solve for every harmonic index ------------------------------ .. GENERATED FROM PYTHON SOURCE LINES 90-105 .. code-block:: Python results = solve_cyclic_modal( K, M, low_dofs, high_dofs, n_sectors=N_SECTORS, n_modes=4, pair_rotation=R_sector, ) print(f"Cyclic modal on {N_SECTORS}-sector rotor:") for res in results: fs = ", ".join(f"{f:.2f}" for f in res.frequency[:4]) print(f" k = {res.harmonic_index}: f [Hz] = [{fs}]") .. rst-class:: sphx-glr-script-out .. code-block:: none Cyclic modal on 4-sector rotor: k = 0: f [Hz] = [0.00, 0.00, 503.32, 503.32] k = 1: f [Hz] = [0.00, 0.00, 262.09, 452.92] k = 2: f [Hz] = [91.75, 91.75, 164.68, 164.68] .. GENERATED FROM PYTHON SOURCE LINES 106-108 Visualise mode 0 of k = 1 ------------------------- .. GENERATED FROM PYTHON SOURCE LINES 108-121 .. code-block:: Python r1 = next(res for res in results if res.harmonic_index == 1) phi = r1.mode_shapes[:, 0].reshape(-1, 3) grid.point_data["mode_k1"] = np.abs(phi) plotter = pv.Plotter(off_screen=True) plotter.add_mesh( grid.warp_by_vector("mode_k1", factor=0.05 / (np.max(np.abs(phi)) + 1e-30)), scalars=np.linalg.norm(phi, axis=1), cmap="viridis", scalar_bar_args={"title": f"|mode| (k=1, f={r1.frequency[0]:.1f} Hz)"}, ) plotter.add_axes() plotter.show() .. image-sg:: /gallery/analyses/cyclic/images/sphx_glr_example_rotor_sector_001.png :alt: example rotor sector :srcset: /gallery/analyses/cyclic/images/sphx_glr_example_rotor_sector_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.506 seconds) .. _sphx_glr_download_gallery_analyses_cyclic_example_rotor_sector.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_rotor_sector.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_rotor_sector.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: example_rotor_sector.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_