Explicit Third-Order Model Reduction Formulas for General Nonlinear Mechanical Systems

TL;DR

This paper derives explicit third-order reduced models for nonlinear mechanical systems using invariant manifold theories, enabling straightforward analysis of system characteristics without complex coordinate transformations.

Contribution

It provides closed-form, cubic-order reduction formulas for general nonlinear systems based on LSM and SSM theories, simplifying model analysis and parameter extraction.

Findings

Reduced models accurately predict backbone and forced-response curves.

Formulas are explicit and depend only on physical and modal parameters.

Validated on complex systems like a nonlinear Timoshenko beam.

Abstract

For general nonlinear mechanical systems, we derive closed-form, reduced-order models up to cubic order based on rigorous invariant manifold results. For conservative systems, the reduction is based on Lyapunov Subcenter Manifold (LSM) theory, whereas for damped-forced systems, we use Spectral Submanifold (SSM) theory. To evaluate our explicit formulas for the reduced model, no coordinate changes are required beyond an initial linear one. The reduced-order models we derive are simple and depend only on physical and modal parameters, allowing us to extract fundamental characteristics, such as backbone curves and forced-response curves, of multi-degree-of-freedom mechanical systems. To numerically verify the accuracy of the reduced models, we test the reduction formulas on several mechanical systems, including a higher-dimensional nonlinear Timoshenko beam.

Tables

Table 1. Table 1: Hardening coefficients of different model reduction methods.

Method	${\omega }_{\_}2$
LSM
Normal forms
MDs

Table 2. Table 2: Geometric and material parameters.

Parameter	Value
$L$	$1200\mathrm{mm}$
$h$	$40\mathrm{mm}$
$b$	$40\mathrm{mm}$
$\rho$	$7850\cdot {10}^{-9}\mathrm{kg}{\mathrm{mm}}^{-3}$
$E$	$90\mathrm{GPa}$
$G$	$34.6\mathrm{GPa}$
$\eta$	$13.4\mathrm{MPa}s$
$\mu$	$8.3\mathrm{MPa}s$

Equations

$$\ddot{x}+{\omega }^{2}x+P(x,\mathbf{y})$$

$$\mathbf{M}\ddot{\mathbf{y}}+\mathbf{Ky}+\mathbf{Q}(x,\mathbf{y})$$

$$\omega =\sqrt{\frac{k}{m}},{\mathbf{\Omega }}_{\_}^2\text{p}={\mathbf{M}}^{-1}\mathbf{K},$$

$$P(x,\mathbf{y})=\sum _{\_}{j+|\mathbf{k}|=2}^l{p}_{\_}j\mathbf{k}{x}^j{y}_{\_}{1}^{k_{\_}1}\dots y_{\_}{n}^{k_{\_}n},\mathbf{Q}(x,\mathbf{y})=\sum _{\_}{j+|\mathbf{k}|=2}^l{\mathbf{q}}_{\_}j\mathbf{k}{x}^j{y}_{\_}{1}^{k_{\_}1}\dots y_{\_}{n}^{k_{\_}n},$$

$$p_{\_}j\mathbf{k}\in \mathbb{R},{\mathbf{q}}_{\_}j\mathbf{k}\in {\mathbb{R}}^n,\mathbf{k}\in {\mathbb{N}}^n.$$

$$0=(0,\dots ,0)\in {\mathbb{N}}^n,{\mathbf{e}}_{\_}j=(0,\dots ,0,\underset{j^{\text{th}}}{1},0,\dots ,0)\in {\mathbb{N}}^n,$$

$${\mathbf{p}}_{\_}j\mathbf{I}={\left[p_{\_}j{\mathbf{e}}_{\_}1,\dots ,p_{\_}j{\mathbf{e}}_{\_}n\right]}^{\top }\in {\mathbb{R}}^n,{\mathbf{Q}}_{\_}j\mathbf{I}=\left[{\mathbf{q}}_{\_}j{\mathbf{e}}_{\_}1,\dots ,{\mathbf{q}}_{\_}j{\mathbf{e}}_{\_}n\right]\in {\mathbb{R}}^{n\times n}.$$

$$\frac{{\omega }_{\_}i}{\omega }\notin \mathbb{Z},i=1,\dots ,n,$$

$$\mathcal{W}=\{(x,\dot{x},\mathbf{y},\dot{\mathbf{y}})\in {\mathbb{R}}^{2n+2},\mathbf{y}=\mathbf{w}(x,\dot{x}),\dot{\mathbf{y}}=\dot{\mathbf{w}}(x,\dot{x})\}$$

$$\mathbf{y}=\mathbf{w}(x,\dot{x})=\sum _{\_}i+j\in \{2,3\}{\mathbf{w}}_{\_}ij{x}^i{\dot{x}}^j+\mathcal{O}(|(x,\dot{x}){|}^4).$$

$$\ddot{x}+{\omega }^{2}x+P\left(x,\mathbf{w}\left(x,\dot{x}\right)\right)+\mathcal{O}\left(|\left(x,\dot{x}\right){|}^4\right)=0.$$

$$\ddot{x}+{\omega }^{2}x+p_{\_}20x^2+\left(p_{\_}30+⟨{\mathbf{p}}_{\_}1\mathbf{I},{\mathbf{w}}_{\_}20⟩\right)x^3+⟨{\mathbf{p}}_{\_}1\mathbf{I},{\mathbf{w}}_{\_}02⟩x{\dot{x}}^2+\mathcal{O}(|(x,\dot{x}){|}^4)=0,$$

$${\mathbf{w}}_{\_}20$$

$${\mathbf{w}}_{\_}02$$

$$\ddot{x}+{\omega }^{2}x+R(x,\mathbf{\eta })$$

$$\ddot{\mathbf{\eta }}+{\mathbf{\Omega }}^2\mathbf{\eta }+\mathbf{S}(x,\mathbf{\eta })$$

$$R(x,\mathbf{\eta })=\sum _{\_}{j+|\mathbf{k}|=2}^l{r}_{\_}j\mathbf{k}{x}^j{\eta }_{\_}{1}^{k_{\_}1}\dots {\eta }_{\_}{n}^{k_{\_}n},\mathbf{S}(x,\mathbf{\eta })=\sum _{\_}{j+|\mathbf{k}|=2}^l{\mathbf{s}}_{\_}j\mathbf{k}{x}^j{\eta }_{\_}{1}^{k_{\_}1}\dots {\eta }_{\_}{n}^{k_{\_}n},$$

$$r_{\_}j\mathbf{k}\in \mathbb{R},{\mathbf{s}}_{\_}j\mathbf{k}\in {\mathbb{R}}^n,\mathbf{k}\in {\mathbb{N}}^n.$$

$${\mathbf{r}}_{\_}j\mathbf{I}:={\left[r_{\_}j{\mathbf{e}}_{\_}1,\dots ,r_{\_}j{\mathbf{e}}_{\_}n\right]}^{\top }\in {\mathbb{R}}^n,{\mathbf{S}}_{\_}j\mathbf{I}:=\left[{\mathbf{s}}_{\_}j{\mathbf{e}}_{\_}1,\dots ,{\mathbf{s}}_{\_}j{\mathbf{e}}_{\_}n\right]\in {\mathbb{R}}^{n\times n}.$$

