Global System Reduction Order Modeling for Localized Feature Inclusion

Quinn, D. Dane; Brink, Adam R.

doi:10.1115/1.4048890

Abstract

The development of reduced-order models remains an active research area, despite advances in computational resources. The present work develops a novel order-reduction approach that is designed to incorporate isolated regions that contain, for example, nonlinearitites or accumulating damage. The approach is designed to use global modes of the overall system response, which are then naturally coupled to the response in the isolated region of interest. Two examples are provided to demonstrate both the accuracy and the computational efficiency of the proposed approach. The performance of this approach is compared to the exact response corresponding to a finite element simulation for the chosen problems. In addition, the accuracy and computational efficiency are shown relative to a standard Galerkin reduction based on the linear normal modes. It is found that the proposed reduction offer computational efficiency comparable to a Galerkin reduction, but more accurately represents the response of the system when both are compared to the finite element simulation.

1 Introduction

Despite recent advances in computing power, large-scale finite element analysis (FEA) can still present computational difficulties simply arising from, for example, the overall size of the problems considered or large collections of analyses required for uncertainty quantification. Therefore, reduction techniques originating when computational resources were a fraction of current capabilities continue to retain their importance and usefulness. Component mode synthesis (CMS), based on the concept of substructuring, is one of the most widely used reduction methods in structural dynamics. Originally based on the the work of Hurty [1] and extended by Craig and Bampton [2], CMS decomposes the response of an individual component of a larger structure into component and constraint modes [3,4]. The dominant modes of the component can be identified and retained, thus significantly reducing the overall order of the system, while the constraint modes can be used to couple adjacent components together. Such an approach is particularly attractive when the response of the overall structure is localized to only a few of the components of the system, but can be used to develop reduced-order models of the response of complex structural systems.

However, there are numerous examples where the response of the overall system can be characterized as global in nature, so that the response of the system is more closely represented by the structural modes of the overall system rather than localized within only a few components. Jointed structures provide a canonical example of such systems [5–10], where the response of the overall system with joints is closely related to the modes of the overall structure in the absence of the interfaces, known as the monolithic system. In such a case, the structural modes of the monolithic system are significantly different from the vibration modes of the individual components. As identified by Segalman [11–13], in such cases, it is desirable to describe the response and ultimately develop reduced-order models in terms of the structural modes of the overall system, rather than a collection of component-level modes. However, in such a system, the dynamics introduced by the presence of the interfaces, which are internal to the structure, nonetheless significantly influence the response of the overall system [6].

More generally, the response of system with isolated features, including nonlinearitites, could in principle be represented by CMS, so that the isolated effects are introduced through coupling between the constraint modes, while the exterior linear regions were decomposed into their individual component modes. However, the same issue remains—the response of the overall system may be better represented by the modes of the overall system rather than those of the individual (linear) components. In the present work, we develop a novel order-reduction method for systems with isolated nonlinearitites. In contrast to CMS, which utilizes the modes of the individual components to achieve the reduction in the overall model order, the present method employs the modes of the global system, obtained from the corresponding monolithic system in the absence of the nonlinearitites.

2 Global Model Reduction

We consider a general structure as illustrated in Fig. 1, composed of two exterior regions, identified as

C_{1}

and

C_{2}

⁠, which are described by their internal variables

{\underline{x}}_{1} \in R^{N_{1}}

and

{\underline{x}}_{2} \in R^{N_{2}}

⁠, respectively. In addition, there exists an internal isolated region

C_{r}

with variables

{\underline{x}}_{r} \in R^{N_{r}}

⁠. The equations of motion for the exterior regions

C_{1}

and

C_{2}

are assumed to be linear, while the localized interior region includes terms represented by

\underline{N} \in R^{N_{r}}

⁠, which could describe nonlinearitites in the isolated region arising from, for example, the presence of joints and interfaces, or accumulating damage in the interior region. The overall linear mass, damping, and stiffness matrices are represented as

M

⁠,

C

⁠, and

K

⁠, respectively, so that the overall equations of motion are

M \ddot{\underline{x}} + C \dot{\underline{x}} + K \underline{x} + {\underline{N}}_{⋆} = \underline{f} (t)

(1)

with

\underline{x} = [{\underline{x}}_{1} {\underline{x}}_{r} {\underline{x}}_{2}]^{T}, {\underline{N}}_{⋆} = [\underline{0} \underline{N} \underline{0}]^{T}

(2)

The mass, damping, and stiffness matrices are defined over the entire structure and in particular include the linear terms in the equations of motion within the region

C_{r}

⁠. Moreover, in this formulation, it is assumed that the additional terms defined by

\underline{N}

do not depend on the response of the exterior regions as described by

{\underline{x}}_{1}

and

{\underline{x}}_{2}

but instead only depend on the response in the isolated region, so that

\underline{N} = \underline{N} ({\underline{x}}_{r})

⁠. In the investigations below, this system is referred to as the reference system. A mixed formulation for the equations of motion, similar to that presented in Quinn [14], is developed here but extended to a form appropriate for the equations of motion arising from a finite element discretization of the equations of motion for a general structural system.

2.1 Mixed Formulation.

The variables

{\underline{x}}_{r}

within the interior region are further decomposed as

{\underline{x}}_{r} = [{\underline{x}}_{α} {\underline{x}}_{c} {\underline{x}}_{β}]^{T}

(3)

where the components described by

{\underline{x}}_{α}

and

{\underline{x}}_{β}

describe degrees-of-freedom (DOF) on the boundary between components

C_{1}

and

C_{r}

⁠, and

C_{2}

and

