Harmonic / frequency-response analysis
======================================

When to use
-----------

Harmonic analysis recovers the **steady-state response of a
structure to sinusoidal forcing at a fixed frequency**, or
across a frequency sweep.  Use it when:

- You have rotating-equipment harmonics (motors, pumps,
  compressors) exciting at known frequencies and want the
  response amplitude / phase at each.
- Your loading is a deterministic sinusoid (mooring loads at
  wave frequency, NVH input at engine orders) and you want
  the frequency-by-frequency gain.
- You need a transfer-function curve (FRF) for control-system
  / vibration-isolation design.

For random-vibration / PSD inputs use a future
random-vibration solver path.  For non-stationary or
shock-type loading use :doc:`/user-guide/solving/transient`.  For free-vibration
modes alone use :doc:`modal`.

Boundary conditions and loads
-----------------------------

- **Dirichlet constraints** — same as :doc:`linear-static`.
- **Complex-valued nodal forces** via the harmonic API
  surface.  Each force has a magnitude *and* phase; the
  solver returns the complex displacement.
- **Damping** — Rayleigh damping
  (:math:`\mathbf{C} = \alpha \mathbf{M} + \beta
  \mathbf{K}`) at the model level, picked via
  :meth:`Model.rayleigh_damping <femorph_solver.Model.rayleigh_damping>`
  or supplied directly to the solve.

The math (one paragraph)
------------------------

Harmonic analysis solves the complex linear system

.. math::

    \bigl(-\omega^{2}\, \mathbf{M}
        + i\, \omega\, \mathbf{C}
        + \mathbf{K}\bigr)\,
    \mathbf{u}(\omega) = \mathbf{f}(\omega)

at each forcing angular frequency :math:`\omega = 2\pi f`.
The complex displacement vector :math:`\mathbf{u}(\omega)`
encodes both amplitude and phase relative to the forcing.
Two solution paths: **direct** (factor the complex stiffness
at each :math:`\omega` — accurate, expensive) and **modal
superposition** (project onto pre-computed mode shapes —
much cheaper, valid when :math:`f \ll f_\text{nyquist of mode set}`).

Running the solve
-----------------

.. code-block:: python

    freqs = np.linspace(0.0, 1000.0, 201)  # 0 to 1 kHz, 5 Hz spacing
    res = m.solve_harmonic(frequencies=freqs)

Returns a :class:`~femorph_solver.solvers.harmonic.HarmonicResult`:

- ``res.frequencies`` — ``(n_freq,)`` forcing frequencies [Hz].
- ``res.displacement`` — ``(n_dof, n_freq)`` complex; row k
  column j is :math:`u_k(\omega_j)`.
- ``res.diagnostics`` — solver-statistics dict (factor-cost,
  iteration count, etc.).

For a single-DOF / single-frequency lookup, see the focused
recipe at
:ref:`sphx_glr_gallery_analyses_harmonic_example_sdof_frf.py`.

Post-processing recipes
-----------------------

- **Frequency-response curves** — index ``res.displacement``
  by the DOF of interest, take ``abs()`` for magnitude and
  ``np.angle(...)`` for phase, plot vs ``res.frequencies``.
- **Resonance detection** — peaks in the magnitude curve
  cluster at the natural frequencies from
  :doc:`modal`.  Cross-check the location of each peak.
- **Mode-superposition mode** — pass ``method="modal"`` and
  a pre-computed ``ModalResult``; the solve projects each
  forcing onto modal coordinates instead of solving the full
  complex system.  Big speedup at the cost of mode-set
  truncation error.

Common gotchas
--------------

- **Excitation near a resonance with no damping** produces
  unbounded amplitude — the factorisation becomes
  ill-conditioned right at :math:`f = f_n`.  Always
  configure Rayleigh (or modal) damping before sweeping
  through a resonance.
- **Frequency-step too coarse** misses sharp resonance peaks.
  At low damping ratios (ζ < 1 %) the half-power bandwidth
  narrows below 1 Hz; resolve with ``∆f`` ≤ ``f_n × ζ``.
- **Mode-superposition truncation.**  Including only the
  lowest 12 modes when the forcing reaches 5 kHz misses the
  contribution of higher modes; use the modal-participation
  curve to choose ``n_modes``.

Verification + gallery
----------------------

- :ref:`sphx_glr_gallery_analyses_harmonic_example_sdof_frf.py` —
  single-DOF spring-mass-damper FRF vs textbook closed form.

Solver mechanics + performance
------------------------------

For complex-stiffness assembly, mode-superposition
mechanics, and frequency-sweep batching strategies, see
:doc:`/user-guide/solving/harmonic` — that page is the solver-engineering
deep-dive on the same ``Model.solve_harmonic()``
invocation.