$$\mathbf{\eta }=\tilde{\mathbf{w}}(x,\dot{x})=\sum _{\_}i+j\in \{2,3\}{\tilde{\mathbf{w}}}_{\_}ij{x}^i{\dot{x}}^j+\mathcal{O}(|(x,\dot{x}){|}^4).$$

$$\ddot{x}+{\omega }^{2}x+R(x,\tilde{\mathbf{w}}(x,\dot{x}))+\mathcal{O}(|(x,\dot{x}){|}^4)=0.$$

$$\begin{array}{rl}& \ddot{x}+{\omega }^{2}x+r_{\_}20x^2+\left(r_{\_}30-⟨{\mathbf{r}}_{\_}1\mathbf{I},({\mathbf{\Omega }}^2{\mathbf{D}}_{\_}4{)}^{-1}{\mathbf{D}}_{\_}2{\mathbf{s}}_{\_}20⟩\right)x^3+2⟨{\mathbf{r}}_{\_}1\mathbf{I},({\mathbf{\Omega }}^2{\mathbf{D}}_{\_}4{)}^{-1}{\mathbf{s}}_{\_}20⟩x{\dot{x}}^2+\mathcal{O}(|(x,\dot{x}){|}^4)=0,\end{array}$$

$${\mathbf{D}}_{\_}2=({\mathbf{\Omega }}^2-2{\omega }^2\mathbf{I}),{\mathbf{D}}_{\_}4=({\mathbf{\Omega }}^2-4{\omega }^2\mathbf{I}).$$

$$r_{\_}20=0,⟨{\mathbf{r}}_{\_}1\mathbf{I},({\mathbf{\Omega }}^2{\mathbf{D}}_{\_}4{)}^{-1}{\mathbf{s}}_{\_}20⟩=0.$$

$$\ddot{x}+{\omega }^{2}x+r_{\_}30x^3+\mathcal{O}(x^4)=0.$$

$$\begin{array}{rl}\dot{x}=& y,\\ \dot{y}=& -{\omega }^{2}x-\alpha x^2-\beta x^3-\gamma x{y}^2.\end{array}$$

$$\omega ={\mathrm{e}}^{-\gamma x^2}dx\wedge dy,H(x,y)=\frac{1}{2}{e}^{\gamma x^2}{y}^2+\frac{{\omega }^2}{2\gamma }{e}^{\gamma x^2}+\int e^{\gamma x^2}\left(\alpha x^2+\beta x^3\right)dx.$$

$$D^{2}H(0,0)={\omega }^2>0,$$

$$\mathbf{\Theta }=-\frac{1}{2}{\partial }_{\_}\mathbf{\eta }\mathbf{S}(x,\mathbf{\eta }{)}^{-1}{{\partial }_{\_}^{2}xx\mathbf{S}(x,\mathbf{\eta })|}_{\_}x=0,\mathbf{\eta }=0,$$

$$\mathbf{S}(x,\mathbf{\eta })={\mathbf{\Omega }}^2\mathbf{\eta }+\sum _{\_}{j+|\mathbf{k}|=2}^3{\mathbf{s}}_{\_}j\mathbf{k}{x}^j{\eta }_{\_}{1}^{k_{\_}1}\dots {\eta }_{\_}{n}^{k_{\_}n}.$$

$$\mathbf{\Theta }=-{\mathbf{\Omega }}^{-2}{\mathbf{s}}_{\_}20,$$

$$\ddot{x}+{\omega }_{\_}{0}^{2}x+p_{\_}20x^2+\left(p_{\_}30-⟨{\mathbf{p}}_{\_}1\mathbf{I},\mathbf{\Phi }{\mathbf{\Omega }}^{-2}{\mathbf{s}}_{\_}20⟩\right)x^3=0.$$

$$q_{\_}i$$

Full text

1 Introduction

The dimensional reduction of nonlinear systems in structural dynamics has two underlying drivers. One is the necessity to decrease computational time, as computing a reduced-order model can be orders of magnitudes faster than solving the full system. The second driver is model identification, which is easier to carry out when a reduced number of parameters needs to be fitted to experimental or simulation data.

The appealing properties of modal decomposition and reduction in linear systems motivate the exploration of a similar decomposition and reduction for nonlinear systems. The latter can be accomplished in a mathematically exact fashion if the spectral subspaces of the linearized system can be continued as invariant manifolds of equal dimension, tangent to those spectral subspaces at a fixed point. Moreover, if these manifolds are unique and analytic, they can be sought in the form of Taylor expansions near the fixed point. As a consequence, the dynamics on these manifolds will serve as a mathematically exact reduced-order model that can be constructed up to any required order of accuracy.

In this work, we provide ready-to-use explicit formulas for third-order reduced models of autonomous and non-autonomous mechanical systems with general second- and third-order polynomial nonlinearities, without transforming the second-order equations of motion to a first-order form. For conservative systems, this third-order model reduction is performed onto Lyapunov Subcenter Manifolds (LSMs) [1], whereas for damped-forced systems, we perform the reduction onto Spectral Submanifolds (SSMs) [2], which are mathematically rigorous versions of the nonlinear normal modes proposed first by Shaw and Pierre [3]. To obtain higher-order LSM- and SSM-reduced models, one may use the automated numerical reduction procedure developed in [4].

Alternative approaches to nonlinear model reduction in mechanical systems include the method of normal forms ([5, 6, 7] ) and the method of modal derivatives [8]. Normal forms involve a recursive simplification of the full system via a sequence of nonlinear coordinate changes near a fixed point. Truncated at any finite order, the normal form admits invariant subspaces to which one may reduce the dynamics. It is, however, generally unclear from this approach what will happen to these invariant subspaces under the addition of the remaining, unnormalized terms. Furthermore, the reduced model is available in coordinates that are nonlinear functions of the original coordinates and hence their physical meaning is not immediate. Finally, the procedure assumes the knowledge of a full linear modal decomposition at the start, which is not necessarily available for large systems.

The method of normal forms has nevertheless been applied to equations of motion in second-order form including periodic forcing (cf. Neild and Wagg [6] and Neild et al. [7]). As an advantage, the method allows for the extraction of the steady-state response in addition to backbone curves in the unforced and undamped limit. As opposed to SSM theory, however, the method does not provide a mathematically exact reduced-order model and requires additional smallness assumptions on the nonlinear and damping terms.

In earlier work, Touzé and Amabili [5] propose a normal form procedure for the construction of reduced-order models for geometrically nonlinear vibrations of thin structures. They introduce the external forcing directly into the normal form along non-physical, curvilinear coordinates. In this work, we prove that at leading-order, the SSM-reduced system of a periodically forced mechanical system can be seen as the autonomous SSM-reduced system with the modal-participation factor (see Géradin and Rixen [9]) of the first mode added to the reduced system, therefore justifying the Touzé–Amabili approach. Our proof utilizes recent work by Breunung and Haller [10] on time-periodic SSMs in mechanical systems.