C_{r}

⁠, respectively. In particular, the coupling between

{\underline{x}}_{1}

and

{\underline{x}}_{α}

is linear, as is that between

{\underline{x}}_{2}

and

{\underline{x}}_{β}

⁠. However, the degrees-of-freedom

{\underline{x}}_{c}

within

C_{r}

are assumed to be nonlinearly coupled to these boundary coordinates

{\underline{x}}_{α}

and

{\underline{x}}_{β}

⁠. Therefore, the general N degree-of-freedom mechanical system given in Eq. can be expanded in the form

\begin{aligned} [\begin{matrix} M_{11} & M_{1 α} & 0 & 0 & 0 \\ M_{1 α}^{T} & M_{α α} & M_{α c} & 0 & 0 \\ 0 & M_{α c}^{T} & M_{c c} & M_{c β} & 0 \\ 0 & 0 & M_{c β}^{T} & M_{β β} & M_{β 2} \\ 0 & 0 & 0 & M_{β 2}^{T} & M_{22} \end{matrix}] [\begin{matrix} {\ddot{\underline{x}}}_{1} \\ {\ddot{\underline{x}}}_{α} \\ {\ddot{\underline{x}}}_{c} \\ {\ddot{\underline{x}}}_{β} \\ {\ddot{\underline{x}}}_{2} \end{matrix}] \\ + [\begin{matrix} C_{11} & C_{1 α} & 0 & 0 & 0 \\ C_{1 α}^{T} & C_{α α} & C_{α c} & 0 & 0 \\ 0 & C_{α c}^{T} & C_{c c} & C_{c β} & 0 \\ 0 & 0 & C_{c β}^{T} & C_{β β} & C_{β 2} \\ 0 & 0 & 0 & C_{β 2}^{T} & C_{22} \end{matrix}] [\begin{matrix} {\dot{\underline{x}}}_{1} \\ {\dot{\underline{x}}}_{α} \\ {\dot{\underline{x}}}_{c} \\ {\dot{\underline{x}}}_{β} \\ {\dot{\underline{x}}}_{2} \end{matrix}] \\ + [\begin{matrix} K_{11} & K_{1 α} & 0 & 0 & 0 \\ K_{1 α}^{T} & K_{α α} & K_{α c} & 0 & 0 \\ 0 & K_{α c}^{T} & K_{c c} & K_{c β} & 0 \\ 0 & 0 & K_{c β}^{T} & K_{β β} & K_{β 2} \\ 0 & 0 & 0 & K_{β 2}^{T} & K_{22} \end{matrix}] [\begin{matrix} {\underline{x}}_{1} \\ {\underline{x}}_{α} \\ {\underline{x}}_{c} \\ {\underline{x}}_{β} \\ {\underline{x}}_{2} \end{matrix}] \\ + [\begin{matrix} \underline{0} \\ {\underline{N}}_{α} ({\underline{x}}_{α}, {\underline{x}}_{c}) \\ {\underline{N}}_{c} ({\underline{x}}_{α}, {\underline{x}}_{c}, {\underline{x}}_{β}) \\ {\underline{N}}_{β} ({\underline{x}}_{c}, {\underline{x}}_{β}) \\ \underline{0} \end{matrix}] = [\begin{matrix} {\underline{f}}_{1} (t) \\ {\underline{f}}_{α} (t) \\ {\underline{f}}_{c} (t) \\ {\underline{f}}_{β} (t) \\ {\underline{f}}_{2} (t) \end{matrix}] \end{aligned}

(4)

Recall that the nonlinearities described by

{\underline{N}}_{(\cdot)}

are localized to the interior region

C_{r}

⁠. As a result they do not depend on the response of the exterior regions, as described by

{\underline{x}}_{1}

and

{\underline{x}}_{2}

⁠, but instead depend only on the response in the interior region, as indicated in Eq. .

Note that this representation is very similar to that taken in the development of the CMS approach. For example, the equations of motion restricted to the adjoining regions

C_{1}

and

C_{α}

can be expressed as

\begin{aligned} [\begin{matrix} M_{11} & M_{1 α} \\ M_{1 α}^{T} & M_{α α} \end{matrix}] [\begin{matrix} {\ddot{\underline{x}}}_{1} \\ {\ddot{\underline{x}}}_{α} \end{matrix}] + [\begin{matrix} C_{11} & C_{1 α} \\ C_{1 α}^{T} & C_{α α} \end{matrix}] [\begin{matrix} {\dot{\underline{x}}}_{1} \\ {\dot{\underline{x}}}_{α} \end{matrix}] + [\begin{matrix} K_{11} & K_{1 α} \\ K_{1 α}^{T} & K_{α α} \end{matrix}] [\begin{matrix} {\underline{x}}_{1} \\ {\underline{x}}_{α} \end{matrix}] \\ = [\begin{matrix} {\underline{f}}_{1} (t) \\ {\underline{f}}_{α} (t) - (M_{α c} {\ddot{\underline{x}}}_{c} + C_{α c} {\dot{\underline{x}}}_{c} + K_{α c} {\underline{x}}_{c} + {\underline{N}}_{α} ({\underline{x}}_{α}, {\underline{x}}_{c})) \end{matrix}] \end{aligned}

(5)

so that

C_{1}

describes a subcomponent and

C_{α}

represents the boundary interval, while the right-hand side describes the external forcing on the appropriate region.

Over the interior region let

\underline{y}

be the response of an ideal system to be subsequently defined, so that the response of the reference system in the isolated region can be expressed as

{\underline{x}}_{r} = \underline{y} + \underline{z}

(6)

Thus,

