.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/analyses/harmonic/example_sdof_frf.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_harmonic_example_sdof_frf.py: SDOF receptance — closed form vs modal superposition ==================================================== A single-DOF mass-spring-damper has the analytical receptance .. math:: H(\omega) \;=\; \frac{1}{k - m\,\omega^{2} + i\,\omega\,c} This example drives the same SDOF system through :func:`~femorph_solver.solvers.harmonic.solve_harmonic_modal` and overlays the modal-superposition output on the closed form. With a single retained mode, the superposition reproduces the closed form to machine precision — the curves overlap exactly. Tuning the example ------------------ * The eigenvalue ``ω_n = √(k / m)`` sets the resonance peak. * The damping ratio ``ζ = c / (2 √(k m))`` sets the height. * Reducing damping makes the peak taller; the closed-form and modal curves still match. References ---------- * Ewins, D. J. *Modal Testing: Theory, Practice and Application*, 2nd ed., Research Studies Press (2000), §2.3 — receptance formulation used throughout this example. * Craig, R. R. and Kurdila, A. J. *Fundamentals of Structural Dynamics*, 2nd ed., Wiley (2006), §11.4 — modal superposition for forced harmonic response. .. GENERATED FROM PYTHON SOURCE LINES 35-44 .. code-block:: Python from __future__ import annotations import matplotlib.pyplot as plt import numpy as np import scipy.sparse as sp from femorph_solver.solvers.harmonic import solve_harmonic_modal .. GENERATED FROM PYTHON SOURCE LINES 45-47 1-DOF mass-spring-damper ------------------------ .. GENERATED FROM PYTHON SOURCE LINES 47-59 .. code-block:: Python m = 1.0 # kg k = 4.0e4 # N/m → ω_n = 200 rad/s, f_n ≈ 31.83 Hz zeta = 0.02 c = 2.0 * zeta * np.sqrt(k * m) omega_n = np.sqrt(k / m) f_n_hz = omega_n / (2.0 * np.pi) print(f"Resonance: ω_n = {omega_n:.2f} rad/s f_n = {f_n_hz:.3f} Hz ζ = {zeta:.3f}") K = sp.csr_array(np.array([[k]], dtype=np.float64)) M = sp.csr_array(np.array([[m]], dtype=np.float64)) .. rst-class:: sphx-glr-script-out .. code-block:: none Resonance: ω_n = 200.00 rad/s f_n = 31.831 Hz ζ = 0.020 .. GENERATED FROM PYTHON SOURCE LINES 60-64 Excitation sweep ---------------- Span 0.1 ω_n through 3 ω_n on a log grid so the resonance peak is well-resolved without piling samples at the high-frequency tail. .. GENERATED FROM PYTHON SOURCE LINES 64-66 .. code-block:: Python freq = np.geomspace(0.1 * f_n_hz, 3.0 * f_n_hz, 400) .. GENERATED FROM PYTHON SOURCE LINES 67-69 Modal-superposition response ---------------------------- .. GENERATED FROM PYTHON SOURCE LINES 69-79 .. code-block:: Python F = np.array([1.0], dtype=np.complex128) # unit forcing on the single DOF res = solve_harmonic_modal( K, M, F, frequencies=freq, n_modes=1, modal_damping_ratio=zeta, ) .. GENERATED FROM PYTHON SOURCE LINES 80-82 Closed-form receptance ---------------------- .. GENERATED FROM PYTHON SOURCE LINES 82-93 .. code-block:: Python omega = 2.0 * np.pi * freq H_closed = 1.0 / (k - m * omega**2 + 1j * omega * c) # Sanity-check: with one retained mode the two should match to floating # point precision. np.testing.assert_allclose(res.displacement[:, 0], H_closed, atol=1e-12, rtol=1e-12) print( "Modal-superposition matches closed form to " f"max |Δ| = {np.max(np.abs(res.displacement[:, 0] - H_closed)):.2e}" ) .. rst-class:: sphx-glr-script-out .. code-block:: none Modal-superposition matches closed form to max |Δ| = 0.00e+00 .. GENERATED FROM PYTHON SOURCE LINES 94-96 Overlay the magnitude / phase ----------------------------- .. GENERATED FROM PYTHON SOURCE LINES 96-118 .. code-block:: Python fig, (ax_mag, ax_phase) = plt.subplots(2, 1, sharex=True, figsize=(8, 6)) ax_mag.loglog(freq, np.abs(H_closed), "k-", lw=2.5, label="closed form") ax_mag.loglog(freq, np.abs(res.displacement[:, 0]), "C1--", lw=1.5, label="modal superposition") ax_mag.axvline(f_n_hz, color="grey", linestyle=":", alpha=0.6, label=f"f_n = {f_n_hz:.2f} Hz") ax_mag.set_ylabel("|H(ω)| [m / N]") ax_mag.set_title("SDOF receptance") ax_mag.grid(True, which="both", alpha=0.3) ax_mag.legend() phase = np.angle(res.displacement[:, 0], deg=True) phase_closed = np.angle(H_closed, deg=True) ax_phase.semilogx(freq, phase_closed, "k-", lw=2.5, label="closed form") ax_phase.semilogx(freq, phase, "C1--", lw=1.5, label="modal superposition") ax_phase.axvline(f_n_hz, color="grey", linestyle=":", alpha=0.6) ax_phase.set_xlabel("frequency [Hz]") ax_phase.set_ylabel("phase [deg]") ax_phase.grid(True, which="both", alpha=0.3) ax_phase.legend() fig.tight_layout() plt.show() .. image-sg:: /gallery/analyses/harmonic/images/sphx_glr_example_sdof_frf_001.png :alt: SDOF receptance :srcset: /gallery/analyses/harmonic/images/sphx_glr_example_sdof_frf_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.423 seconds) .. _sphx_glr_download_gallery_analyses_harmonic_example_sdof_frf.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_sdof_frf.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_sdof_frf.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: example_sdof_frf.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_