A more recent application of the method of normal forms for model-order reduction purposes can be found in the work of Denis et al. [11]. They observe that, in the absence of internal resonances, a reduced-order model for each mode can be approximated by a Duffing oscillator. The method requires full knowledge of the linear modal decomposition and excludes resonances up to any order in order to make the unnormalized tail of the asymptotic expansion arbitrarily small. In this work, we confirm this idea mathematically using LSM theory and show that for a conservative system with symmetric cubic nonlinearities, a third-order reduced model can indeed be sought as a Duffing oscillator, without using a modal transformation.

The other alternative to nonlinear model reduction is the method of modal derivatives [8], which assumes the existence of quadratic invariant manifolds in the configuration space (space of positions). While such a manifold can formally be sought, it generally does not actually exist in the configuration space but in the phase space. Accordingly, modal derivatives only give an accurate reduced-order model in slow-fast systems in which spectral submanifolds have a weak dependence on the velocities (see [12]).

In our present development of explicit, cubic reduction formulas, we circumvent the issues with normal forms and modal derivatives and provide generally applicable, explicit formulas for nonlinear reduced-order models. We also illustrate the use of our reduction formulas on different mechanical systems varying from simple to complex.

2 Conservative systems

2.1 Setup

We consider $(n+1)$-degree-of-freedom, undamped nonlinear autonomous mechanical systems of the form

[TABLE]

where $x$ is a scalar, $\mathbf{y}$ is an $n$-vector; $P$ and $\mathbf{Q}$ are nonlinear, polynomial functions of $x$ and $\mathbf{y}$ up to order $2\leq l\leq 3$; $\mathbf{M}$ and $\mathbf{K}$ are the mass and stiffness matrices corresponding to the $n$-degree-of-freedom generalized coordinate vector $\mathbf{y}$, while $m$ and $k$ are the modal mass and stiffness coefficients corresponding to $x$. We refer to the $x$ variable as the modeling variable and $\mathbf{y}$ as the vector of non-modeling variables. Obtaining the form (1)-(2) for a mechanical system requires the use of a linear transformation that decouples the modeling variable from the rest. To this end, only the mode shape associated with the $x$-mode has to be explicitly known; the remaining $n$ columns of the matrix of the linear transformation can be chosen as an arbitrary basis from the $n$-dimensional plane orthogonal to the known mode shape in $\mathbb{R}^{n+1}$.

We will use the notation

[TABLE]

where $\omega$ is the natural frequency belonging to the modeling mode $x$ with the imaginary pair of eigenvalues $\pm\mathrm{i}\omega$. The non-modeling modes belong to the purely imaginary eigenvalues $\pm\mathrm{i}\omega_{\_}1,\ldots,\pm\mathrm{i}\omega_{\_}n$, with $\omega_{\_}k$ denoting the $k^{\text{th}}$ natural frequency. Our notation for the coefficients of the $l^{\text{th}}$-order polynomials $P$ and $\mathbf{Q}$ is

[TABLE]

For the integer vectors $\mathbf{k}\in\mathbb{N}^{n}$, specific values we will often need are

[TABLE]

denoting the identically zero index vector and the $j^{\text{th}}$ unit vector, respectively. Additionally, we introduce the notation

[TABLE]

2.2 The LSM-reduced model in the general case

As shown by Kelley [1], under the non-resonance conditions

[TABLE]

a unique, analytic, two-dimensional invariant manifold

[TABLE]

exists for system (1)-(2), tangent to the two-dimensional subspace $\mathbf{y}=\dot{\mathbf{y}}=\mathbf{0}$ at the trivial fixed point $x=\dot{x}=0$, $\mathbf{y}=\dot{\mathbf{y}}=\mathbf{0}$. Often called the Lyapunov subcenter manifold (LSM), $\mathcal{W}$ is known to be filled with periodic orbits of (1)-(2). In the nonlinear vibrations literature, the amplitude-frequency plot of these periodic orbits is called the (conservative) backbone curve of system (1)-(2) (see, e.g., Londono et al. [13]).

Reducing the full dynamics to the invariant manifold $\mathcal{W}$ gives an exact two-dimensional model for the nonlinear dynamics associated with the $x$-mode. Such a model has apparently not been computed for the general class of systems (1)-(2) in the literature before. We will carry out this computation explicitly below.

As a first step, we derive a universally valid third-order approximation to $\mathcal{W}$ in the form

[TABLE]

We restrict (1)-(2) to the LSM by enslaving the $\mathbf{y}$ variable to $(x,\dot{x})$ in equation (1) as

[TABLE]

More specifically, we obtain the following result.

Theorem 2.1**.**

Under assumption (6), the exact, third-order reduced model for system (1)-(2) on the LSM, $\mathcal{W}$, takes the form

[TABLE]

where

[TABLE]

Proof.

This result follows directly upon substitution of the coefficients derived in Lemma A.1 of Appendix A. ∎

When system (1)-(2) is lower dimensional and full modal decomposition can be performed in reasonable time, one can also derive explicit formulas for the LSM coefficients, as we show below in Section 2.3. These formulas enable a direct comparison of our results with other methods, such as the method of normal forms [14], that unavoidably require a full modal transformation of the system.

2.3 Explicit form of the LSM-reduced model when the non-modeling modes are available

After a modal transformation of the form $\mathbf{y}=\bm{\Phi\eta}$, system (1)-(2) can be written as

[TABLE]

where $\omega$ is the natural frequency of the modeling variable $x$, and the diagonal matrix $\bm{\Omega}$ contains the natural frequencies of the non-modeling variables $\mathbf{y}$. The $l^{\text{th}}$-order polynomials $R$ and $\mathbf{S}$ are given by

[TABLE]

Similarly to Section 2.2, we introduce the notation

[TABLE]

We now derive the third-order approximation to $\mathcal{W}$ in modal coordinates in the form

[TABLE]

We restrict (11)-(12) to the LSM by enslaving the $\bm{\eta}$ variable to $(x,\dot{x})$ in eq. (11), which gives

[TABLE]

We then obtain the following result.

Theorem 2.2**.**

Under assumption (6), the exact, third-order reduced model of (11)-(12) on the LSM takes the form

[TABLE]

where

[TABLE]

Proof.

This result is directly obtained upon substitution of the coefficients derived in Lemma B.1 of Appendix B. ∎

Although similar results can be obtained using normal form theory, see e.g. Denis et al. [11], a clear advantage of the LSM reduction method is that it requires no nonlinear change of coordinates. As a consequence, the reduced model (17) keeps the original modeling coordinate $x$ as a physically meaningful modal coordinate. Additionally, the coefficients of the reduced model (17) are simple to compute numerically.

Note that equation (17) takes the form of a Duffing-oscillator if

[TABLE]