\underline{z}

is the difference between the response of the ideal and actual system, referred to as the deviatoric response. Provided that

\underline{N}

is small these are expected to be in some sense close. Using this expansion the equations of motion for

C_{1}

for example can be expressed as

\begin{aligned} M_{11} {\ddot{\underline{x}}}_{1} & + C_{11} {\dot{\underline{x}}}_{1} + K_{11} {\underline{x}}_{1} + (M_{1 α} ({\ddot{\underline{y}}}_{α} + {\ddot{\underline{z}}}_{α}) \\ + C_{1 α} ({\dot{\underline{y}}}_{α} + {\dot{\underline{z}}}_{α}) + K_{1 α} ({\underline{y}}_{α} + {\underline{z}}_{α})) = {\underline{f}}_{1} (t) \end{aligned}

(7)

This equation of motion exactly represents the response of the original structure but requires the response for the deviatoric displacement

{\underline{z}}_{α}

⁠. A similar expression holds in the exterior region

C_{2}

⁠. Nonetheless proceeding, a mixed displacement vector can be defined as

\underline{w} \equiv {[\begin{matrix} {\underline{x}}_{1} & {\underline{y}}_{α} & {\underline{y}}_{c} & {\underline{y}}_{β} & {\underline{x}}_{2} \end{matrix}]}^{T}

(8)

where

({\underline{x}}_{1}, {\underline{x}}_{2})

represent the exact response in the exterior regions and

({\underline{y}}_{α}, {\underline{y}}_{c}, {\underline{y}}_{β})

described the response of the ideal system in the interior region. The resulting equations of motion for

\underline{w}

can be written as

\begin{aligned} [\begin{matrix} M_{11} & M_{1 α} & 0 & 0 & 0 \\ M_{1 α}^{T} & M_{α α} & M_{α c} & 0 & 0 \\ 0 & M_{α c}^{T} & M_{c c} & M_{c β} & 0 \\ 0 & 0 & M_{c β}^{T} & M_{β β} & M_{β 2} \\ 0 & 0 & 0 & M_{β 2}^{T} & M_{22} \end{matrix}] [\begin{matrix} {\ddot{\underline{x}}}_{1} \\ {\ddot{\underline{y}}}_{α} \\ {\ddot{\underline{y}}}_{c} \\ {\ddot{\underline{y}}}_{β} \\ {\ddot{\underline{x}}}_{2} \end{matrix}] \\ + [\begin{matrix} C_{11} & C_{1 α} & 0 & 0 & 0 \\ C_{1 α}^{T} & C_{α α} & C_{α c} & 0 & 0 \\ 0 & C_{α c}^{T} & C_{c c} & C_{c β} & 0 \\ 0 & 0 & C_{c β}^{T} & C_{β β} & C_{β 2} \\ 0 & 0 & 0 & C_{β 2}^{T} & C_{22} \end{matrix}] [\begin{matrix} {\dot{\underline{x}}}_{1} \\ {\dot{\underline{y}}}_{α} \\ {\dot{\underline{y}}}_{c} \\ {\dot{\underline{y}}}_{β} \\ {\dot{\underline{x}}}_{2} \end{matrix}] \\ + [\begin{matrix} K_{11} & K_{1 α} & 0 & 0 & 0 \\ K_{1 α}^{T} & K_{α α} & K_{α c} & 0 & 0 \\ 0 & K_{α c}^{T} & K_{c c} & K_{c β} & 0 \\ 0 & 0 & K_{c β}^{T} & K_{β β} & K_{β 2} \\ 0 & 0 & 0 & K_{β 2}^{T} & K_{22} \end{matrix}] [\begin{matrix} {\underline{x}}_{1} \\ {\underline{y}}_{α} \\ {\underline{y}}_{c} \\ {\underline{y}}_{β} \\ {\underline{x}}_{2} \end{matrix}] \\ = [\begin{matrix} {\underline{f}}_{1} (t) - (M_{1 α} {\ddot{\underline{z}}}_{α} + C_{1 α} {\dot{\underline{z}}}_{α} + K_{1 α} {\underline{z}}_{α}) \\ {\underline{f}}_{α} (t) \\ {\underline{f}}_{c} (t) \\ {\underline{f}}_{β} (t) \\ {\underline{f}}_{2} (t) - (M_{2 β}^{T} {\ddot{\underline{z}}}_{β} + C_{2 β}^{T} {\dot{\underline{z}}}_{β} + K_{2 β}^{T} {\underline{z}}_{β}) \end{matrix}] \end{aligned}

(9)

or

M \ddot{\underline{w}} + C \dot{\underline{w}} + K \underline{w} = \underline{f} (t) + \underline{Q}

(10)

with

\underline{Q}

defined as

\begin{aligned} \underline{Q} \equiv [\begin{matrix} - (M_{1 α} {\ddot{\underline{z}}}_{α} + C_{1 α} {\dot{\underline{z}}}_{α} + K_{1 α} {\underline{z}}_{α}) & \underline{0} & \underline{0} & \underline{0} \end{matrix} \\ {\begin{matrix} - (M_{β 2}^{T} {\ddot{\underline{z}}}_{β} + C_{β 2}^{T} {\dot{\underline{z}}}_{β} + K_{β 2}^{T} {\underline{z}}_{β}) \end{matrix}]}^{T} \end{aligned}

(11)

Similar to the definition of

\underline{z}

⁠,

\underline{Q}

can be interpreted as the deviatoric force acting on the external regions from the internal region. The displacement vector

\underline{w}

mixes the response of the original system in the exterior regions

C_{1}

and

