Algebraic Harmonic Balance Method

This project provides a Python 3.10 implementation of the algebraic harmonic balance method (algebraic HBM, aHBM) as proposed here.

Theoretic background

We are considering second order ordinary differential equations (ODEs) with polynomial coefficients in the state x : \mathbb T \subset \mathbb R \to \mathbb R, that is ODEs of the form

r(t,x;u) = \rho x''(t) + \delta x'(t) + \sum_{i=1}^q \alpha_i x^i(t) - u(t) = 0 \,, \quad u(t) = \hat u \cos(\Omega t) \,.

The idea of the HBM is to yield approximations x_n(t) = c_0 + \sum_{i=1}^n c_{2i-1} \cos(i \Omega t) + c_{2i} \sin(i \Omega t) of stationary periodic solutions x of the ODE. Given a excitation frequency \Omega, the algebraic HBM of order n yields a system of multivariate polynomials R_i, i=0,1,\ldots,2n, in the variables c_0,c_1,\ldots,c_{2n} that solve the algebraic system

F_n(\mathbf c; \Omega) = [R_i(\mathbf c; \Omega)]_{i=0}^{2n} = 0

where \mathbf c = [c_0,c_1,\ldots,c_{2n}] \in \mathbb R^{2n+1}. A solution \mathbf c \leftrightarrow x_n of F_n(\mathbf c; \Omega) = 0 is also a solution of the (original) HBM defining system of integral equations

\langle r(x_n), \phi_j\rangle = \frac{1}{T} \int_0^T r(t,x_n(t)) \phi_j(t) \mathrm d t \,, \quad j = 0,1,\ldots,2n \,,

with basis functions \phi_0(t) = 1, \phi_{2i-1}(t) = \cos(i \Omega t) and \phi_{2i}(t) = \sin(i \Omega t), i=1,\ldots,n. Note that building and evaluating F_n does not require the computation of integrals as e.g. in the classical or Alternating Frequency-Time HBM [Krack, Gross].

Minimal working example

A MWE can be found here. Below is a slighlty more extensive version of said example.

Required imports:

from ode import ODE_2nd_Order_Poly_Coeffs
from harmonic_balance_method import Algebraic_HBM

Define the classical softening Duffing oscillator

r(t,x) = x''(t) + 0.4 x'(t) + -0.4 x^3(t) - 0.3 \cos(\Omega t) = 0

as a second order ODE with polynomial coefficients via

ode = ODE_2nd_Order_Poly_Coeffs(mass=1, damping=.4, stiffness=1, excitation=(0,.3), monomials={3: -.4})

Then, initialize the algebraic HBM for ansatz order n.

HBM = Algebraic_HBM(ODE=ode, order=n)

Now generate the multivariate polynomials that define the algebraic equation system of the algebraic HBM.

HBM.generate_multivariate_polynomials()

Compile multivariate polynomials into excecutable functions F_n and \mathrm DF_n=(\frac{\mathrm d F_n}{\mathrm d \mathbf c}, \frac{\mathrm d F_n}{\mathrm d a}).

F, DF = HBM.compile()

The system F_n and its Jacobian \mathrm d F_n may now be used to perform a bifurcation analysis of the system. For example, computing the frequency response of the system for different initial guesses [Dänschel, Lentz, von Wagner].