This condition is satisfied, for instance, for an even potential function that leads to vanishing second-order nonlinearities. In that case, the reduced model, accurate up to $\mathcal{O}(x^{3})$, is given by

[TABLE]

Consequently, for a conservative system with symmetric cubic nonlinearities, a third-order reduced model can immediately be sought as a Duffing-oscillator. For such systems, therefore, the numerical and experimental procedure of Olivier et al. [11] is justified when they define an equivalent Duffing-oscillator for reduced models obtained from a normal form analysis.

2.4 Reduced and truncated dynamics on the LSM

The full reduced dynamics on the LSM, a symplectic invariant manifold of a canonical Hamiltonian system, is necessarily Hamiltonian. However, system (10) is a truncation of the full reduced dynamics. Therefore, it is a priori unclear if (10) is Hamiltonian or not.

In the following, we show that the truncated reduced model on the LSM is a conservative nonlinear oscillator with closed orbits, and therefore accurately captures the long-term qualitative dynamics on the LSM. As a first step, we rewrite the third-order reduced model (10) as a two-dimensional first-order system

[TABLE]

Note that system (20) is invariant under the transformation $x\rightarrow x$, $y\rightarrow-y$, $t\rightarrow-t$, and therefore is reversible with respect to the line $y=0$ [15].

Theorem 2.3**.**

Let (20) be nondegenerate, i.e., assume that the determinant of its Jacobian at $x=y=0$ is positive, and assume that (20) is reversible. Then, the $x=y=0$ equilibrium of (20) is a center, i.e., system (20) has only periodic orbits in a vicinity of the origin.

Proof.

See Stepanov and Nemytskii [15] who prove this result for general reversible systems. ∎

Theorem 2.3 guarantees that locally around the origin, the truncated system qualitatively captures the long-term dynamics on the LSM. For a complete global understanding of the phase space of system (20), we now prove that the third-order reduced model is a generalized Hamiltonian system.

Theorem 2.4**.**

System (20) is Hamiltonian with the symplectic form $\omega$ and the Hamiltonian $H$ defined as

[TABLE]

Proof.

We derive this result in Appendix C. ∎

The linearization of system (20) reveals that its fixed point, $(x,y)=(0,0)$, is a center in linear approximation. The Hessian of $H(x,y)$, evaluated at $(x,y)=(0,0)$, is

[TABLE]

thus we conclude that the truncated reduced-order model is Lyapunov stable. This proves that, beyond being a close pointwise approximation, the truncated model also accurately represents the long-term, qualitative dynamics on the LSM around the fixed point.

2.5 Comparison with reduced models obtained from the method of modal derivatives

Although the application of modal derivatives (MDs) is based on an idea similar to the LSM reduction, MDs lack the underlying invariant manifold theory that would ensure robustness or validity of the reduced-order model [8]. The manifold proposed by MDs is sought in the quadratic form $\bm{\eta}=\bm{\Theta}x^{2}$, where $\bm{\Theta}$ is expressed, according to [12], as

[TABLE]

where

[TABLE]

This leads to the coefficient vector

[TABLE]

and ultimately to the reduced model

[TABLE]

Theorem 2.5**.**

In the limit of $\Omega_{\_}i\rightarrow\infty$, i.e., for large modal gaps between the modeling and non-modeling frequencies, the MD-reduced model approaches the LSM-reduced model. The two methods, however, do not agree in any other nontrivial case. Near resonances among the natural frequencies, the error of the MD-based reduction becomes unbounded.

Proof.

We derive this result in Appendix D. ∎

By Theorem 2.5, the error of a model reduction using modal derivatives, instead of the LSM-reduced model, will be small for extremely large spectral gaps. Generally, however, the error introduced by the formal MD approach is substantial and grows unbounded near external resonances.

2.6 Comparison with reduced models obtained from the method of normal forms

The method of normal forms, described in detail in [14], offers an alternative way to reduce a nonlinear mechanical system to a single mode and infer its hardening or softening behavior. The decoupling of the modeling mode from the rest is achieved through a nonlinear, near-identity change of variables. In that setting, the system is assumed to be in modal coordinates to begin with, therefore throughout this section, we also assume the input system for reduction in the form of equations (11)-(12).

The key step in the derivation of a normal-form-based reduced model is to introduce a near-identity transformation of the modal coordinate system $\mathbf{q}=[x;\bm{\eta}]$ to a curvilinear coordinate system $\mathbf{s}=[u;\mathbf{v}]$ of the form

[TABLE]

which constitutes a decoupling up to a given order between the modeling and non-modeling modes. The reduced model obtained by this approach, with the same notation as in (4) and (13), reads as

[TABLE]

Although (29) is similar to the LSM reduced model, it can only be computed if the system is available in modal coordinates. In that case, one has to compute only a small subset of the coefficients in the near-identity transformation (27)-(28) to obtain (29). The reduced model correctly identifies hardening or softening behavior up to second order (see Section 4).

The method of normal forms, however, does not guarantee that an actual LSM exists in the phase space under addition of higher-order, unnormalized terms in system (11)-(12). In that sense, a truncated normal-form based reduced model is not justified at the same level of mathematical rigor as the LSM-based reduced model.

3 Dissipative systems

3.1 Setup

We now consider $(n+1)$-degree-of-freedom, periodically forced, damped nonlinear mechanical systems of the form

[TABLE]

where $x$ is a scalar, $\mathbf{y}$ is an $n$-vector; $P$ and $\mathbf{Q}$ are now nonlinear, polynomial functions of $x$, $\dot{x}$, $\mathbf{y}$ and $\dot{\mathbf{y}}$ up to order $2\leq l\leq 3$; $\mathbf{M}$, $\mathbf{C}$ and $\mathbf{K}$ are the mass, damping and stiffness matrices corresponding to the $n$-degree-of-freedom generalized coordinate vector $\mathbf{y}$, while $m$, $c$ and $k$ are the modal mass, damping and stiffness coefficients, corresponding to $x$. We assume Rayleigh-damping, which implies that the modal matrices simultaneously diagonalize the mass, the stiffness and the damping matrices. The mono-harmonic forcing vectors $F_{\_}1\sin(\phi)$ and $\mathbf{F}_{\_}2\sin(\phi)$, do not depend on the positions and velocities. In this setting, our notation of the coefficients of the $l^{\text{th}}$-order polynomials $P$ and $\mathbf{Q}$ is

[TABLE]

Additionally, we introduce the notation

[TABLE]