C_{2}

with that of the ideal system in the interior region. Note that the ideal system is not the system in the absence of the nonlinearities as described by

\underline{N} = \underline{0}

⁠. Instead, the coordinates

\underline{y}

describe the response of the system with the additional forcing arising from

\underline{Q}

⁠.

The closure of this system requires the deviatoric response $({\underline{z}}_{α}, {\underline{z}}_{β})$ in the boundary regions be determined outside of this system of equations. In addition, Eq. (10) cannot be used to determine the response within the interior region defined by the boundary and isolated regions $(C_{α}, C_{c}, C_{β})$ ⁠. However, if $({\underline{z}}_{α}, {\underline{z}}_{β})$ can be determined then the resulting response in the exterior regions is exact. In Ref. [14], the deviatoric force $\underline{Q}$ was specified as a constitutive equation, based on a reduced-order model for an interface as developed by Ref. [15].

In terms of

\underline{y}

and

\underline{z}

⁠, the equations of motion in the interior region can be reduced to

\begin{aligned} [\begin{matrix} M_{α α} & M_{α c} & 0 \\ M_{α c}^{T} & M_{c c} & M_{c β} \\ 0 & M_{c β}^{T} & M_{β β} \end{matrix}] [\begin{matrix} {\ddot{\underline{z}}}_{α} \\ {\ddot{\underline{z}}}_{c} \\ {\ddot{\underline{z}}}_{β} \end{matrix}] + [\begin{matrix} C_{α α} & C_{α c} & 0 \\ C_{α c}^{T} & C_{c c} & C_{c β} \\ 0 & C_{c β}^{T} & C_{β β} \end{matrix}] [\begin{matrix} {\dot{\underline{z}}}_{α} \\ {\dot{\underline{z}}}_{c} \\ {\dot{\underline{z}}}_{β} \end{matrix}] \\ + [\begin{matrix} K_{α α} & K_{α c} & 0 \\ K_{α c}^{T} & K_{c c} & K_{c β} \\ 0 & K_{c β}^{T} & K_{β β} \end{matrix}] [\begin{matrix} {\underline{z}}_{α} \\ {\underline{z}}_{c} \\ {\underline{z}}_{β} \end{matrix}] \\ + [\begin{matrix} {\underline{N}}_{α} ({\underline{y}}_{α} + {\underline{z}}_{α}, {\underline{y}}_{c} + {\underline{z}}_{c}) \\ {\underline{N}}_{c} ({\underline{y}}_{α} + {\underline{z}}_{α}, {\underline{y}}_{c} + {\underline{z}}_{c}, {\underline{y}}_{β} + {\underline{z}}_{β}) \\ {\underline{N}}_{β} ({\underline{y}}_{c} + {\underline{z}}_{c}, {\underline{y}}_{β} + {\underline{z}}_{β}) \end{matrix}] = \underline{0} \end{aligned}

(12)

or

M_{r} \ddot{\underline{z}} + C_{r} \dot{\underline{z}} + K_{r} \underline{z} + \underline{N} (\underline{y} + \underline{z}) = \underline{0}

(13)

Thus, the deviatoric response as identified by

\underline{z}

is driven through the nonlinearities

\underline{N}

in the interior region by the response of the ideal displacements

\underline{y}

⁠. The equations of motion in terms of

\underline{w}

and

\underline{z}

are referred to below as the unreduced isolated system.

In summary, the nonlinear equations of motion for the reference system given in Eq. (1) are written in terms of $\underline{x} \in R^{N}$ ⁠, with N = N₁ + N_r + N₂. These have been decomposed into a set of linear equations in terms of the mixed variable $\underline{w} \in R^{N}$ given in Eq. (10) together with a system of nonlinear equations given in Eq. (13), written in terms of $\underline{z} \in R^{N_{r}}$ ⁠. At first glance, this appears to be a step backwards, as the original N degree-of-freedom system has been replaced by a system with N + N_r degrees-of-freedom. However, the linear equations for $\underline{w}$ may be able to be reduced using existing model reduction techniques, so that if N ≫ N_r, a significant computational reduction may be achieved and, as illustrated below, a significant increase in accuracy as compared to traditional Galerkin reduction in terms of the eigenmodes of the linear system.

2.2 Modal Analysis.

With

\underline{N} \neq \underline{0}

⁠, Eq. has no simple representation and in particular no modal decomposition is assumed to be available^². However, such a eigenanalysis can be performed in the linear system with

\underline{N} \equiv \underline{0}

⁠. More generally, various order reduction techniques could be applied to this linear system, including projection based methods [16]. Moreover, assuming proportional damping, the undamped modes

\underline{ϕ}

defined as

ω^{2} M \underline{ϕ} = K \underline{ϕ}

(14)

can be used to decouple the mixed coordinates

\underline{w}

⁠. In particular, these identified modes are global to the system and are not localized to a single exterior region, in contrast to a typical reduction based on component mode synthesis. Using the resulting eigenvectors, the modal decomposition of Eq. reduces to

{\ddot{q}}_{k} + 2 ζ_{k} ω_{k} {\dot{q}}_{k} + ω_{k}^{2} q_{k} = {\underline{ϕ}}_{k}^{T} \underline{f} (t) + {\underline{ϕ}}_{k}^{T} \underline{Q}, k = 1, \dots, N

(15)

with q_k(t) the k^th modal coordinate, defined as

q_{k} = {\underline{ϕ}}_{k}^{T} \underline{w}, {\underline{ϕ}}_{m}^{T} M {\underline{ϕ}}_{n} = δ_{m n}

(16)

