Modal analysis (free vibration) ================================ When to use ----------- Modal analysis recovers the **natural frequencies and mode shapes** of an unconstrained or constrained linear-elastic structure. Use it when: - You need to know whether the structure's natural frequencies fall inside any operating-band excitation (the Campbell-diagram clearance check). - You're sizing modes for a downstream **response-spectrum** or **modal-superposition transient** analysis (cumulative effective mass ≥ 90 % per direction is the standard gate). - You want a pre-stress baseline for buckling or a sanity check for BC completeness (the rigid-body-mode diagnostic in :doc:`/best-practices/bc-pitfalls`). For static deflection under fixed loads, use :doc:`linear-static`. For frequency-domain forced response, use :doc:`/user-guide/solving/harmonic`. For time-domain transient response, use :doc:`/user-guide/solving/transient`. Boundary conditions and loads ----------------------------- - **Dirichlet constraints** via :meth:`Model.fix(nodes=..., dof=...) ` — identical to :doc:`linear-static`. Most modal solves are run on a fully-constrained structure; for free-free modal analysis (a structure floating in space — a spacecraft, an unconstrained rod), see the ``sigma`` discussion below. - **No applied loads.** Modal analysis is the *unforced* free-vibration response — nodal-force / pressure / gravity are ignored. - **Pre-stress** (geometric stiffness from a prior static load) is on the planned roadmap for buckling / pre-stressed modal; not yet shipped. The math (one paragraph) ------------------------ Modal analysis solves the **generalised symmetric-definite eigenvalue problem** (:doc:`/reference/theory/eigenproblem`) .. math:: \mathbf{K}\, \boldsymbol{\phi}_i = \omega_i^{2}\, \mathbf{M}\, \boldsymbol{\phi}_i, returning eigenpairs :math:`(\omega_i^{2}, \boldsymbol{\phi}_i)` where :math:`\omega_i` is the angular natural frequency in rad/s and :math:`\boldsymbol{\phi}_i` is the mode shape. femorph-solver dispatches the lowest ``n_modes`` via shift-invert Lanczos (``scipy.sparse.linalg.eigsh``); see the theory chapter for the full algorithm. Running the solve ----------------- The high-level path: .. code-block:: python res = m.solve_modal(n_modes=12) Returns a :class:`~femorph_solver.solvers.modal.ModalResult`: - ``res.frequency`` — ``(n_modes,)`` natural frequencies in Hz, ascending. - ``res.omega_sq`` — ``(n_modes,)`` squared angular frequencies in (rad/s)². - ``res.mode_shapes`` — ``(n_dof, n_modes)`` mass- orthonormalised eigenvectors. Mode shapes are mass-orthonormalised by default (:math:`\boldsymbol{\phi}_i^{T}\, \mathbf{M}\, \boldsymbol{\phi}_j = \delta_{ij}`), which means post-processing recipes (modal participation, mode-shape animation) work directly without re-normalising. Free-free analysis ------------------ A free-free structure has six rigid-body modes at :math:`\omega = 0`, which makes :math:`\mathbf{K}` rank-deficient. Pass a non-zero shift to move the spectrum off the zero-eigenvalue ridge: .. code-block:: python res = m.solve_modal(n_modes=10, sigma=-1.0) The shift-invert path then factors :math:`(\mathbf{K} - \sigma \mathbf{M})` which is positive-definite for any :math:`\sigma < \omega_1^{2}`. See :ref:`sphx_glr_gallery_verification_example_verify_free_free_rod_nf.py` for the canonical free-free pattern. Targeting a frequency band -------------------------- For an off-zero shift (modes near a target frequency rather than the lowest), set ``sigma`` to the squared angular frequency of interest: .. code-block:: python # Modes near 2 kHz: sigma = (2 * np.pi * 2_000.0) ** 2 res = m.solve_modal(n_modes=10, sigma=sigma) Post-processing recipes ----------------------- - **Modal participation factors + effective mass** — :ref:`sphx_glr_gallery_post-processing_example_modal_participation.py`. The recipe used by Tutorial 2 to pick ``n_modes`` for a 90 %-cumulative-mass response-spectrum analysis. - **Mode-shape rendering** — :ref:`sphx_glr_gallery_post-processing_example_mode_shape_plot.py`. - **Mode-shape animation** (one mode swept over a period) — :ref:`sphx_glr_gallery_post-processing_example_mode_shape_animation.py`. Common gotchas -------------- - **Surviving rigid-body modes** when the structure should be constrained. Run the :doc:`/best-practices/bc-pitfalls` rigid-body-mode check on the unloaded model — every near-zero frequency is an unconstrained DOF. - **Bi-degenerate transverse pairs** on symmetric cross- sections. A square-section cantilever has ``f_y = f_z`` for every transverse-bending mode; the eigensolver returns them as two modes at the same frequency. The translation-energy mode-classification recipe in Tutorial 2 distinguishes them. - **Stress-recovery vs displacement accuracy.** Nodal averaging stress on a recovered mode shape is meaningless — mode shapes carry no absolute scale. Use ``stress`` only on a static or transient solve, or on a participation- factor-scaled mode for response-spectrum analysis. End-to-end tutorial ------------------- Tutorial 2 walks the modal-survey workflow on a cantilever bracket: :ref:`sphx_glr_gallery_tutorials_tutorial_02_modal_survey.py`. Verification examples that exercise this analysis type: - :ref:`sphx_glr_gallery_verification_example_verify_cantilever_nf.py` — first cantilever NF - :ref:`sphx_glr_gallery_verification_example_verify_cantilever_higher_modes.py` — first three cantilever modes - :ref:`sphx_glr_gallery_verification_example_verify_cc_beam_modes.py` — clamped-clamped beam modes - :ref:`sphx_glr_gallery_verification_example_verify_ss_beam_modes.py` — simply-supported beam modes - :ref:`sphx_glr_gallery_verification_example_verify_axial_rod_nf.py` — fixed-free axial rod - :ref:`sphx_glr_gallery_verification_example_verify_free_free_rod_nf.py` — free-free axial rod - :ref:`sphx_glr_gallery_verification_example_verify_ss_plate_modes.py` — simply-supported plate Solver mechanics + performance ------------------------------ For shift-invert Lanczos algorithm details, dense vs sparse path selection, eigen-solver backend tradeoffs, and the ``tol`` parameter, see :doc:`modal` — that page is the solver-engineering deep-dive on the same ``Model.solve_modal()`` invocation.