We denote the solutions of the eigenproblem $\det(m\lambda^{2}+c\lambda+k)=0$ by $\lambda_{\_}{1}^{x}=\bar{\lambda}_{\_}{2}^{x}$, corresponding to the underdamped $x$-mode. We denote the solutions of the eigenproblem $\det(\lambda^{2}\mathbf{M}+\lambda\mathbf{C}+\mathbf{K})=0$ by $\lambda_{\_}{i}^{y}$, $i\in 1,\dots,2n$. We order the eigenvalues by decreasing real parts so that $\text{Re}\lambda_{\_}{2n}^{y}\leq\dots\leq\text{Re}\lambda_{\_}1^{y}<0$. We define $\bm{\Phi}_{\_}2$ as the modal matrix, which consists of the eigenvectors corresponding to each eigenvalue $\lambda_{\_}i^{y}$, for $i\in 1,\dots,2n$. The modeling subspace is fixed to be the slow spectral subspace, i.e., the subspace spanned by the slowest decaying mode, which, in turn, implies that $\text{Re}\lambda_{\_}{2n}^{y}\leq\dots\leq\text{Re}\lambda_{\_}1^{y}<\text{Re}\lambda_{\_}{2}^{x}=\text{Re}\lambda_{\_}{1}^{x}<0$.

3.2 The SSM-reduced model

Haller and Ponsioen [2] give a detailed explanation of different classes of ODEs to which the theory of spectral submanifolds (SSMs) applies. The case relevant for the above setup is a time-periodic slow SSM. To discuss the conditions under which such manifolds exist, we define the absolute spectral quotient of the slow spectral subspace (aligned with the x-degree of freedom) in the form

[TABLE]

with the function $\text{Int}(\cdot)$ referring to the integer part of a real number.

Next, we assume that the low-order non-resonance condition

[TABLE]

holds. Then a two-dimensional, analytic, invariant and attracting time-periodic SSM exists that is unique in the differentiability class $C^{\Sigma+1}$. For a more detailed formulation of absolute and relative spectral quotients of a general spectral subspace, see Haller and Ponsioen [2].

Theorem 3.1**.**

Under assumption (36), the single degree-of-freedom, third-order SSM-reduced model of system (30)-(31) is of the form

[TABLE]

with $\mathbf{z}^{\otimes 2}:=\left[x^{2},x\dot{x},\dot{x}x,\dot{x}^{2}\right]^{\top}$ and with

[TABLE]

where $w_{\_}{11,i}$, $w_{\_}{12,i}$, $w_{\_}{22,i}$ and $\tilde{w}_{\_}{11,i}$, $\tilde{w}_{\_}{12,i}$, $\tilde{w}_{\_}{22,i}$ are defined in (137) and (141), respectively.

Proof.

We derive this result in Appendix F. ∎

In the present dissipative case, we could also use a higher-dimensional modeling variable, $\mathbf{x}\in\mathbb{R}^{\nu}$, $\nu\geq 1$. As an advantage, the corresponding SSM exists under less restrictive non-resonance conditions. By increasing the number of modeling variables, resonances between the modeling and the non-modeling modes can be avoided. We discuss the construction of the SSM-reduced model for $\mathbf{x}\in\mathbb{R}^{\nu}$ in Appendix E.

4 Hardening or softening behavior

By smoothness of all terms involved, finite Taylor expansions for the SSM-reduced dynamics necessarily converge to finite expansions of the LSM-reduced dynamics over the same modeling variable $x$, as the damping matrix $\mathbf{C}$ tends to zero. This follows because both SSM and LSM satisfy the same invariance equations and these equations depend smoothly on the dissipative perturbation to their conservative limit. As a consequence, we can determine the hardening or softening nature of a dissipative system through the periodic orbits of the conservative limit.

4.1 Extracting the backbone curve from the LSM

Under assumption (6), we define the backbone curve of system (1)-(2) as the set of ordered pairs $\left(r,\omega\left(r\right)\right)$ for a periodic orbit family on the LSM, starting from the initial condition

[TABLE]

To obtain an approximation to this curve, we first rewrite (10) as

[TABLE]

where we introduced the short-hand notation

[TABLE]

With the help of this notation, we obtain the following result

Theorem 4.1**.**

The backbone curve of system (1)-(2) can be approximated by

[TABLE]

where,

[TABLE]

Proof.

We derive this result in Appendix G. ∎

We use Theorem 4.1 to generate frequency-energy plots for a conservative system in section 5.2. Additionally, we observe that the second-order LSM approximation (42) takes the form of a diagonal quadratic form (see $\mathbf{w}_{\_}{11}=\mathbf{0}$ in Lemma B.1) for each non-modeling mode, which gives that the eigendirections are aligned with the $x$ and $\dot{x}$ axes. Therefore the periodic orbits on the LSM take their positional extrema in the $x$ direction for each non-modeling mode. To compute this extremum we merely need to take the values that this second order approximation takes in the first eigendirection ($x$ axis). This gives the second-order amplitude estimate for the non-modeling modes

[TABLE]

along periodic orbits on the LSM. Knowing the potential for the forces in a conservative mechanical system together with the estimate (44) gives a fourth-order estimate of the the energy of the periodic orbits of the LSM.

4.2 Comparison of the backbone curves obtained from different model reduction methods

Here we compute the backbone curves of system (11)-(12) obtained under different reduction methods, using the approach explained in the proof of Theorem 4.1. This approach leads to the general approximation formula

[TABLE]

where $\omega_{\_}2$ varies with the reduction method used, as summarized in Table 1.

Table 1 shows that SSM theory and normal forms give the same result for $\omega_{\_}2$, while the prediction for $\omega_{\_}2$ from modal derivatives differs substantially.

5 Numerical simulations

We have implemented the model reduction formulas obtained in the previous sections in Julia [16]. In this section, we illustrate the power of these formulas on several examples. In each case, the full reduced and linearized models are solved numerically using the DifferentialEquations.jl [17] package.

5.1 An analytic example

Haller and Ponsioen [12] gave a minimal nontrivial example of a system on which SSM reduction can be simply demonstrated. The system is of the form

[TABLE]

which is already in modal coordinates, and hence we have $\bm{\Phi}=\mathbf{I}$. The $\dot{x}y$, $x\dot{y}$ and $\dot{x}\dot{y}$ terms in equation (46) are absent. All second-order nonlinearities in $x,\dot{x}$ are zero, therefore $p_{\_}{20\mathbf{0}\mathbf{0}}=p_{\_}{11\mathbf{0}\mathbf{0}}=p_{\_}{02\mathbf{0}\mathbf{0}}=0$. In addition, $p_{\_}{03\mathbf{0}\mathbf{0}}=p_{\_}{12\mathbf{0}\mathbf{0}}=0$ and $p_{\_}{21\mathbf{0}\mathbf{0}}=\mu_{\_}1$ and $p_{\_}{30\mathbf{0}\mathbf{0}}=b$. The linear part in the modeling mode $x$ is characterized by $2\zeta\omega=c_{\_}1$ and $\omega^{2}=k_{\_}1$.

Using equation (37), we obtain the reduced-order model in the form

[TABLE]

in agreement with [12]. We show a simulation of the reduced dynamics on the SSM in Figure 1. Additionally, we show a trajectory of the full system, illustrating its convergence to the reduced dynamics on the SSM.