where δ_mn is the Kronecker delta function. Other expansions, including those based on Proper Orthogonal Decompositions [17–19], could be employed as well. In this,

{\underline{ϕ}}_{k}^{T} \underline{f} (t)

describes the modal forcing while

{\underline{ϕ}}_{k}^{T} \underline{Q}

describes the effect of the isolated nonlinearitites, both acting on the k^th mode. However, as described above, the deviatoric forcing

\underline{Q}

must be determined outside of Eq. . In particular, the deviatoric displacement

\underline{z}

is driven by the mixed response of the interior region as defined by

\underline{y} = {[\begin{matrix} {\underline{y}}_{α} & {\underline{y}}_{c} & {\underline{y}}_{β} \end{matrix}]}^{T}

(17)

which is expressed in terms of the modal coordinates through the mixed coordinates

{[\begin{matrix} {\underline{x}}_{1} & \underline{y} & {\underline{x}}_{2} \end{matrix}]}^{T} \equiv \underline{w} (t) = \sum_{i = 1}^{N} {\underline{ϕ}}_{i} q_{i} (t)

(18)

If the response of the interior region in the mixed formulation can be well approximated by an M mode truncation then the deviatoric response in the interior region as defined by

\underline{z}

can be determined from Eq. , which can then be used to determine the reduced-order modal response of the system as given in Eq. . With the modal reduction, the resulting system has M + N_α + N_c + N_β degrees-of-freedom. While M ≤ N, typically the retained modes offer a significant reduction in the overall degree-of-freedom for the system, so that M ≪ N₁ + N₂. This formulation nonetheless couples together all of the retained modes through the reconstruction of the mixed response in the interior region from the modal expansion, as determined from Eq. . The resulting response of the system can then be reconstructed from the modal coordinates

\underline{q} (t)

and the interior coordinates

\underline{z} (t)

as

{\underline{x}}^{reduced} = (\sum_{i = 1}^{M} {\underline{ϕ}}_{i} q_{i}) + [\underline{0} \underline{z} \underline{0}]^{T}

(19)

3 Example: Finite Element Analysis Rod Elements

As an illustrative example for this approach, consider the finite discretization of a one-dimensional rod with fixed-free boundary conditions. If overall mass and length of the rod are m and ℓ, respectively, while the axial rigidity is EA, then for the the uniform rod with N elements the corresponding mass and stiffness matrices are

\begin{aligned} M = \frac{ρ ℓ}{N} [\begin{matrix} 1 & 0 & \dots & 0 \\ 0 & 1 & ⋱ \\ ⋮ & ⋱ & ⋱ & ⋱ & ⋮ \\ ⋮ & ⋱ & 1 & 0 \\ 0 & \dots & 0 & 1 \end{matrix}], \\ K = \frac{E A (N + 1)}{ℓ} [\begin{matrix} 2 & - 1 & \dots & 0 \\ - 1 & 2 & ⋱ \\ ⋮ & ⋱ & ⋱ & ⋱ & ⋮ \\ ⋮ & ⋱ & 2 & - 1 \\ 0 & \dots & - 1 & 2 \end{matrix}] \end{aligned}

(20)

where the elements are equally spaced, so that the distance between nodes is Δ = ℓ/N. The spatial position of the rod, denoted s lies in the interval 0 ≤ s ≤ ℓ, and if the interior region is located in the interval s ∈ (s₁, s₂), then the boundary nodes are assumed to correspond to

α = ⌊ \frac{s_{1} N}{ℓ} ⌋

and

β = ⌈ \frac{s_{2} N}{ℓ} ⌉

⁠, so that

\begin{aligned} {\underline{x}}_{1} (t) \equiv {x_{i} (t), i = 1, \dots, α - 1}, {\underline{x}}_{α} (t) \equiv {x_{α} (t)} \\ {\underline{x}}_{c} (t) \equiv {x_{i} (t), i = α + 1, \dots, β - 1} \\ {\underline{x}}_{β} (t) \equiv {x_{β} (t)}, {\underline{x}}_{2} (t) \equiv {x_{i} (t), i = β + 1, \dots, N} \end{aligned}

(21)

The nonlinearities between nodes i and j in the isolated region are assumed to depend on the relative displacement between node pairs, defined as v_i,j(t) ≡ x_i(t) − x_j(t), so that individual forces are assumed to take the form

f (v (t)) = c_{1} v (t) + c_{2} v^{3} (t) + c_{3} \dot{v} (t) + c_{4} {\dot{v}}^{3} (t) + h [v (t)]

(22)

where h[v(t)] is a hysteretic operator acting on the relative displacement v(t). These additional terms used in the interior region are representative of nonlinearities that are often seen in structural systems, including linear detuning and cubic stiffness and damping terms, as well as a hysteretic force commonly arising from an interface [20]. These terms were chosen to highlight the ability of the proposed approach to incorporate such common terms. Therefore, the nonlinear force acting between nodes i and j is defined as

f_{i, j} (t) \equiv f (v_{i, j} (t))

(23)

The resulting nonlinearities in the interior region are therefore identified as

\begin{aligned} {\underline{N}}_{α} \equiv {f_{α + 1, α}} \\ {\underline{N}}_{c} \equiv {- f_{i, i - 1} + f_{i + 1, i}, i = α + 1, \dots, β - 1} \\ {\underline{N}}_{β} \equiv {f_{β, β - 1}} \end{aligned}

(24)

In this description, h[v(t)] is a force component arising from hysteretic damping and is determined by a differential constitutive equation. In the case of a Bouc model for the hysteretic damping, the equations of motion are augmented by constitutive relations of the form

