The effect of nonlinear cross-coupling on reduced-order modelling.

The use of reduced-order models (ROMs) for nonlinear systems has received significant attention due to their potential to greatly reduce computational cost, compared to full nonlinear finite-element models. Here, we consider and compare two indirect methods; the applied modal force and enforced modal displacement techniques, paying particular attention to the effect of nonlinear cross-coupling terms. The analysis we present shows that the applied modal force technique is able to account for some effects arising from modes that are not retained in the ROMs, but the resulting accuracy of the ROM depends on the amplitudes selected for the set of forces used to estimate the coefficients of the ROMs. This analysis also shows that the enforced modal displacement technique does not compensate for the effect of modal interactions with modes that are not included in the ROM, but its accuracy is independent of the amplitude of the forces used to estimate the coefficients. The mechanisms that lead to the differences between these techniques is firstly demonstrated using a two conceptually-simple, discrete systems, before a nonlinear beam model is considered.


Introduction
The numerical models used by engineers are becoming increasingly complex. Even with the dramatic increase in computational power, procedures such as optimisation or sensitivity analysis come at great computational expense -an expense that increases significantly with the complexity of the model. Due to this, there has been significant interest in developing reduced-order models (ROMs) that capture, with reasonable accuracy, the dynamics predicted by nonlinear finiteelement (FE) models. Such models are needed to assist in designing lighter, and hence more flexible, structures, such as turbine blades, aircraft and bridges [1,2] and assessing the effect of nonlinearities on hypersonic vehicles [3][4][5].
Reduced-order models have been extensively adopted for linear systems to reduce their dimensions and, more recently, reduction techniques have been extended to account for some nonlinear behaviour. For instance, Craig-Bampton [6] methods can deal with local nonlinearities such as nonlinear springs [7,8], however they cannot be used for deriving reduced-order models of systems with global nonlinearities. In order to deal with global nonlinearities, with particular emphasis on geometrically nonlinear effects [9], direct and indirect methods have been developed. The differences between the direct and indirect methods can be seen in the way in which the nonlinear stiffness forces are described and the modal bases adopted. In direct methods, the nonlinear stiffness forces and their modal counterparts are assumed to be known and the relationship between them is used to determine the nonlinear modal * Corresponding author.
space [10]. Indirect methods assume that the nonlinear terms featuring the finite-element tensors are not known and assumptions on how the resulting nonlinear modal terms are described in the ROM are made [9]. In most commercial software the nonlinear tensors are not known explicitly, instead the nonlinearity is captured via an iterative scheme, thus the indirect methods have attracted considerable attention and are considered here. Two common indirect methods are the Applied Modal Force (AMF) and Enforced Modal Displacement (EMD); both use a set of nonlinear static solutions to estimate a set of nonlinear stiffness coefficients adopted in the ROMs. The main difference between the methods is in the inputs applied to the full FE model and the corresponding outputs used to determine their nonlinear coefficients. The AMF applies a set of static loads to the full FE model and calculates the resulting nonlinear displacement responses. This approach was first introduced by McEwan [2] (referred to as the applied loads procedure in [3] or as Implicit Condensation (IC) in [5]). In contrast, the EMD, introduced by Muravyov and Rizzi [11], and sometimes termed the enforced displacements procedure [3], applies a set of static displacements to the full FE model and uses the resultant forces to determine the stiffness coefficients for the reduced-order model.
With either reduction method, the final common step is the validation of the ROMs. An accurate, but computationally-expensive, solution is the use of nonlinear time-domain simulations for multiple sets of input forces/initial conditions [2,5,12]. An alternative solution is the use of nonlinear normal modes (NNMs) or backbone curves [13,14], which describe the nonlinear undamped-unforced frequency response -an approach that has been adopted in many fields [1,13,[15][16][17][18][19][20]. The NNMs can be analytically determined for some systems by adopting an approximate perturbation method such as the harmonic balance method or normal forms [21]. Alternatively, NNMs can be found numerically, or experimentally, using a Newton-Raphson method; however, computing the Jacobian can be computationally-expensive [7,8]. Considering FE modelling, full-modal continuation methods have been used to compare AMF and EMD techniques applied to full-order continuous systems [1]. Here, we adopt a different approach to understanding the difference between the two techniques, namely, to consider them analytically when applied to a simple system. Using the insights this brings we then consider the results they produce for more complex systems.
Section 2 presents a motivating example, and briefly discusses the two indirect methods -the AMF and the EMD -before considering them analytically. Specifically, the accuracy of the techniques in capturing the nonlinear behaviour of systems is considered, and their ability to account for any effects arising from modal interaction with modes that are not retained in the ROM is assessed. Emphasis is placed on the effect of selecting the input scaling factor on the resulting nonlinear coefficients in the reduced-order model and the resulting curve fit. Although these techniques are typically applied to FE models, the use of analytical models in the current work allows for a direct mathematical comparison between the full-order models and resulting ROMs. In Section 3, further analysis of the motivating example is presented, before a lumped-mass system, exhibiting membrane-like coupling, is considered. In Section 4, a continuous system is used to provide further insight into the applications of this approach, before conclusions are drawn in Section 5.

Reduced-order modelling
A model reduction has the main aim of reducing the dimension of a system in order to allow analysis to be performed in a faster and less computationally-expensive way. In this section, a conceptually-simple, two-degree-of-freedom (2-DoF) model is considered to highlight the issues faced during the reduction process. The observations resulting from this are further discussed in Section 3, where a 10-DoF model is also considered, and then tested on a model of a continuous structure in Section 4.

Motivating example
A 2-DoF spring-mass system, shown schematically in Fig. 1, consists of two masses, both with = 1. These masses are both grounded via linear springs (both with rate = 1) and are connect via another linear spring (also with rate = 1). A nonlinear cubic spring, with coefficient = 0.25, connects the first mass to the ground. For analytical simplicity, the 2-DoF system is treated here as the full-order model. As ROMs are typically described in terms of linear modal coordinates, prior to the reduction the system is transformed using the linear mode shapes, resulting in decoupling of the linear stiffness terms (but typically retaining coupling between linear modes through nonlinear stiffness terms). One advantage of using a linear modal basis for the ROMs is that these are easily computed using commercial finite-element packages [2,5,11,12,22]. Here we consider generating a ROM consisting of the first mode, 1 , with a linear and a cubic stiffness term -reducing the 2-DoF system to a 1-DoF one is artificial, but we believe instructive here due to the relative simplicity of the resulting analysis. Fig. 2 shows the variation of the nonlinear coefficient, which captures cubic stiffness effects via the force = 3 1 , as the input to the full-order model is scaled (this scale factor, , is defined formally in Section 2.2). The cubic coefficient value determined using the EMD (×) is constant and equal to the actual coefficient value in the full model (line) while the value obtained with the AMF (+) is a function of the scale factor. We now consider why this is the case.

The AMF and EMD methods
Let us express the full model in the form where , ,̈, , , ( ) and are the matrices of inertia and stiffness, vector of generalised accelerations and displacements, vector of nonlinear forces and external forcing (both defined in physical coordinateshence the subscript ), respectively.
For reduction purposes, we consider the full model in its modal coordinates . This is achieved using the transform = Φ , where Φ is the mass-normalised matrix of eigenvectors, allowing Eq. (1) to be written ( Here, Λ is the diagonal matrix of the squares of the natural frequencies, is the vector of nonlinear forces projected into the modal coordinates and = Φ is the modal forcing vector. All the linear coefficients in Eq. (2) can be obtained directly from the model; however, the nonlinear coefficients are not always available as, typically, finite-element software iteratively determines the nonlinear effects on the system responses.
Consider a ROM in which modes are retained from the full model (of modes), wherêare the estimates of the modal coordinates retained in the model. Note, if the reduced model perfectly captured the response of the full model then̂= with being a × reduction matrix in which the th row consists of zeros except for a unity element at the location of the th retained mode. Here the forcing is reduced usinĝ = (the effect of forcing in non-retained modes is not captured in the ROM) and the linear properties usinĝ= Λ . As ( ) is not typically available from the full model, cannot be used to reduce it. Instead these nonlinearities are assumed to be a sum of cubic and quadratic terms that combine all the retained modes [2,3,5,15,23], such that the th term in vector̂is given by where the summations in , and span the retained modes. The coefficients , , and , characterise the nonlinear terms and must be estimated in the reduction process.
This identification process differentiates the AMF and EMD techniques. Static responses from the full model are used for both the techniques: the AMF applies a set of static modal forces ( ) as inputs to the full model and uses the resulting displacement responses to estimate the ROM parameters; whereas, for the EMD, the inputs to the full model are static physical deformations ( ) with the resulting forces used in the parameter estimation. Fig. 3 shows flow charts summarising the identification of the coefficients , , and , for the two techniques. The identification is typically achieved using the least squares method with the number of static responses needing to be at least equal to the number of coefficients. However, often the number of coefficients is reduced using physical constraints [3][4][5]23,24].
These physical constraints arise from ensuring conservation of energy, as discussed in [11]. Considering the 2-DoF mass-spring model with the generic cubic stiffness terms discussed earlier (but now, for simplicity, ignoring external forces -= 0), the Lagrangian may be written as Now comparing these equations with the general form for the nonlinearity in the ROM (Eq. (4)), it can be seen that the following constraints exist 1,2,2 = 3 (1) 2,2,2 .
These constrains, along with the equivalent ones for any quadratic terms in the ROM, are needed to ensure conservation of energy in the model. The identification of the coefficients , , and , requires a set of forces/displacements defined using a reference scale factor for the structure of interest ( ). Similarly to the formulation in [5] we define the full model input force as where is a column vector of length . As multiple force (or displacement) inputs are needed for the AMF (or EMD) technique, a family of vectors are used. Each consist of a simple combination of up to 3 modes being non-zero [3,4,23], with the non-zero values based on for the th retained mode. In this paper, only single-mode models are derived, and hence the set of vectors, , take the form

Mathematical insights
Using the 2-mass oscillator, we now consider why the coefficients of the ROM found using the AMF are a function of scale factor, as seen in Fig. 2. In contrast, as observed numerically by Rizzi and Przekop [25], the EMD-based coefficients are not sensitive to scale factor -also shown in the figure. The full model, Eqs. (7) with external forcing added as in Eq. (2), may be written as Retaining the first mode, as before, we can write the ROM equation, Eq. (3), aŝ where we have dropped the (1) 11 term for simplicity. With the AMF method, (1) 1,1,1 can be approximated by applying the static input force = [ 1 , 0] to the full model, Eqs. (12), i.e. we remove the dynamic terms and set 2 = 0. Solving the resulting expression, to find the static modal displacements, provides the data used to identify the ROM coefficients. Here, we analytically investigate how this identification is performed. This analysis is based on a further simplification, namely that as the reduction is centred around the first mode, hence it is reasonable to assume that only a small amount of 2 present, otherwise such a reduction would be ill-advised. Based on this, we neglect terms containing 2 2 and 3 2 , allowing Eqs. (12) to be written as Here we can see that the identified ROM coefficient is a function of the modal displacement, 1 and hence is also a function of the scaling factor. If the EMD is used, the input to the full model is a modal displacement vector = [ 1 , 0] , and the static force response of the full model, Eqs. (12), is Comparing the first of these equations to the static version of the ROM, Eq. (13) witḧ1 = 0, it can be seen that (1) 1,1,1 = ∕(4 2 ), in other words, the identified coefficient in the ROM exactly matches that of the same nonlinear term in the full model. Note that this also means that the identified term is independent of the input amplitude, or scaling factor. Furthermore, this shows that the effect of any coupling with the second, unmodelled, mode is not captured using the EMD method. Fig. 4 shows the cubic coefficient of the full model, AMF method and EMD method, as shown previously in Fig. 2, along with the approximation to the identified nonlinear coefficient using the AMF techniques, i.e. Eq. (15). These predictions, represented by black circles, show a good agreement with numerically-calculated values.
In summary, this analysis shows that the two methods result in different ROMs. The coefficients in the ROM identified using AMF are  a function of the scale factor used in the input forces to the full model. This is due to the fact that all the modal coordinates can be triggered by this force input and that the scale affects the degree of the resulting cross-coupling. This can be useful when strong cross-coupling between modes is present and needs to be captured in an approximate fashion, as is explained in the next section. In contrast, the EMD method predicts the coefficient accurately for this well-defined problem but does not have the ability to capture any cross-coupling between the modes, as 2 is not excited in the static displacement inputs applied to the full model.
These features are examined further in the following sections, where static stiffness and backbone curves for discrete and continuous systems are discussed.

Discrete system analysis
Before considering a continuous system, we extend the discussion of the two-mass system before considering a five-mass system that exhibits coupling between axial-and membrane-like modes.

Two-DoF spring-mass system
The 2-DoF spring-mass system, presented in Section 2.1, is considered further here.
The red line in Fig. 5 shows the relationship between the nonlinear force and the modal displacement for the full model when 2 = 0, and hence is representative of the locus of forces used for the AMF fit. The fit using the AMF-based ROM (red dots) matches this well for both scale factors shown -see inset panel of Fig. 5(a) which shows the fit over the range relevant to the lower scale factor. In contrast, the EMD-based ROM does not capture the full model well when the higher scale factor is considered. Instead, it exactly matches the stiffness relationship for the full model when it is constrained such that 2 = 0 through the appropriate choice of non-zero 2 . This may be considered somewhat artificial, as the second mode is unlikely to be constrained to zero in the full system. For the lower scale factor the match is close, as the 2 content in the full model (with input = [ 1 , 0] ) is small at these amplitudes.
This improvement in force-modal displacement with the AMF method translates, for this simple example, into improved dynamic performance of the AMF-based ROM over that of the EMD-based one. Fig. 6( ) shows the two backbone curves for the full model, projected into the modal space. We see the modal contribution to the first and second nonlinear normal modes (NNMs) -here represents the amplitude of modal coordinate . As the ROM is based on only the first mode, ideally it should match the 1 line for the first NNM. For this system it can be seen that there is a sizeable 2 component to the first NNM and so generating an accurate single-mode ROM is challenging. Fig. 6( ) shows the backbone curves predicted by the EMD-and AMF-based ROMs -note that, as the ROM is a single-mode model, the second modal amplitude is not captured. These backbone curves have been computed from the ROM and the full-order model using the numerical continuation software COCO [26]. It can be seen that, for this example, neither technique produces a ROM capable of fully capturing the nonlinear behaviour of the system using just the first mode, over the region shown here. The main reason is the strong coupling between the two modes, as is apparent in Fig. 6( ). This modal interaction increases with amplitude, which corresponds to a decrease in the accuracy of the ROMs. Considering the techniques in turn, the backbone curves obtained using the EMD do not change with scale factor as the input to the full model does not excite the second mode, as seen in Eq. (16b) and in Fig. 4. Here, the identification process results in a perfect fitting of the cubic coefficient of the first mode. However, as no information regarding the influence of the second mode is generated, no cross-coupling effects are captured, resulting in a poor fit to the full model backbone curves at higher amplitudes (where the cross-coupling influence of the second mode increases). For the AMF method, the backbone curves are a function of the scale factor used. Here, horizontal lines are used to indicate the maximum modal displacement of the full model due to the force inputs for each scale factor. It can be seen that, for the higher scale factors, the ROM curves seem to soften such that they more closely follow the full model at higher amplitudes, at the cost of a poorer fit at lower amplitudes. We suggest that this is due to the fact that the influence of the second mode is partially captured in the AMF technique as seen in the analytical analysis, Eq. (15). A consequence of this is that the coefficient value tends away from that of the cubic term in the full model, ∕(4 2 ). Note that, as the scale factor reduces, the AMF-based ROM tends to that of the EMD-based one, this is because the effect of the cross-coupling reduces as the displacements of the full model are reduced.

A five-mass, 10-DoF system
This section considers the 5-mass, 10-DoF system shown in Fig. 7. These masses are free to move in two directions, and , and are connected via linear springs of stiffness , with an additional two springs grounding the end masses. The length of the springs when at equilibrium is denoted 0 . As such, to achieve a tension of in all springs when the system is stationary, the separation of the grounded points is given by 6 ( 0 + ∕ . Throughout this section, the following parameter values are considered: = 0.1 kg, = 100 kN m −1 , 0 = 10 cm and = 100 N. With these values, the first five linear modes (i.e. the five modes with the lowest natural frequencies) consist of purely vertical ( -direction) motion, whilst the following five linear modes are purely horizontal (in the -direction). These sets of modes are analogous to sets of axial and membrane modes in a continuous structure.
Although the springs in this system are linear, this system exhibits nonlinear behaviour due to the stretching of the springs that results from vertical displacement of the masses. However, in contrast to the system considered in Section 2, this nonlinearity is not perfectly captured by quadratic and cubic nonlinear terms.

Single-mode ROM
In order to investigate the behaviour of the reduced-order modelling techniques when applied to this system, the AMF and EMD methods have been used to compute a ROM consisting of only the first mode of the system. The resulting values of the coefficient of the 3 1 term are shown in Fig. 8. As in the previous case, the EMD method predicts  the same coefficient as the full model 1 at low scale factors. However, in contrast to the case shown in Fig. 4, the EMD method varies with scale factor. This is due to the higher-order terms that are present in the full-order model -i.e. the force-displacement relationship cannot be perfectly captured by a function with only quadratic and cubic terms. The coefficients predicted using the AMF method appear be more robust to scale factor; however, these do not converge to the fullmodel coefficient at low scale factors. To investigate this, the coupling between the modes is now considered.
When applying loads to the first mode, during the AMF method, a strong coupling with the seventh mode (the second axial mode) is observed -i.e. the seventh mode exhibits a relatively large displacement. As such, we now consider the behaviour of the first and seventh modes when a static force is applied to the first mode. For simplicity, it is assumed that the nonlinearity may be approximated up to the cubic order and that all other modes are negligible; hence the equations of motion of the first and seventh modes may be written where 1 and 7 represent the static forces applied to the two modes. Note that some nonlinear parameters (e.g. 2 ) appear in both equations of motion to satisfy the Lagrangian, as discussed in Section 2.2. When the AMF method is used to generate a ROM of the first mode, a force is only applied to the first mode, i.e. 7 = 0. If it is assumed that 7 is small, and hence removing terms containing 2 7 or 3 7 , Eqs. (17) may be written  1 The coefficient of the full model has been found using a Taylor approximation of the full system, up to the third order. (18b) Substituting Eq. (18b) into Eq. (18a), and assuming any terms containing 4 1 or above are negligible, leads to where: (1) .
This shows that, following the assumptions stated above, the quadratic coefficient in the ROM, (1) 1,1 , is equal to the coefficient of 2 1 in the fullorder model, Eq. (17a), and is invariant with scale factor. Note that, for the system and parameters considered here, 1 is zero, and hence it is expected that (1) 1,1 = 0. As seen in the 2-DoF system considered in Section 2, the cubic parameter, (1) 1,1,1 , is a function of 1 and hence will vary with scale factor. However, unlike the previous case, when the scale factor is low (i.e. when 1 is small), the cubic parameter will tend towards which is not equal to the coefficient of 3 1 in the full-order model, Eq. (17a). This difference is highlighted in Fig. 8 where, for low scale factors, the coefficients predicted using the AMF method (blue crosses) tend towards a different value to that of the full-order model and those predicted by the EMD method. The dashed-black line shows the value to which the AMF method is expected to converge, given by Eq. (21). The black circles indicate the predicted value of the AMF method found using Eq. (20), and agree very well with the values found using the AMF method.
As previously depicted in Fig. 5, the red lines in Fig. 9 show the relationship between the nonlinear force and the modal displacement for the full model when a force is applied to only the first modeequivalent to the AMF method. The blue lines show this relationship when the displacements of all but the first mode are constrained to zero -equivalent to the EMD method. Note that, due to the symmetry of the nonlinearity in this system, only the positive forces and displacements are shown. The fit using the AMF-based ROM (red dots) shows a good agreement with the red line at low amplitudes (although it diverges at higher amplitudes). Furthermore, the AMF-based ROM appears to exhibit little change between the scale factor of = 0.01 -in Fig. 9( ) -and = 0.1 -in Fig. 9( ). This reflects the observation in Fig. 8, where the cubic coefficient predicted by the AMF method is relatively invariant to the scale factor. The fit using the EMD-based ROM (blue dots) appears to diverge quickly, with amplitude, for the low scale factor fit (in Fig. 9( )) and shows a different curve, with a poor fit at low amplitudes, for the higher scale factor (in Fig. 9( )). This demonstrates that the constrained full model exhibits more complex force-displacement relationships than can be captured by a quadratic and cubic nonlinearity. This illustrates the reason for the sensitivity of the EMD fit to the scale factor, as seen in Fig. 8. Fig. 10 compares the backbone curves found using the 1-mode ROMs of the AMF and EMD methods, to the first backbone curve of the full-order model, computed using COCO [26]. Three different scale factors are considered for each of the methods: = 0.01, = 0.05 and = 0.1. For all scale factors, the backbone curves of the AMF-based ROMs (dashed lines) show an excellent agreement with the backbone curve of the full-order model (solid blue line), within the region of amplitudes used for fitting (the horizontal lines). Beyond these regions, the AMF backbone curves only diverge slightly. Furthermore, the results of the three different scale factors are almost indistinguishable, due to the robustness of the coefficients to the scale factor. The EMDbased backbone curves (solid lines) show a much greater variation between scale factors (due to the sensitivity of the parameters) and a poor fit to the backbone curve of the full model. Additionally, due to the strong hardening nonlinearity in this structure, the load cases used for the EMD method contain significantly higher displacements than those used for the AMF method (note that the maximum displacement at = 0.1 for the EMD method is beyond the axis limits in Fig. 10). These backbone curves demonstrate that, whilst the EMD method is able to correctly estimate the coefficients of the true model (at low scale factors), it does so by disregarding the coupling between the first mode and the remaining modes of the system (namely the seventh mode, which exhibits a strong coupling). As such, the backbone curves that result from the EMD method can be inaccurate. The AMF method, meanwhile is able to capture this coupling, but this results in a different set of coefficients to the full-order model. Fig. 11. The variation of the cubic coefficient in the first equation of motion, with scale factor, for the 5-mass system. A two-mode model is derived using the EMD and AMF methods.

Two-mode ROM
As it has been observed that the seventh mode couples strongly with the first, a ROM including both of these modes in now derived using the AMF and EMD methods. Note that the dynamics of the first mode (i.e. the first backbone curve) are still of primary interest, and the inclusion of the seventh mode is intended as a means of improving the accuracy of the model.
As described in Eq. (10), the forces (in the case of the AMF method) or displacements (in the case of the EMD method) are given by Note that, due to the extremely high stiffness of the seventh mode, an arbitrary scaling of 0.01 has been applied to 2 . Without this scaling, the forces applied to this mode become extremely large. As described in [5], the loading or displacement combinations are given by  full model is represented by a red line, and the red and blue crosses represent the values predicted by the EMD and AMF methods for the 2mode ROM, respectively. As the coupling with the seventh mode is now captured in this two-mode model, the coefficient predicted by the AMF method now matches that of the full-order model for low scale factors (in contrast with the single-mode model, Fig. 8, where a different value is predicted). As with the single-mode model, the value of the coefficient predicted by the EMD method varies significantly with the scale factor, due to the higher-order terms in the full model. A similar trend is now seen in the AMF method, where the predicted coefficient increases with scale factor, rather than remaining relatively stable as in the single-mode case. This may also be due to the higher-order coupling terms between 1 and 7 , which were considered negligible previously, but which are more significant here as 7 is forced to a higher amplitude.
The backbone curves of the two-mode ROMs, derived using the EMD and AMF methods at different scale factors, are compared to the first backbone curve of the full model in Fig. 12. Comparing these to the backbone curves of the single-mode ROMs, in Fig. 10, clearly shows an improvement in the EMD method. This is due to the coupling between the first-and seventh-mode that is now captured. The backbone curves of the AMF method, however, appear to falsely predict an internal resonance between these two modes (indicated by the sudden change in direction of the backbone curve). Additionally, the AMF backbone curves are significantly less accurate at low amplitudes, before this internal resonance.
These results indicate that the AMF method is more robust when low-frequency modes are included in the ROM, and is able to capture the coupling with higher-frequency mode. However, it does so by predicting coefficient values that are not equal to those in the full model. In contrast, the EMD method requires that higher-frequency modes are included in the ROM to achieve a good fit. However, even without including these modes, it is able to correctly estimate the coefficient values in the full model. It is clear that the accuracy of both methods are affected if the order of the nonlinearity in the ROMs is not sufficiently high.

A continuous system
To further investigate the accuracy of the AMF and EMD methods, a continuous system is now considered. To allow us to assess the actual value of the 2 1 and 3 1 coefficients in the first mode of the full model, rather than using an FE model, we select a Galerkin-based model. The system considered is a pinned-pinned beam with a linear rotational spring at one end, shown in Fig. 13. This system is asymmetric, and is similar to that in [27].
Based on the Euler-Bernoulli beam theory, the equations of motion for the free vibration are given by (24) where defines the distance along the beam, is the transverse displacement of the beam, = 0.5 m is the beam length, and = 100 N m rad −1 denotes the rotational spring constant. Here, is the Dirac delta function, ( ) is the rotation of the beam at the sprungend ( = ), is its second moment of area and is its cross-sectional area. The nonlinearity arises from the dynamic tension in the beam, which is captured by the third term in Eq. (24). Additionally, the beam has a width of 0.03 m and a thickness of 0.001 m, with a density of = 7800 kg m −3 and a Young's modulus of = 200 GPa. We use the separation of variables substitution where ( ) and ( ) are the th mass-normalised modeshape and modal displacement respectively. After applying the boundary conditions, the linear modeshapes are written witĥ and where is found by solving the expression A Galerkin decomposition is now applied by multiplying the equation of motion by and integrating over the length of the beam. Noting the orthogonality between and when ≠ , and between 4 ∕ 4 and when ≠ , leads to the modal equation for the th modë This system is simulated using the first 15 modes, assuming that the contribution of higher modes is negligible. Note that similar expressions are derived in [27,28]. Fig. 14 shows the first backbone curve of the beam, simulated using the first 15 modes. The amplitude of the first 5 modes of this model are shown in terms of the response frequency, . This demonstrates that, in this response region, the first mode dominates the response, and the second modal contribution is significant at higher frequencies, but that the other modes are negligible. This also shows the onset of an internal resonance between the first and second modes, where the second modal amplitude grows, whilst the first mode decreases. It is not expected that a single-mode ROM will be able to capture this resonance. Fig. 15( ) presents the variation of the cubic coefficients adopting the EMD and AMF techniques, varying the values for the scale factor, having retained just one mode. As observed when considering the discrete spring-mass model, the accuracy of the AMF in identifying the cubic coefficient decreases as the scale factor increases, due to the triggering of other modes. As with the two-mass system considered in Section 3, where the nonlinearity of the full-order model is cubic, the EMD method is robust to the scale factor, and predicts the same coefficient as the full model.   Fig. 15( ) compares the backbone curves of the full, 15-mode beam model to those of the ROMs. The backbone curve of the ROM found using the EMD method shows a very close agreement with the full model before the internal resonance. Similarly, the AMF backbone curve, found using a low scale factor of = 10 −4 , shows a good agreement. The horizontal line shows the maximum amplitude of 1 reached when fitting the parameters using the AMF method. As the scale factor is increased, to = 10 −3 and = 10 −2 , the backbone curves of the AMF method begin to deviate, as seen in the previous examples. However, in this example, this deviation is small, given that the scale factor is varying significantly. Note that, due to the large variations in scale factor, the maximum values of 1 that are reached for scale factors of = 10 −3 and = 10 −2 are beyond the limits of this plot. Fig. 16 shows the nonlinear force against the modal displacement for the first mode of the 15-mode beam model. As previously, the red line, showing the full model, represents the response of the system when a static force is applied to just the first mode, whilst the blue line shows the case where the displacement of the first mode is nonzero and other displacements are constrained to zero. These correspond to the loci of forces used for the AMF and EMD methods respectively. The force-displacement relationships predicted by the AMF and EMD methods (red and blue dots respectively) show a decrease in stiffness as the scale factor increases, as reflected in the backbone curves in Fig. 15.

Conclusions
In this paper, two methods for the identification of reduced-order models (ROMs) have been discussed and an investigation of their capability in capturing nonlinear behaviour of systems has been performed; particular emphasis has been placed on the nonlinear coupling between modes. The methods are the enforced modal displacement (EMD) and the applied modal force (AMF). The analysis has been conducted by applying these techniques to two discrete nonlinear analytical models (a two-mass and five-mass model) and a continuous nonlinear analytical model (a beam with pinned ends). As analytical models are used, the true values of the nonlinear coefficients in the full-order models are known, allowing for direct and detailed comparison with the values predicted by the reduced-order modelling methods. The difference between these techniques has been demonstrated by considering the magnitude of the forces (in the case of the AMF method) or displacements (for the EMD) used to calibrate the ROMs.
The coefficients determined by the EMD method are invariant to the magnitude of the enforced displacements when the order of nonlinearity in the full-and reduced-order models match; however, the constraints required to calibrate the EMD could be considered unrepresentative of the true dynamics. The AMF identifies different coefficients as it changes the maximum displacement obtained given a set of forces, even when the correct order of nonlinearity is employed. It has been demonstrated that this variation is due to the capability of the AMF method in capturing the nonlinear coupling with the modes that are not retained in the ROM. These coupling terms may be captured by the EMD method, however this requires that the coupled modes are included in the ROM, resulting in a larger model. At this point, it could be difficult for the user to decide upon which technique is likely to give  the stronger results. While the EMD method can provide a consistent set of parameters, there is potential for the AMF results to capture the dynamic behaviour of the system more accurately. However, it is important to consider that the triggering of other modes using the AMF can badly affect the accuracy of the ROM if the scale of the forcing used for calibration is incorrect. Unfortunately, the optimal choice of this scaling is not typically known a-priori. In such a case, the EMD is a sensible option, particularly if the nonlinearities the system presents are well-modelled. As shown in the present discussion, the EMD is able to provide the correct values for the coefficients linked to the nonlinear terms defined just by the retained modes, assuming the nonlinearities are well-described by the nonlinear terms in the ROM. If the order of the nonlinearity in the ROMs is not sufficiently high, the accuracy of both methods is compromised. This suggests that, although additional data would be required for calibration, ROMs with higherorders of nonlinearity could provide more accurate results. Regardless of which technique is used, the importance of accurately capturing the cross-coupling terms in ROM is a key factor in predicting the correct nonlinear behaviour.