All of our reduced models are defined locally around a fixed point. Their domain of validity reduces as the parameters approach a configuration in which the non-resonance conditions (36) no longer hold. In the resonant case, the reduced models are not defined anymore. Using the same example, we demonstrate the breakdown of the reduced model in a near-resonant case, where we show that the response of the reduced model is inaccurate even when the full response of the system is in the linear regime (see Figure 2).

5.2 Nonlinear oscillator chain

We continue with a twelve-degree-of-freedom, third-order nonlinear mechanical system of the form

[TABLE]

describing the mechanical model shown in Figure 3. In this system, all the mass, spring and damping coefficients are identical, i.e., $m=k=c=1$. The springs admit a potential function of the form

[TABLE]

As system (49) is not in the form of (30)-(31), a linear transformation is needed first to linearly decouple the modeling mode from the non-modeling ones. Since the mass matrix is diagonal, while the stiffness and damping matrices are tridiagonal Töplitz matrices, we can compute the modal matrix by using the closed form formula (cf. [18])

[TABLE]

For the conservative limit of system (49), we compute the second-order approximated backbone curves and verify our result via numerical continuation using the po toolbox of coco [19] (see Figure 4). To quantify the region of validity of these approximations, we refer to Figure 5, where we observe an increase in energy of two orders of magnitude into the nonlinear regime for both the hardening and the softening systems.

Additionally, we numerically simulate the original damped system, the reduced system based on modal derivatives and the reduced system based on SSM theory. The mismatch between the modal derivatives and the SSM-reduced model is clearly visible in Figure 6.

5.3 Nonlinear Timoshenko beam

In this section, we construct an SSM-reduced model for a discretized nonlinear Timoshenko beam consisting of three elements, resulting in a 32-dimensional phase space (cf. Ponsioen et al. [4]). The chosen beam parameter values are listed in Table 2.

Here $L$ is the length of the beam; $h$ the height of the beam; $b$ the with of the beam; $\rho$ the density; $E$ the Young’s modulus; $G$ the shear modulus; $\eta$ the axial material damping constant and $\mu$ the shear material damping constant.

For the simulation of the full stiff system, we used a fourth-order [20] ESDRIK integrator with a relative and absolute tolerance of $10^{-8}$. The reduced system was integrated with an adaptive $7/6$ Runge-Kutta scheme [21] where we set the tolerances equal to $10^{-12}$. With these settings a speed up factor of 2820 was obtained (from $197$s to $0.07$s) between the integration of the full system and the reduced model, see Figure 7.

In Figure 8, we compare the backbone curve extracted from our LSM reduced model with the numerically continued backbone curve using the po toolbox of coco. In this case the damping parameters have been set to zero to obtain a conservative system, required for LSM reduction and the numerical backbone curve continuation.

6 Conclusions

We have shown that the third-order LSM- and SSM-reduced models accurately capture backbone curves and forced-response curves, respectively, of multi-degree-of-freedom mechanical systems, including a higher-dimensional nonlinear Timoshenko beam, while remaining simple to implement numerically. An advantage of our reduction method is that no near-identity transformations are needed, as opposed to normal-form-based reduction methods. The reduction method keeps the original modeling coordinate as a physically meaningful modal coordinate, resulting in a reduced model that only depends on physical and modal parameters.

Our LSM-reduction method can also be applied when the modal transformation matrix is unfeasible to compute. This is often the case for large systems in which such computations are expensive and often numerically inaccurate. In the conservative case, we have identified conditions under which the formal modal-derivatives-based reduction gives a reasonable approximation of the exact LSM reduction. Additionally, we have shown under which conditions the third-order LSM reduced model is actually a Duffing-oscillator, supporting experimental observations in the literature.

We rigorously justified that in a leading-order approximation, the SSM-reduced system of a periodically forced mechanical system can be seen as the autonomous SSM-reduced system with the modal-participation factor of the first mode added to the reduced system, which justifies the proposed normal form method by Touzé and Amabili [5]. An important further development of the present results will allow for the inclusion of quasi-periodic forcing in system (1)-(2), as the results of Haller and Ponsioen [2] are general enough to allow for such forcing.

Appendix A Derivation of the LSM coefficients in the general case

Lemma A.1**.**

Under assumption (6), the unknown coefficients in (8) satisfy $\mathbf{w}_{\_}{11}=\mathbf{w}_{\_}{21}=\mathbf{w}_{\_}{03}=\mathbf{0}$ and the linear equations

[TABLE]

[TABLE]

with

[TABLE]

Proof.

In the following, we use the invariance of $\mathcal{W}$ to derive the LSM coefficients listed in lemma B.1. Differentiating equation (8) with respect to time, we obtain, for the $i^{\text{th}}$ element of $\dot{\mathbf{y}}$,

[TABLE]

where the gradient of $w_{\_}i$ is

[TABLE]

Similarly, the second time derivative can be written as

[TABLE]

where the Hessian of $w_{\_}i$ is written as

[TABLE]

We express the higher order derivatives, $\ddot{x}$ and $\dddot{x}$, in equation (57) as a function of $x$ and $\dot{x}$, by using equation (1), restricted to the manifold,

[TABLE]

Differentiating equation (59) with respect to time yields,

[TABLE]

Substitution of equations (59) and (60) in (57) gives

[TABLE]

Rewriting equation (2) and substituting $\bm{y}=\mathbf{w}(x,\dot{x})$, we obtain the following expression for the $i^{\text{th}}$ element of $\ddot{\mathbf{y}}$, truncated at $\mathcal{O}(|(x,\dot{x})|^{4})$

[TABLE]

Comparing the coefficients of each monomial term in (61) and (62) leads to a system of algebraic equations for the second order coefficients

[TABLE]

From equation (64), we conclude that $w_{\_}{11,i}=0$, leading to the third order coefficient equations

[TABLE]

where

[TABLE]

Equations (67) and (69) can be solved for $\mathbf{w}_{\_}{21}$ and $\mathbf{w}_{\_}{03}$ leading to $\mathbf{w}_{\_}{21}=\mathbf{w}_{\_}{03}=\mathbf{0}$ Note that the algebraic equations (66)-(69), related to the third order LSM coefficients, depend on the second order LSM coefficients. Therefore one has to solve for the quadratic LSM coefficients first, and substitute the result into the algebraic equations related to the third order. Note that the coupling between different modes of the system only occurs at third-order. Introducing the matrix notation $\mathbf{\Omega}^{2}_{\_}\text{p}:={\mathbf{M}^{-1}\mathbf{K}}$ and $\tilde{\mathbf{Q}}_{\_}{1\mathbf{I}}:=\left[\mathbf{q}_{\_}{1\mathbf{e}_{\_}1},\ldots,\mathbf{q}_{\_}{1\mathbf{e}_{\_}n}\right]$, yields the result stated in Lemma A.1. ∎

Appendix B Derivation of the LSM coefficients when the non-modeling modes are available

Lemma B.1**.**

Under assumption (6), the unknown coefficients in (15) satisfy

[TABLE]

where

[TABLE]

Proof.

This result directly follows from transforming equations (63-69) in modal coordinates. This step makes $\mathbf{M}^{-1}\mathbf{K}$ diagonal and decouples the equations from each other for every mode:

[TABLE]

where

[TABLE]

The solution of this system gives the statement of Lemma B.1. ∎

The LSM theory does not guarantee the existence of a reduced-order model in case of resonance. The LSM coefficients will blow up exactly when the matrices $\bm{\Omega}^{2}\mathbf{D}_{\_}4$ and $\mathbf{D}_{\_}1\mathbf{D}_{\_}9$ become singular. These matrices are diagonal, where their $i^{\text{th}}$ diagonal element can be written as

[TABLE]

As can be seen, 1:1, 1:2 and 1:3 lower-order resonances are explicitly causing a blow-up in the expressions. Higher-order resonances will have the same effect in higher-order approximation. Even though the lower-order approximations are formally computable, there seems to be no awareness of this in existing literature.

Appendix C Proof of Theorem 2.4

For convenience, we restate our two-dimensional first-order system (20)

[TABLE]

Assuming that system (93) is Hamiltonian, then it must be of the general form

[TABLE]

For simplicity we assume that the scalar function $a(x,y)$ only depends on $x$. Equating equations (93) and (94), we obtain

[TABLE]

Integrating equation (95) over $y$ yields

[TABLE]

where $F(x)$ is a scalar function, arising from the integration, that can solely depend on $x$. Differentiating equation (97) with respect to $x$ and equating the result with equation (96), gives

[TABLE]

We observe that $\partial_{\_}xF(x)$ only depends on $x$, therefore we must have that

[TABLE]

having the solutions

[TABLE]

Substituting the solutions for $a(x)$ and $F(x)$ into equation (97), proves Theorem 2.4. ∎

Appendix D Proof of Theorem 2.5

For algebraic equivalence of the MDs and the LSM reduced models, we must have that $\tilde{\mathbf{w}}_{\_}{02}=\mathbf{0}$ because the MDs reduction misses the $\dot{x}^{2}$ term. This condition leads to the equation

[TABLE]

which is satisfied only in the case of $\mathbf{s}_{\_}{2\mathbf{0}}=\mathbf{0}$, which in turn leads to flat manifolds from both the MDs and the LSM reductions.

The main statement of the theorem follows from the following limits for $\mathbf{s}_{\_}{2\mathbf{0}}\neq\mathbf{0}$, written in coordinates:

[TABLE]

Observe that in the resonant limit cases, the error of the MDs reduction is unbounded

[TABLE]

as the MDs reduction method does not take the possible resonances into consideration. This concludes the proof of Theorem 2.5. ∎

Appendix E Construction of the time-periodic SSM reduced model in the general case

E.1 Setup

We now consider $(n+\nu)$-degree-of-freedom, periodically forced, damped nonlinear mechanical systems with a linearly independent partition of the variables so that $\mathbf{x}\in\mathbb{R}^{\nu}$ (modeling variables) and $\mathbf{y}\in\mathbb{R}^{n}$ (non-modeling variables):

[TABLE]

where the mono-harmonic external forcing does not depend on positions and velocities.

We assume Rayleigh-damping, which implies that the modal matrices simultaneously diagonalize the mass, the stiffness and the damping matrices. We denote the solutions of the eigenproblem $\det(\lambda^{2}\mathbf{M}_{\_}1+\lambda\mathbf{C}_{\_}1+\mathbf{K}_{\_}1)=0$ by $\lambda_{\_}{i}^{x}$, $i\in 1,\dots,2\nu$. The eigenvalues are ordered by their real parts so that $\text{Re}\lambda_{\_}{2\nu}^{x}\leq\dots\leq\text{Re}\lambda_{\_}1^{x}<0$. The solutions of the eigenproblem $\det(\lambda^{2}\mathbf{M}_{\_}2+\lambda\mathbf{C}_{\_}2+\mathbf{K}_{\_}2)=0$ are given by $\lambda_{\_}{i}^{y}$, $i\in 1,\dots,2n$. Again, we order the eigenvalues by their real parts so that $\text{Re}\lambda_{\_}{2n}^{y}\leq\dots\leq\text{Re}\lambda_{\_}1^{y}<0$. The modeling subspace is declared to be a slow spectral subspace, i.e., the subspace spanned by the $\nu$ slowest decaying modes. This, in turn, implies that $\text{Re}\lambda_{\_}{2n}^{y}\leq\dots\leq\text{Re}\lambda_{\_}1^{y}<\text{Re}\lambda_{\_}{2\nu}^{x}\leq\dots\leq\text{Re}\lambda_{\_}1^{x}<0$.

E.2 Existence of the reduced model

In the current setting, we seek a time-periodic slow SSM. We assume that the low-order non-resonance conditions

[TABLE]

hold for every $j$ and nonnegative set of integers $\nu_{\_}i$ such that $2\leq\sum_{\_}{k=1}^{\nu}\nu_{\_}k\leq\Sigma$. Then a $2\nu$-dimensional, analytic, invariant and attracting time-periodic SSM exists that is unique in the differentiability class $C^{\Sigma+1}$, and higher, up to analytic.

E.3 Construction of the time-periodic SSM reduced model

We express the damped mechanical system (105)-(106) in modal coordinates by substituting $\mathbf{x}=\bm{\Phi}_{\_}1\bm{\xi}$, $\mathbf{y}=\bm{\Phi}_{\_}2\bm{\eta}$. Subsequently, we project the equations of motion onto the modal directions by left multiplying equations (105) and (106) with the transposed modal matrices $\bm{\Phi}_{\_}1^{\top}$ and $\bm{\Phi}_{\_}2^{\top}$:

[TABLE]

Using the results from Haller and Ponsioen [2], we make use the fact that the SSM perturbs smoothly from the modeling subspace of the linear unperturbed part of system (108)-(109) under the addition of $\mathcal{O}(\varepsilon)$ terms in system (108)-(109). We restrict ourselves to the setting of Breunung and Haller [10] and expand each non-modeling variables $\eta_{\_}i$ in $\varepsilon$. Subsequently, for different orders of $\varepsilon$, every non-modeling variable $\eta_{\_}i$ is expanded in $\bm{\xi}$ and $\dot{\bm{\xi}}$. After a truncating at $\mathcal{O}(|(\bm{\xi},\bm{\dot{\xi}})|^{3},\varepsilon|(\bm{\xi},\bm{\dot{\xi}})|,\varepsilon^{2})$, which is justified if $\bm{\xi}$ and $\dot{\bm{\xi}}$ are of $\mathcal{O}(\varepsilon^{\frac{1}{4}})$, we obtain

[TABLE]

Here, $\mathbf{W}_{\_}{11,i},\mathbf{W}_{\_}{12,i}$ and $\mathbf{W}_{\_}{22,i}$ are $m\times m$ matrices. We specifically aim for a third-order reduced model, which implies that the manifold only has to be approximated up to second order. We uniquely define the quadratic form in equation (110), by requiring that $\mathbf{W}_{\_}{11,i}=\mathbf{W}_{\_}{11,i}^{\top}$ and $\mathbf{W}_{\_}{22,i}=\mathbf{W}_{\_}{22,i}^{\top}$.