{\dot{h}}_{i, j} = (k_{1} - k_{2} sgn ({\dot{v}}_{i, j}) h_{i, j}) {\dot{v}}_{i, j}

(25)

Note that any representation for the hysteretic damping could be used in the above formulation, including Bouc-Wen or Iwan models [20]. In the present work, Bouc models are chosen for simplicity [21].

For the corresponding undamped linear system, the energy in each mode is defined as

E_{i} = \frac{1}{2} ({\underline{ϕ}}_{i}^{T} M {\underline{ϕ}}_{i}) {\dot{q}}_{i}^{2} + \frac{1}{2} ({\underline{ϕ}}_{i}^{T} K {\underline{ϕ}}_{i}) q_{i}^{2}

(26)

so that the total linearized energy in the system is defined as

E (t) = \sum_{i = 1}^{N} E_{i} (t)

(27)

The error between the reference system and the reduced order approximation can be defined as

E = \frac{1}{N} ‖ {\underline{E}}^{reference} - {\underline{E}}^{reduced} ‖

(28)

where

\underline{E}

is the vector of modal energies, so that

[\underline{E}]_{i} = E_{i}

⁠.

For a system with N degrees-of-freedom, the system is simulated with the nominal parameters

\begin{aligned} ρ = 1, ℓ = 1, E A = 1 \\ c_{1} = 0.1 N, c_{2} = 0.1 N^{3}, c_{3} = 0.1 N, c_{4} = 0.1 N^{3} \\ k_{1} = 0.5 N, k_{2} = 5 N \end{aligned}

(29)

The underlying linear model for the rod could in principle be nondimensionalized so that the above choices for (ρ, ℓ, EA) do not restrict the system dynamics. In the internal region the parameters associated with the nonlinearities are scaled by the number of elements N, so that the equivalent stiffness and damping across the internal interval s ∈ (s₁, s₂) = (0.25, 0.35) remains (approximately) constant as the number of elements varies.

The original equations of motion are simulated with N = 80, and for this discretization, the number of elements in each region is

\begin{aligned} N = 80, ⟶ N_{1} = 19, N_{α} = 1, N_{c} = 7, \\ N_{β} = 1, N_{2} = 52 \end{aligned}

(30)

To illustrate the performance of the above model reduction approach, the model described in Sec. 2.2 is simulated and the response is compared against both the response of the original equations of motion and the response from a standard Galerkin reduction. The Galerkin reduction is obtained by assuming that the response is expanded as

\underline{x} (t) = \sum_{i = 1}^{DOF} {\underline{ϕ}}_{i} p_{i} (t)

(31)

analogous to the expansion given in Eq. , but substituted directly into the original coordinates without isolating the response of the interior region.

For this comparison, the initial conditions of the system are chosen so that the rod is statically displaced by a unit force applied at s = 0.30 and then released from rest. Note that this initial state excites every mode in the structure, and the response of the system at s = 0.50 is shown in Fig. 2. In Fig. 2(a), the response of the system is shown for the original equations of motion, while in Fig. 2(b), the response is shown at the same location for the modally reduced isolated model, retaining M = 10 structural modes. Note that this system has nine degrees-of-freedom associated with the interior region (c.f., Eq. (12)), so that the total degree-of-freedom of the reduced model is 19. Finally, for comparison, the response for a 19 degrees-of-freedom Galerkin reduction is shown in Fig. 2(c).

In Fig. 3(a) the error as defined by Eq. (28) is shown for both the localized reduced-order model as well as the corresponding Galerkin reduction. Although not shown, when the system is simulated with the unreduced isolated formulation given in Eqs. (10) and (13), the response is identical to the simulation with the original equations up to numerical accuracy of the integrator used. The nonlinearities are localized in the interior region s ∈ (s₁, s₂) = (0.25, 0.35). As a result, accurate resolution of the dynamics in this region is critical to capture the effect of the nonlinearities, including the energy loss due to dissipation. In Fig. 3(b), the error in the response of the interior region is shown as a function of time for the modally reduced isolated model as well as for the Galerkin reduction, both retaining 19 total degrees-of-freedom. As illustrated, the isolated model is more accurate in the interior region as compared to the Galerkin reduction.

The performance of this isolated reduction technique is further investigated for both the modally reduced isolated model and the Galerkin reduction as the total degree-of-freedom varies. In Fig. 4(a), the energy error given in Eq. is shown as M, the number of structural modes retained, varies, while in Figs. 4(b), the computational efficiency is shown, defined as the ratio of the computational time of the reduced system relative to that of the original equations of motion, that is

η \equiv \frac{T_{reduced}}{T_{exact}}

(32)

In each simulation, the computational time T is measured using the matlab’s built-in toc function. As observed in Fig. 4(a), there is a minimum degree-of-freedom required for the modally reduced isolated model due to the fixed number of degrees-of-freedom associated with the interior region, which with N = 80 is 9. The error in the reduced isolated system quickly decreases as the number of retained spatial modes increases, while the error in the Galerkin reduced model decreases more slowly. However, because of the additional degrees-of-freedom associated with the interior region in the reduced isolated system, it requires a minimum number of spatial modes to match the accuracy of the Galerkin approximation, in this case with M = 6. Above this spatial resolution, however, the reduced isolated model with retaining M > 6 spatial modes is more accurate than the corresponding Galerkin approximation with M + 9 spatial modes. Moreover, as seen in Fig. 4(b), the computational effort associated with the reduced isolated model is comparable to that of the Galerkin reduction of the same size. Note that the simulation uses the matlab adaptive integrator ode23t, so that the computational time depends on a number of factors in additional to model size, including internal integrator tolerances and estimated numerical error.