Lemma E.1**.**

The coefficient matrices $\mathbf{W}_{\_}{11,i}$, $\mathbf{W}_{\_}{12,i}$ and $\mathbf{W}_{\_}{22,i}$ of the second-order SSM approximation (110) solve the linear matrix equations

[TABLE]

where $\mathbf{S}_{\_}{11,i}$, $\mathbf{S}_{\_}{12,i}$ and $\mathbf{S}_{\_}{22,i}$ are the $m\times m$ matrices filled with the coefficients of the $\bm{\xi},\bm{\dot{\xi}}$ dependent terms from the quadratic part of $S_{\_}{i}(\bm{\xi},\bm{\dot{\xi}},\bm{\eta},\dot{\bm{\eta}})$.

Lemma E.2**.**

The time-dependent coefficients $\bar{W}_{\_}i(\phi)$, of the second order SSM approximation (110), are the periodic solutions of the second-order non-homogeneous ordinary differential equations

[TABLE]

and can be expressed as

[TABLE]

Proof.

We prove Lemmas E.1 and E.2 simultaneously using a direct invariance computation. Differentiating equation (110) with respect to time twice yields

[TABLE]

We express the higher derivatives $\ddot{\bm{\xi}}$ and $\bm{\dddot{\xi}}$ using equation (108) and its derivative with respect to time

[TABLE]

Substituting expressions (118) and (119) into equation (117) gives

[TABLE]

Equation (120) is still quadratic in $\bm{\xi}$ and $\bm{\dot{\xi}}$, and can be rearranged as

[TABLE]

where the $\mathbf{B}_{\_}{kl,i}$ matrices are equal to

[TABLE]

The $i^{\text{th}}$ element of equation (109) can be written as

[TABLE]

which, truncated at third-order in $\bm{\xi}$ and $\bm{\dot{\xi}}$ and taking into consideration that $\eta_{\_}i=\mathcal{O}(|(\bm{\xi},\bm{\dot{\xi}})|^{2})$, leads to

[TABLE]

where $\mathbf{S}_{\_}{11,i}$, $\mathbf{S}_{\_}{12,i}$ and $\mathbf{S}_{\_}{22,i}$ are $\nu\times\nu$ matrices filled with coefficients of the $\bm{\xi},\bm{\dot{\xi}}$ dependent terms from the quadratic part of $S_{\_}{i}(\bm{\xi},\bm{\dot{\xi}},\bm{\eta},\dot{\bm{\eta}})$. Substituting $\eta_{\_}i$ from equation (110) and its time derivative $\dot{\eta}_{\_}i$ into equation (127), leads, after some rearrangement, to the quadratic approximation:

[TABLE]

where $\mathbf{C}_{\_}i$ denotes the hypermatrix

[TABLE]

At leading order, we are now able to equate the quadratic forms in equations (121) and (128), which implies that $\mathbf{B}_{\_}i$ and $\mathbf{C}_{\_}i$ must have equal symmetric parts, i.e.,

[TABLE]

which returns the linear system stated in Lemma E.1. Equating the $\mathcal{O}(\varepsilon|(\bm{\xi},\bm{\dot{\xi}})|^{0})$ terms will return a non-homogeneous second-order ordinary differential equation

[TABLE]

where we introduced the notation $\partial_{\_}\phi(\cdot)=(\cdot)^{\prime}$. Using the method of undetermined coefficients, we obtain the periodic solution

[TABLE]

which proves Lemma E.2. ∎

Appendix F Proof of Theorem 3.1

A general solution of equations (111)-(113), without specifying the dimension $\nu$, cannot be obtained. However, the solution for the practically most relevant case of reducing the system to a single degree of freedom ($\nu=1$) can be written in closed form,

[TABLE]

The determinant of $\mathbf{L}_{\_}i$ can be expressed as

[TABLE]

The construction of the SSM breaks down exactly when $\mathbf{L}_{\_}i$ becomes singular. This implies that the determinant $D$ is zero, which will be the case for a 1:2 resonance ($\zeta_{\_}i=\zeta$, $\omega_{\_}i=2\omega$). Note that a 1:1 resonance is allowed by the SSM theory.

The modal non-modeling coordinate $\eta_{\_}{i}$ (110), for $\nu=1$, can be written as

[TABLE]

Differentiating equation (139) with respect to time and plugging in $\bm{\ddot{\xi}}=\ddot{x}$ from equation (118) gives

[TABLE]

with

[TABLE]

We rewrite equations (139) and (140) using the Kronecker product notation

[TABLE]

where

[TABLE]

The modal non-modeling positions and velocities can be written in vector form as

[TABLE]

The reduced model is obtained by substituting $\mathbf{y}=\bm{\Phi}_{\_}2\bm{\eta}$ and $\mathbf{\dot{y}}=\bm{\Phi}_{\_}2\bm{\dot{\eta}}$ into equation (30)

[TABLE]

Note that the $\mathbf{y}$ and $\mathbf{\dot{y}}$ dependency in the cubic part of $P$ leads to at least fourth-order terms in ($x,\dot{x}$), therefore they are omitted from the third-order reduced model. The polynomial function $P$, can explicitly be written as

[TABLE]

Substituting (152) into equation (151) yields the SSM-reduced model stated in Theorem 3.1.∎

Appendix G Proof of Theorem 4.1

We rescale equation (40) for small vibrations,

[TABLE]

Theorem 2.4 guarantees the existence of periodic orbits in the third order reduced dynamics around the fixed point. We extract the backbone curve from the periodic orbits by continuing these orbits from the linear equation as an unperturbed limit of equation (153).

We denote the, $\epsilon$ dependent, angular frequency with

[TABLE]

where $\omega_{\_}0$ is the linear angular frequency. Introducing a rescaling in time as $\tau=\omega(\epsilon)t$, we obtain

[TABLE]

where the prime denotes differentiation with respect to the rescaled time. We expand $x(\tau)$ in $\epsilon$,

[TABLE]

and substitute the result into equation (155), after which we set $\phi(0)=r,\phi^{\prime}(0)=0$ to define an initial value problem. Collecting terms of $\mathcal{O}(1)$, satisfying the initial conditions, leads to

[TABLE]

Collecting terms of $\mathcal{O}(\epsilon)$, and substituting the solution for $\phi_{\_}0(\tau)$, leads to the following equation

[TABLE]

We set $\omega_{\_}1$ in equation (158) equal to [math] in order to avoid having an external resonance, which would lead to an aperiodic solution. This results in the following solution for $\phi_{\_}1(\tau)$

[TABLE]

Finally, collecting terms of $\mathcal{O}(\epsilon^{2})$, we obtain the resulting equation for $\phi_{\_}2(\tau)$, where we immediately observe that in order to avoid having an external resonance, we must have that

[TABLE]

which completes the second order approximation of the backbone curve and proves Theorem 4.1. ∎