4 Example: Finite Element Analysis Beam Elements

Extending on the previous example, the bending deformation of a cantilever beam of unit length is considered. In addition, the equations of motion can be nondimensionalized so that both EI and ρA are considered to be unity and the system is modeled using Euler-Bernoulli beam elements. Each node contains both rotational and transverse degrees-of-freedom, so an P node discretization has N = 2P degrees-of-freedom. In addition, proportional damping is used with ζ = 0.02. The system is assumed to have contain cubic stiffness and damping nonlinearities over the interval s ∈ (0.40, 0.50). External transverse loading is applied to the beam locations at s = s₁ = 0.33 and s = s₂ = 0.67, of the form

f (s_{1}, t) = 15.00 \sin (3.50 t) + 7.50 \sin (5.00 t)

(33a)

f (s_{2}, t) = - 5.00 \sin (3.50 t) - 2.50 \sin (5.00 t)

(33b)

This particular excitation was chosen somewhat arbitrarily, but contains two frequencies that are near the lowest modal frequency of ω = 3.516 of the beam. Finally, the beam is initially at rest in its undeformed position.

With P = 35 nodes for the system, the resulting system has N = 70 degrees-of-freedom. As illustrated in Fig. 5(a), the average error as defined in Eq. (28) depends on the number of modes retained in the expansion of the structural system, and as M increases the error decreases, while as shown in Fig. 5(b), the computational time ratio is reduced for M < N, and approaches unity as M → N.

5 Conclusions

A novel order-reduction approach has been developed specifically for systems that contain isolated features, including nonlinearitites. The approach assumes linear exterior regions connected by localized regions that could, for example, include nonlinearities such as those arising from joints and interfaces. In contrast to a familiar component mode synthesis, the overall system is reduced using the global modes of the system, which then couple to the deviatoric effects introduced by the isolated region. The result is a reduced-order modeling approach that accurately represents the response in the isolated region and therefore can accurately describe the nonlinearities, yet is nonetheless written in terms of the global modes of the system. Two examples are considered, the transient response of a longitudinal rod and the forced response of an Euler-Bernoulli beam. The former is compared against both a standard FEA analysis as well as a Galerkin expansion using the linear mode shapes; the latter is compared against FEA analysis. Note however, that the response within the isolated region as described by $\underline{z}$ has not been reduced. Future efforts can seek to develop appropriate reduction methods focused on this interior region, that can then be combined with the framework presented here.

Footnote

2

Note that $\underline{N}$ need not necessarily be nonlinear, but instead could, for example, include some linear loss of stiffness due to increasing damage. While in this instance a modal decomposition could in principle be undertaken in the damaged state, it is not utilized here. In the present work any additional terms in $\underline{N}$ will be referred to as “nonlinear” regardless of the specific form of the terms.

Acknowledgment

Sandia National Laboratories is a multimission laboratory managed and operated by National Technology and Engineering Solutions of Sandia, LLC, a wholly owned subsidiary of Honeywell International, Inc., for the U.S. Department of Energy’s National Nuclear Security Administration under contract DE–NA0003525. This paper describes objective technical results and analysis. Any subjective views or opinions that might be expressed in the paper do not necessarily represent the views of the U.S. Department of Energy or the United States Government.

Conflict of Interest

There are no conflicts of interest.

Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request. The authors attest that all data for this study are included in the paper. Data are provided by a third party listed in Acknowledgment.

References

1.

Hurty

,

W. C.

,

1965

, “

Dynamic Analysis of Structural Systems Using Component Modes

,”

AIAA J.

,

3

(

4

), pp.

678

–

685

. 10.2514/3.2947

Google Scholar

Crossref

2.

Craig.

,

R. R.

, and

Bampton, Jr

,

M. C. C.

,

1968

, “

Coupling of Substructures for Dynamic Analyses

,”

AIAA J.

,

6

(

7

), pp.

1311

–

1319

. 10.2514/3.4741

Google Scholar

Crossref

3.

Craig, Jr

,

R. R.

, and

Kurdila

,

A. J.

,

2006

,

Fundamentals of Structural Dynamics

. 2nd ed.,

John Wiley & Sons, Inc.

,

Hoboken, NJ

.

4.

Farhat

,

C.

, and

Geradin

,

M.

,

1994

, “

On a Component Mode Synthesis Method and Its Application to Incompatible Substructures

,”

Comput. Struct.

,

51

(

5

), pp.

459

–

473

. 10.1016/0045-7949(94)90053-1

Google Scholar

Crossref

5.

Brake

,

M. R. W.

, editor,

2018

,

The Mechanics of Jointed Structures: Recent Research and Open Challenges for Developing Predictive Models for Structural Dynamics

,

Springer

,

Cham, Switzerland

.

6.

Gaul

,

L.

, and

Lenz

,

J.

,

1997

, “

Nonlinear Dynamics of Structures Assembled by Bolted Joints

,”

Acta Mech.

,

125

(

1

), pp.

169

–

181

. 10.1007/BF01177306

Google Scholar

Crossref

7.

Khattak

,

A. R.

,

Garvey

,

S.

, and

Popov

,

A.

,

2010

, “

Proper Orthogonal Decomposition of the Dynamics in Bolted Joints

,”

J. Sound. Vib.

,

329

, pp.

1480

–

1498

. 10.1016/j.jsv.2009.11.016

Google Scholar

Crossref

8.

Festjens

,

H.

,

Chevallier

,

G.

, and

Dion

,

J. L.

,

2014

, “

Nonlinear Model Order Reduction of Jointed Structures for Dynamic Analysis

,”

J. Sound. Vib.

,

333

, pp.

2100

–

2113

. 10.1016/j.jsv.2013.11.039

Google Scholar

Crossref

9.

Roettgen

,

D. R.

, and

Allen

,

M. S.

,

2017

, “

Nonlinear Characterization of a Bolted, Industrial Structure Using a Modal Framework

,”

Mech. Syst. Signal Proc.

,

84

, pp.

152

–

170

. 10.1016/j.ymssp.2015.11.010

Google Scholar

Crossref

10.

Lacayo

,

R. M.

, and

Allen

,

M. S.

,

2019

, “

Updating Structural Models Containing Nonlinear Iwan Joints Using Quasi-static Modal Analysis

,”

Mech. Syst. Signal Proc.

,

118

, pp.

133

–

157

. 10.1016/j.ymssp.2018.08.034

Google Scholar

Crossref

11.

Segalman

,

D. J.

,

A modal approach to modeling spatially distributed vibration energy dissipation

.

Technical Report SAND2010–4763

,

Sandia National Laboratories

,

August

2010

.

12.

Segalman

,

D. J.

,

2006

, “

Modelling Joint Friction in Structural Dynamics

,”

Struct. Control Health Monitoring

,

13

, pp.

430

–

453

. 10.1002/stc.119

Google Scholar

Crossref

13.

Segalman

,

D. J.

,

2007

, “

Model Reduction of Systems With Localized Nonlinearities

,”

ASME J. Comput. Nonlinear. Dyn.

,

2

(

3

), pp.

249

–

266

. 10.1115/1.2727495

Google Scholar

Crossref

14.

Quinn

,

D. D.

,

2012

, “

Modal Analysis of Jointed Structures

,”

J. Sound. Vib.

,

331

(

1

), pp.

81

–

93

. 10.1016/j.jsv.2011.08.017

Google Scholar

Crossref

15.

Miller

,

J. D.

, and

Quinn

,

D. D.

,

2009

, “

A Two-sided Interface Model for Dissipation in Structural Systems With Frictional Joints

,”

J. Sound. Vib.

,

321

(

1–2

), pp.

201

–

219

. 10.1016/j.jsv.2008.09.037

Google Scholar

Crossref

16.

Benner

,

P.

,

Gugercin

,

S.

, and

Willcox

,

K.

,

2015

, “

A Survey of Projection-based Model Reduction Methods for Parametric Dynamical Systems

,”

SIAM Rev.

,

57

(

4

), pp.

483

–

531

. 10.1137/130932715

Google Scholar

Crossref

17.

Feeny

,

B. F.

, and

Kappagantu

,

R.

,

1998

, “

On the Physical Interpretation of Proper Orthogonal Modes in Vibrations

,”

J. Sound. Vib.

,

211

(

4

), pp.

607

–

616

. 10.1006/jsvi.1997.1386

Google Scholar

Crossref

18.

Chatterjee

,

A.

,

2000

, “

An Introduction to the Proper Orthogonal Decomposition

,”

Curr. Sci.

,

78

(

7

), pp.

808

–

817

.

19.

Kerschen

,

G.

,

Golinval

,

J.-C.

,

Vakakis

,

A. F.

, and

Bergman

,

L. A.

,

2005

, “

The Method of Proper Orthogonal Decomposition for Dynamical Characterization and Order Reduction of Mechanical Systems: An Overview

,”

Nonlinear Dyn. Vol.

,

41

, pp.

147

–

169

. 10.1007/s11071-005-2803-2

Google Scholar

Crossref

20.

Mathis

,

A. T.

,

Balaji

,

N. N.

,

Kuether

,

R. J.

,

Brink

,

A. R.

,

Brake

,

M. R. W.

, and

Quinn

,

D. D.

,

2020

, “

A Review of Damping Models for Structures with Mechanical Joints

,”

ASME Appl. Mech. Rev.

,

72

, p.

040802

. 10.1115/1.4047707

Google Scholar

Crossref

21.

Leine

,

R.

, and

Nijmeijer

,

H.

,

2014

,

Dynamics and Bifurcations of Non-Smooth Mechanical Systems

.

Lecture Notes in Applied and Computational Mechanics

.

Springer

Berlin, Heidelberg

,

ISBN 9783642535796

.

Global System Reduction Order Modeling for Localized Feature Inclusion

Abstract

1 Introduction

2 Global Model Reduction

2.1 Mixed Formulation.

2.2 Modal Analysis.

3 Example: Finite Element Analysis Rod Elements

4 Example: Finite Element Analysis Beam Elements

5 Conclusions

Footnote

Acknowledgment

Conflict of Interest

Data Availability Statement

References

Contents

Data & Figures

Supplements

References

Get Email Alerts

Cited By

ASME Journals

ASME Conference Proceedings

ASME eBooks

Resources

Opportunities

Global System Reduction Order Modeling for Localized Feature Inclusion

Abstract

1 Introduction

2 Global Model Reduction

2.1 Mixed Formulation.

2.2 Modal Analysis.

3 Example: Finite Element Analysis Rod Elements

4 Example: Finite Element Analysis Beam Elements

5 Conclusions

Footnote

Acknowledgment

Conflict of Interest

Data Availability Statement

References

Contents

Data & Figures

Supplements

References

Related

Get Email Alerts

Cited By

Related Articles

Related Proceedings Papers

Related Chapters

ASME Journals

ASME Conference Proceedings

ASME eBooks

Resources

Opportunities

This Feature Is Available To Subscribers Only