Formulation and application of SMU: An open-source MATLAB package for structural model updating

Abstract

This paper describes a MATLAB package for structural model updating, named SMU. The SMU package updates parameter values of a finite element model by solving optimization problems utilizing modal properties obtained from sensor measurements. In particular, the package offers three model updating formulations, namely, (1) the modal assurance criterion value approach, (2) the eigenvector difference approach, and (3) the modal dynamic residual formulation. The first two belong to the family of modal property difference formulations. For each formulation, the analytical Jacobian derivative of the objective function is derived and implemented in SMU. Since the formulated optimization problems are generally nonconvex, the global optimality of the solution cannot be guaranteed using off-the-shelf optimization algorithms. In order to increase the chance of finding a better local minimum, the SMU package can perform gradient search from randomly generated starting points. Several examples for the model updating of as-built structures are included in the GitHub package. This paper demonstrates the SMU functionality through model updating of an 18-DOF model and a concrete building frame model.

Keywords

finite element model updating open-source software nonconvex optimization optimization algorithms modal analysis

Introduction

In modern structural analysis, a great amount of effort has been devoted to finite element (FE) modeling towards simulating the behavior of an as-built structure. In general, predictions by FE models often differ from in-situ experimental measurements. Various inaccuracies in the models can contribute to the discrepancy. For example, idealized connections and support conditions are commonly used in structural analysis and design, while these conditions with infinite or zero stiffness do not exist in reality. In addition, material properties of an as-built structure are usually different from the nominal values, particularly for concrete. To obtain a more accurate FE model that truly reflects behaviors of an as-built structure, data collected from the in-situ experiments can be used to update the values of selected model parameters (e.g., support stiffness or mass parameters) in the FE model. This process is known as FE model updating. An updated FE model can more accurately predict structural behavior under various loading conditions and also serve as the baseline for identifying structural property changes over time, potentially due to deterioration. In the meantime, benefiting from the development of low-cost wireless sensing systems (Dong et al., 2016; Lynch and Loh, 2006; Kane et al., 2014), more and more structural sensors are available for measuring structural responses. As a result, large amounts of sensor data collected from as-built structures are becoming available for FE model updating, which on the other hand, poses computational challenges for performing model updating.

Numerous FE model updating algorithms have been developed and practically applied in the past few decades (Friswell and Mottershead, 1995). Most algorithms can be categorized into two groups, that is, time-domain and frequency-domain approaches. Time-domain approaches deal with time history data collected from the actual structure directly, without the requirement for extracting modal properties (Chatzi and Smyth, 2009; Doucet and Tadić, 2003; Sato and Qi, 1998; Yang et al., 2006). Overall, the time-domain approaches suffer convergence difficulties and high computational cost when applied to large-scale FE models. Different from the time-domain approaches, the frequency-domain approaches can update an FE model using frequency-domain modal properties extracted from experimental measurements. These include resonant frequencies, vibration mode shapes, and damping ratios. In particular, early researchers started by minimizing an objective function consisting of the differences between measured and simulated resonant frequencies. Naturally, this is formulated as a mathematical optimization problem with the selected structural parameters as optimization variables. This category of model updating approaches is named the modal property difference formulation. For example, Zhang et al. (2000) proposed an eigenvalue sensitivity-based model updating approach that was applied on a scaled suspension bridge model, and the updated FE model shows resonant frequencies closer to the experimental measurements. Salawu (1997) reviewed various model updating algorithms using resonant frequency differences and concluded that differences in frequencies may not be sufficient enough for accurately identifying structural parameter values. As a result, other modal properties, for example, mode shapes or modal flexibility, were investigated for model updating. For example, Moller and Friberg (1998) adopted the modal assurance criterion (MAC)-related function for updating the model of an industrial structure, in attempt to make the updated FE model generate mode shapes that are closer to those extracted from experimental measurements. FE model updating using differences in simulated and experimental mode shapes and frequencies was also applied to damage assessment of a reinforced concrete beam (Teughels et al., 2002). Yuen (2012) developed an efficient model updating algorithm using frequencies and mode shapes at only a few selected degrees-of-freedom (DOFs) for the first few modes. Another formulation widely applied in frequency domain model updating is the modal dynamic residual formulation, which minimizes the residuals of generalized eigenvalue problems involving stiffness and mass matrices (Farhat and Hemez, 1993; Kosmatka and Ricles, 1999; Li et al., 2018). The model updating formulations are in general nonconvex, and global optimality of the solution cannot be guaranteed using off-the-shelf local optimization algorithms (Wang et al., 2019a). Attempts have been made to obtain the global optimum of the optimization problems in FE model updating (Bakir et al., 2008; Hofmeister et al., 2019; Li et al., 2018; Otsuki et al., 2021a). Furthermore, model updating of large-scale structures is computationally challenging. Substructure-based approaches have been proposed to reduce the computational cost (Weng et al., 2011; Zhu et al., 2016, 2021).

From a practical point of view, the implementation of structural model updating can be challenging especially for those who are not familiar with mathematical optimization. Difficulties arise in how to construct an objective function and variable constraints, how to select optimization algorithms, and how to deal with nonconvex problems. Optimization algorithms perform gradient search using the Jacobian derivative at every iteration. It has been known that the analytical gradient can achieve better accuracy in model updating compared to the numerical gradient (Dong and Wang, 2018). However, the derivation and implementation of analytical gradients for various model updating formulations can be challenging.

SMU, an open-source MATLAB package for structural model updating, is developed and shared with the research community (Wang et al., 2019b). The package implements three model updating formulations in the frequency domain with the corresponding analytical gradients. Currently, SMU supports optimization algorithms available in the MATLAB optimization toolbox (MathWorks Inc., 2019), such as the Levenberg–Marquardt algorithm, the trust-region-reflective algorithm, and the interior-point method. In order to increase the chance of finding a better local minimum for the nonconvex optimization problem, the package performs gradient search from randomly generated starting points. Several examples for the model updating of as-built structures are included in the GitHub package. SMU has been readily utilized to facilitate research in structural model updating (Dong and Wang, 2019; Otsuki et al., 2021a; Zhang and Sun, 2021).

The rest of the paper is organized as follows. The Structural model updating formulations section presents the three model updating formulations adopted in SMU. The Jacobian derivatives of the formulations section shows the mathematical derivation of the analytical Jacobian derivatives for each model updating formulation. The Introduction of SMU package section describes the key components and user interface of SMU. The Example 1: 18-story steel frame section and the Example 2: concrete building frame section demonstrate the model updating of an 18-story steel frame and a concrete building frame using SMU. Finally, the Conclusions section provides a summary and conclusions.

Structural model updating formulations

The formulations adopted in the current version of SMU update stiffnesses of structural models, while if needed, the functionality can be extended for mass and damping updating. For model updating of a linear structure with N DOFs, the stiffness matrices can be expressed as

K (α) = K_{0} + \sum_{j = 1}^{n_{α}} α_{j} K_{j}

(1)

where

K_{0} \in ℝ^{N \times N}

is the initial stiffness matrix;

n_{α}

is the number of stiffness updating variables;

α_{j}

is the j-th entry of the stiffness updating vector variable

α \in ℝ^{n_{α}}

, which represents the relative change of a stiffness parameter from the initial value; and

K_{j} \in ℝ^{N \times N}

is the influence matrix for the j-th stiffness updating variable

α_{j}

The formulations adopted in SMU are based on the generalized eigenvalue problem in structural dynamics

[K (α) - λ_{i} M] {ψ_{i}} = 0, i = 1 \dots n_{modes}

(2)

where

M \in ℝ^{N \times N}

denotes the mass matrix;

λ_{i} \in ℝ

and

ψ_{i} \in ℝ^{N}

are the i-th eigenvalue and eigenvector, respectively; and

n_{modes}

denotes the number of modes. Note that

λ_{i}

and

ψ_{i}

implicitly depend on

α

, and thus can be denoted as

λ_{i} (α)

and

ψ_{i} (α)

. In general, eigenvectors obtained through field testing are limited to the DOFs measured by sensors. To reflect the measured DOFs, we assume the DOFs are rearranged properly and define

ψ_{i} (α) = [ψ_{i}^{ℳ} (α); ψ_{i}^{U} (α)] \in ℝ^{N}

, where

ψ_{i}^{ℳ} (α) \in ℝ^{n_{ℳ}}

represents the eigenvector entries corresponding to the measured DOFs and

ψ_{i}^{U} \in ℝ^{n_{U}}

corresponds to the unmeasured DOFs. In this paper, we use a semicolon “;” to concatenate vectors/matrices by columns. Note that

n_{ℳ} + n_{U} = N

, the total number of DOFs. In this paper, experimentally obtained eigenvalues and eigenvectors are denoted as

λ_{i}^{EXP}

and

ψ_{i}^{EXP, ℳ}

(at measured DOFs), respectively.

Modal property difference formulation #1: MAC value approach

The SMU package supports the MAC (modal assurance criterion) value approach, which is categorized as a modal property difference formulation

\begin{array}{l} \underset{α}{minimize} \sum_{i = 1}^{n_{modes}} {{(\frac{λ_{i}^{EXP} - λ_{i} (α)}{λ_{i}^{EXP}} \cdot w_{λ_{i}})}^{2} + {(\frac{1 - \sqrt{{MAC}_{i} (α)}}{\sqrt{{MAC}_{i} (α)}} \cdot w_{ψ_{i}})}^{2}} \\ subject to L_{α} \leq α \leq U_{α} \end{array}

(3)

where

w_{λ_{i}}

represents the weighting factor of the i-th eigenvalue difference;

w_{ψ_{i}}

represents the weighting factor of the i-th MAC value difference; and

L_{α}

and

U_{α} \in ℝ^{n_{α}}

denote the lower and upper bounds for the variable

α

, respectively.

{MAC}_{i} (α)

represents the modal assurance criterion between the i-th experimental and simulated eigenvectors/mode shapes at measured DOFs, that is,

ψ_{i}^{EXP, ℳ}

and

ψ_{i}^{ℳ} (α)

{MAC}_{i} (α) = \frac{{({(ψ_{i}^{EXP, ℳ})}^{T} ψ_{i}^{ℳ} (α))}^{2}}{{‖ ψ_{i}^{EXP, ℳ} ‖}_{2}^{2} {‖ ψ_{i}^{ℳ} (α) ‖}_{2}^{2}}, i = 1 \dots n_{modes}

(4)

Here ${‖ \cdot ‖}_{2}$ denotes the $ℒ_{2}$ -norm of a vector. The MAC value represents the similarity between two vectors. When the two vectors are collinear, the MAC value is 1. When two vectors are orthogonal, the MAC value is 0.

At every iteration of an optimization gradient search, the generalized eigenvalue problem in equation (2) is solved using values of $α$ at the current iteration, which produces simulated eigenvalues $λ_{i} (α)$ and eigenvectors $ψ_{i} (α)$ . This process determines that with implicit functions $λ_{i} (α)$ and ${MAC}_{i} (α)$ , the objective function in equation (3) has to remain implicit. For $α$ at the current iteration, the optimization algorithm then evaluates the objective function value in equation (3) and calculates the Jacobian derivative to find the next search gradient. It then moves to the next iteration and repeats this process until convergence criteria are satisfied. Besides having an implicit objective function, the MAC value approach provides a nonconvex optimization problem for which the global optimality of the solution cannot be guaranteed using local optimization algorithms.

Modal property difference formulation #2: eigenvector difference approach

Another modal property difference formulation supported in SMU is the eigenvector difference approach

\begin{array}{l} \underset{α}{minimize} \sum_{i = 1}^{n_{modes}} {{(\frac{λ_{i}^{EXP} - λ_{i} (α)}{λ_{i}^{EXP}} \cdot w_{λ_{i}})}^{2} + {‖ {ψ_{- q_{i}, i}^{EXP, ℳ} - ψ_{- q_{i}, i}^{ℳ} (α)} \cdot w_{ψ_{i}} ‖}_{2}^{2}} \\ subject to L_{α} \leq α \leq U_{α} \end{array}

(5)

where

w_{λ_{i}}

represents the weighting factor of the i-th eigenvalue difference;

w_{ψ_{i}}

represents the weighting factor of the i-th eigenvector difference. The experimental eigenvector at measured DOFs

ψ_{i}^{EXP, ℳ}

is normalized so that the largest entry of

ψ_{i}^{EXP, ℳ}

, denoted as the

q_{i}

-th entry, equals to 1. Accordingly, the simulated eigenvector at measured DOFs

ψ_{i}^{ℳ}

is also normalized so that the

q_{i}

-th entry equals to 1.

ψ_{- q_{i}, i}^{EXP, ℳ}

and

ψ_{- q_{i}, i}^{ℳ} \in ℝ^{n_{ℳ} - 1}

represent the eigenvectors at the measured DOFs with the

q_{i}

-th entry removed. Compared to the MAC value approach, the eigenvector difference approach directly minimizes the differences between the entries in the experimental and simulated eigenvectors. Note that the eigenvector difference approach also has an implicit objective function and is a nonconvex optimization problem.

Modal dynamic residual formulation

The last formulation currently implemented in SMU is the modal dynamic residual formulation. The formulation minimizes the residual of generalized eigenvalue equations in structural dynamics. The residuals are calculated based on the stiffness and mass matrices in combination with experimentally obtained eigenvalues and eigenvectors

\begin{array}{l} \underset{α, ψ^{U}}{minimize} \sum_{i = 1}^{n_{modes}} {‖ [K (α) - λ_{i}^{EXP} M] {\begin{matrix} ψ_{i}^{EXP, ℳ} \\ ψ_{i}^{U} \end{matrix}} \cdot w_{i} ‖}_{2}^{2} \\ subject to L_{α} \leq α \leq U_{α} L_{ψ_{i}^{U}} \leq ψ_{i}^{U} \leq U_{ψ_{i}^{U}}, i = 1 \dots n_{modes} \end{array}

(6)

where

w_{i}

represents the weighting factor for the i-th mode, and

ψ_{i}^{U} \in ℝ^{n_{U}}

represents the unmeasured entries of i-th eigenvector. Compared to the modal property difference formulations, the modal dynamic residual formulation has more updating variables

ψ_{i}^{U}, i = 1 \dots n_{modes}

. Note that the objective function in equation (6) is explicitly expressed by optimization variables

α

and

ψ_{i}^{U}

, while equations (3) and (5) for the modal property difference formulations are not. The modal dynamic residual formulation is convex in the rare case when all DOFs are measured (i.e.,

n_{ℳ} = N

n_{U} = 0

, and

ψ_{i}^{U}

is null). However, when only some DOFs are measured, as often the case in practice, the formulation is nonconvex.

Summary of model updating formulations

Each formulation implemented in SMU has pros and cons in terms of computational efficiency and updating performance. In general, by not involving the MAC value term, the eigenvector difference approach (equation (5)) is computationally more efficient than the MAC value approach (equation (3)). Numerical studies have shown that the eigenvector difference approach usually achieves higher accuracy than the MAC value approach (Dong and Wang, 2018). The modal dynamic residual formulation (equation (6)) is based on a different concept from the MAC value and the eigenvector difference approaches. The modal dynamic residual formulation does not directly minimize the discrepancy of modal properties, and thus, the updated model may have modal properties that are more different from experimental values. It is highly recommended to investigate and compare model updating results using different formulations especially when dealing with experimental data.

In general, when using experimental data, the system identification results of eigenvectors have more variability and can be less accurate compared to eigenvalues. Moreover, it is more challenging to accurately identify the eigenvalues and eigenvectors for higher modes than those for lower modes. Therefore, users can take into account the reliability of the identified modal parameters through the use of weighting factors, $w_{λ_{i}}$ and $w_{ψ_{i}}$ in equations (3) and (5) or $w_{i}$ in equation (6). Furthermore, if several vibration records and system identification results are available, statistical properties of identified results, such as standard deviation, can be reflected in weighting factors to consider uncertainties in measurements and system identification results.

Jacobian derivatives of the formulations

The numerical optimization algorithms solving equations (3), (5), and (6) are iterative. At every iteration step, the Jacobian derivative (short-named as the Jacobian) of the objective function is often used to determine the search direction. To facilitate the Jacobian derivation, the optimization problems in equations (3), (5), and (6) can be rewritten as

\underset{x}{minimize} f (x) = r {(x)}^{T} r (x) = {‖ r (x) ‖}_{2}^{2} subject to L_{x} \leq x \leq U_{x}

(7)

where

x \in ℝ^{n_{x}}

is a vector variable;

r (x) : ℝ^{n_{x}} \to ℝ^{n_{r}}

is a residual vector function;

f (x) : ℝ^{n_{x}} \to ℝ

is an objective function;

n_{x}

is the length of the vector variable

x

; and

n_{r}

is the length of the residual vector function

r

. The Jacobian of the objective function

f (x)

is defined as

D_{x} f = [\begin{matrix} \partial f / \partial x_{1} & \partial f / \partial x_{2} & \dots & \partial f / \partial x_{n_{x}} \end{matrix}] \in ℝ^{1 \times n_{x}}

. In the following subsections, the Jacobians of the objective functions shown in equations (3), (5), and (6) will be derived.

Jacobian of modal property difference formulation #1: MAC value approach

When using the MAC value approach in equation (3), the vector variable x in equation (7) is the stiffness updating vector variable $α$ . Hence, the MAC value formulation in equation (3) is rewritten in least squares form as $f (α) = r {(α)}^{T} r (α)$ with a residual vector function $r (α) : ℝ^{n_{α}} \to ℝ^{2 \cdot n_{m o d e s}}$ as

r (α) = [\begin{matrix} r_{1} (α) \\ ⋮ \\ r_{n_{modes}} (α) \end{matrix}] = [\begin{matrix} (λ_{1}^{EXP} - λ_{1} (α)) / λ_{1}^{EXP} \cdot w_{λ_{1}} \\ (1 - \sqrt{{MAC}_{1} (α)}) / \sqrt{{MAC}_{1} (α)} \cdot w_{ψ_{1}} \\ ⋮ \\ (λ_{n_{modes}}^{EXP} - λ_{n_{modes}} (α)) / λ_{n_{modes}}^{EXP} \cdot w_{λ_{n_{modes}}} \\ (1 - \sqrt{{MAC}_{n_{modes}} (α)}) / \sqrt{{MAC}_{n_{modes}} (α)} \cdot w_{ψ_{n_{modes}}} \end{matrix}]

(8)

where

r_{i} (α) : ℝ^{n_{α}} \to ℝ^{2}

equals

[\begin{matrix} (λ_{i}^{EXP} - λ_{i} (α)) / λ_{i}^{EXP} \cdot w_{λ_{i}}; & (1 - \sqrt{{MAC}_{i} (α)}) / \sqrt{{MAC}_{i} (α)} \cdot w_{ψ_{i}} \end{matrix}], i = 1 \dots n_{modes}

. The Jacobian for

f (α)

D_{α} f \in ℝ^{1 \times n_{α}}

, equals

D_{r} f \cdot D_{α} r

by using the chain rule. The first term is

D_{r} f = 2 r^{T} \in ℝ^{1 \times (2 \cdot n_{modes})}

. The second term is

D_{α} r = [\begin{matrix} D_{α} r_{1}; & D_{α} r_{2}; & \dots & D_{α} r_{n_{modes}} \end{matrix}] \in ℝ^{(2 \cdot n_{modes}) \times n_{α}}

. Recall the definition of the MAC value in equation (4), each

D_{α} r_{i} \in ℝ^{2 \times n_{α}}

can be formed as

D_{α} r_{i} = [\begin{matrix} - \frac{w_{λ_{i}}}{λ_{i}^{EXP}} \cdot D_{α} (λ_{i} (α)) \\ (\frac{- w_{ψ_{i}}}{\sqrt{{MAC}_{i} (α)}}) (\frac{{(ψ_{i}^{EXP, ℳ})}^{T}}{{(ψ_{i}^{EXP, ℳ})}^{T} ψ_{i}^{ℳ} (α)} - \frac{{(ψ_{i}^{ℳ} (α))}^{T}}{{‖ ψ_{i}^{ℳ} (α) ‖}_{2}^{2}}) D_{α} (ψ_{i}^{ℳ} (α)) \end{matrix}], i = 1 \dots n_{modes}

(9)

The formulation for $D_{α} (λ_{i} (α)) \in ℝ^{1 \times n_{α}}$ and $D_{α} (ψ_{i}^{ℳ} (α)) \in ℝ^{n_{ℳ} \times n_{α}}$ have been well studied by researchers (Fox and Kapoor, 1968; Nelson, 1976). Nevertheless, a simplified way of obtaining $D_{α} (ψ_{i}^{ℳ} (α))$ based on the normalization of $ψ_{i}^{ℳ} (α)$ is presented, without expressing the derivative as a linear combination of all the eigenvectors (as in Nelson (1976)). Recall the generalized eigenvalue equation for the i-th mode

[K (α) - λ_{i} M] {ψ_{i}} = 0

(10)

By differentiating equation (10) with respect to the j-th updating variable, $α_{j}$ , the following equation can be obtained

[K (α) - λ_{i} M] \frac{\partial ψ_{i}}{\partial α_{j}} = \frac{\partial λ_{i}}{\partial α_{j}} M ψ_{i} - K_{j} ψ_{i}

(11)

Assume that the eigenvalues are distinct and define the modal mass of the i-th mode as $m_{i} = {(ψ_{i})}^{T} M ψ_{i}$ . Then, pre-multiply equation (11) by ${(ψ_{i})}^{T}$ , and note ${(ψ_{i})}^{T} [K (α) - λ_{i} M] = 0$ . Equation (11) can be simplified as follows

{(ψ_{i})}^{T} K_{j} ψ_{i} = \frac{\partial λ_{i}}{\partial α_{j}} {(ψ_{i})}^{T} M ψ_{i}

(12)

and

\partial λ_{i} / \partial α_{j} \in ℝ

can be obtained

\frac{\partial λ_{i}}{\partial α_{j}} = \frac{{(ψ_{i})}^{T} K_{j} ψ_{i}}{m_{i}}

(13)

As a result, the Jacobian of the i-th simulated eigenvalue, $D_{α} (λ_{i}) \in ℝ^{1 \times n_{α}}$ , with respect to updating vector, α, can be found as follows

D_{α} (λ_{i}) = [\begin{matrix} \frac{\partial λ_{i}}{\partial α_{1}} & \frac{\partial λ_{i}}{\partial α_{2}} & \dots & \frac{\partial λ_{i}}{\partial α_{n_{α}}} \end{matrix}] = [\begin{matrix} \frac{{(ψ_{i})}^{T} K_{1} ψ_{i}}{m_{i}} & \frac{{(ψ_{i})}^{T} K_{2} ψ_{i}}{m_{i}} & \dots & \frac{{(ψ_{i})}^{T} K_{n_{α}} ψ_{i}}{m_{i}} \end{matrix}]

(14)

After obtaining $\partial λ_{i} / \partial α_{j}$ , equation (11) is reused to find the only remaining unknown term, $\partial ψ_{i} / \partial α_{j} \in ℝ^{N}$ . However, $\partial ψ_{i} / \partial α_{j}$ cannot be directly obtained from equation (11) because $[K (α) - λ_{i} M]$ is rank deficient by one assuming that the eigenvalue $λ_{i}$ is distinct. Nevertheless, as previously mentioned, $ψ_{i}^{ℳ}$ is normalized so that the $q_{i}$ -th entry always equals a constant 1. As a result, the $q_{i}$ -th entry in vector $\partial ψ_{i}^{ℳ} / \partial α_{j}$ is zero, that is, $\partial ψ_{q_{i}, i}^{ℳ} / \partial α_{j} = 0$ . Because of the separation by measured and unmeasured DOFs, $ψ_{i} (α) = [\begin{matrix} ψ_{i}^{ℳ} (α); & ψ_{i}^{U} (α) \end{matrix}]$ , the $q_{i}$ -th entry in $\partial ψ_{i} / \partial α_{j}$ is also zero, that is, $\partial ψ_{q_{i}, i} / \partial α_{j} = 0$ . This is utilized to resolve the rank deficiency issue of $[K (α) - λ_{i} M] .$ Specifically, define $P_{i} = [\begin{matrix} Q_{i} & 0_{(n_{ℳ} - 1) \times n_{U}} \\ 0_{n_{U} \times n_{ℳ}} & I_{n_{U}} \end{matrix}] \in ℝ^{(N - 1) \times N}$ where $Q_{i} \in ℝ^{(n_{ℳ} - 1) \times n_{ℳ}}$ selects the entries of measured DOFs except for the $q_{i}$ -th entry as:

Q_{i} = [\begin{matrix} I_{q_{i} - 1} & 0_{(q_{i} - 1) \times 1} & 0_{(q_{i} - 1) \times (n_{ℳ} - q_{i})} \\ 0_{(n_{ℳ} - q_{i}) \times (q_{i} - 1)} & 0_{(n_{ℳ} - q_{i}) \times 1} & I_{n_{ℳ} - q_{i}} \end{matrix}]

(15)

where

I_{q_{i} - 1}

and

I_{n_{ℳ} - q_{i}}

denote identity matrices with size of

q_{i} - 1

and

n_{ℳ} - q_{i}

, respectively. Then, pre-multiplying and post-multiplying

[K (α) - λ_{i} M]

in equation (11) by

P_{i}

and

P_{i}^{T}

to cross out the

q_{i}

-th row and

q_{i}

-th column,

B_{i} \in ℝ^{(N - 1) \times (N - 1)}

is generated

B_{i} = P_{i} [K (α) - λ_{i} M] P_{i}^{T}

(16)

Next, pre-multiply $(\frac{\partial λ_{i}}{\partial α_{j}} M ψ_{i} - K_{j} ψ_{i})$ in equation (11) by $P_{i}$ to eliminate the $q_{i}$ -th row and obtain $b_{i j} \in ℝ^{N - 1}$

b_{i j} = P_{i} \cdot (\frac{\partial λ_{i}}{\partial α_{j}} M ψ_{i} - K_{j} ψ_{i})

(17)

Finally, recalling $\partial ψ_{q_{i}, i} / \partial α_{j} = 0$ , the elimination of the $q_{i}$ -th row in equation (11) is equivalent to the following

B_{i} {\begin{matrix} \frac{\partial (Q_{i} ψ_{i}^{ℳ})}{\partial α_{j}} \\ \frac{\partial (ψ_{i}^{U})}{\partial α_{j}} \end{matrix}} = b_{i j}

(18)

Thus, the Jacobian of the i-th simulated eigenvector with respect to the updating variables can be shown as

{\begin{matrix} \frac{\partial (Q_{i} ψ_{i}^{ℳ})}{\partial α_{j}} \\ \frac{\partial (ψ_{i}^{U})}{\partial α_{j}} \end{matrix}} = B_{i}^{- 1} b_{i j}

(19)

In summary, $\partial ψ_{i}^{ℳ} / \partial α_{j}$ is 0 at the $q_{i}$ -th entry, and other entries are provided by the equation above. The Jacobian of the simulated eigenvector at measured DOFs, $D_{α} (ψ_{i}^{ℳ}) \in ℝ^{n_{ℳ} \times n_{α}}$ in equation (9), can be obtained as follows

\begin{array}{l} D_{α} (ψ_{i}^{ℳ}) = [\begin{matrix} \frac{\partial ψ_{i}^{ℳ}}{\partial α_{1}} & \frac{\partial ψ_{i}^{ℳ}}{\partial α_{2}} & \dots & \frac{\partial ψ_{i}^{ℳ}}{\partial α_{n_{α}}} \end{matrix}] \\ = [\begin{matrix} \frac{\partial ψ_{1, i}^{ℳ}}{\partial α_{1}} & \frac{\partial ψ_{1, i}^{ℳ}}{\partial α_{2}} & \dots & \frac{\partial ψ_{1, i}^{ℳ}}{\partial α_{n_{α}}} \\ ⋮ & ⋮ & \dots & ⋮ \\ \frac{\partial ψ_{q_{i} - 1, i}^{ℳ}}{\partial α_{1}} & \frac{\partial ψ_{q_{i} - 1, i}^{ℳ}}{\partial α_{2}} & \dots & \frac{\partial ψ_{q_{i} - 1, i}^{ℳ}}{\partial α_{n_{α}}} \\ 0 & 0 & \dots & 0 \\ \frac{ψ_{q_{i} + 1, i}^{ℳ}}{\partial α_{1}} & \frac{\partial ψ_{q_{i} + 1, i}^{ℳ}}{\partial α_{2}} & \dots & \frac{\partial ψ_{q_{i} + 1, i}^{ℳ}}{\partial α_{n_{α}}} \\ ⋮ & ⋮ & \dots & ⋮ \\ \frac{ψ_{n_{ℳ,} i}^{ℳ}}{\partial α_{1}} & \frac{ψ_{n_{ℳ,} i}^{ℳ}}{\partial α_{2}} & \dots & \frac{ψ_{n_{ℳ,} i}^{ℳ}}{\partial α_{n_{α}}} \end{matrix}] | \begin{matrix} 1 \\ ⋮ \\ q_{i} - 1 \\ q_{i} \\ q_{i} + 1 \\ ⋮ \\ n_{ℳ} \end{matrix} \end{array}

(20)

After obtaining the Jacobian of the simulated eigenvalue and eigenvector at measured DOFs, $D_{α} (λ_{i})$ and $D_{α} (ψ_{i}^{ℳ})$ , the analytical Jacobian in equation (9) can be calculated.

Jacobian of modal property difference formulation #2: eigenvector difference approach

For the eigenvector difference approach in equation (5), the optimization variable x in equation (7) is the stiffness updating vector variable α. In order to derive the Jacobian of the eigenvector difference formulation in equation (5), a residual vector function $r (α) : ℝ^{n_{α}} \to ℝ^{n_{ℳ} \cdot n_{m o d e s}}$ is defined as follows

r (α) = [\begin{matrix} r_{1} (α) \\ ⋮ \\ r_{n_{m o d e s}} (α) \end{matrix}] = [\begin{matrix} \frac{(λ_{1}^{EXP} - λ_{1} (α))}{λ_{1}^{EXP}} \cdot w_{λ_{1}} \\ {ψ_{- q_{1}, 1}^{EXP, ℳ} - ψ_{- q_{1}, 1}^{ℳ} (α)} \cdot w_{ψ_{1}} \\ ⋮ \\ \frac{(λ_{n_{modes}}^{EXP} - λ_{n_{modes}} (α))}{λ_{n_{modes}}^{EXP}} \cdot w_{λ_{n_{modes}}} \\ {ψ_{- q_{n_{modes}}, n_{modes}}^{EXP, ℳ} - ψ_{- q_{n_{modes}}, n_{modes}}^{ℳ} (α)} \cdot w_{ψ_{n_{modes}}} \end{matrix}]

(21)

where

r_{i} (α) : ℝ^{n_{α}} \to ℝ^{n_{ℳ}}

equals

[\begin{matrix} (λ_{i}^{EXP} - λ_{i} (α)) / λ_{i}^{EXP} \cdot w_{λ_{i}}; & {ψ_{- q_{i}, i}^{EXP, ℳ} - ψ_{- q_{i}, i}^{ℳ} (α)} \cdot w_{ψ_{i}} \end{matrix}], i = 1 \dots n_{modes}

. Using

r (α)

, the optimization problem in equation (5) for the eigenvector difference formulation can also be rewritten the same as equation (7), with the objective function

f (α) = r {(α)}^{T} r (α)

. Again, the Jacobian for

f (α)

D_{α} f \in ℝ^{1 \times n_{α}}

equals

D_{r} f \cdot D_{α} r

from the chain rule. However, the residual vector

r

has a different dimension for the MAC value formulation. For the eigenvector difference formulation, the first Jacobian term is

D_{r} f = 2 r^{T} \in ℝ^{1 \times (n_{ℳ} \cdot n_{m o d e s})}

. Meanwhile, the second term is

D_{α} r = [\begin{matrix} D_{α} r_{1}; & D_{α} r_{2}; & \dots & D_{α} r_{n_{modes}} \end{matrix}] \in ℝ^{(n_{ℳ} \cdot n_{modes}) \times n_{α}}

, where each

D_{α} r_{i} \in ℝ^{n_{ℳ} \times n_{α}}

can be formed as follows

D_{α} r_{i} = [\begin{matrix} - \frac{D_{α} (λ_{i} (α))}{λ_{i}^{EXP}} \cdot w_{λ_{i}} \\ - D_{α} (ψ_{- q_{i}, i}^{ℳ} (α)) \cdot w_{ψ_{i}} \end{matrix}], i = 1 \dots n_{modes}

(22)

The Jacobian of the i-th simulated eigenvalue and eigenvector at measured DOFs, $D_{α} (λ_{i} (α))$ and $D_{α} (ψ_{- q_{i}, i}^{ℳ} (α))$ , have been provided in equations (13) and (20).

Jacobian of modal dynamic residual formulation

To derive the Jacobian of the modal dynamic residual formulation in equation (6), define a vector variable $x = [\begin{matrix} α; & ψ_{1}^{U}; & \begin{matrix} \dots \end{matrix} & ψ_{n_{modes}}^{U} \end{matrix}] \in ℝ^{n_{x}}$ with $n_{x} = n_{α} + n_{U} \cdot n_{modes}$ and a residual vector function $r (x)$ $: ℝ^{n_{x}} \to ℝ^{N \cdot n_{m o d e s}}$ as

r (x) = [\begin{matrix} r_{1} (x) \\ ⋮ \\ r_{n_{modes}} (x) \end{matrix}] = [\begin{matrix} [K (α) - λ_{1}^{EXP} M] {\begin{matrix} ψ_{1}^{EXP, ℳ} \\ ψ_{1}^{U} \end{matrix}} \cdot w_{1} \\ ⋮ \\ [K (α) - λ_{n_{modes}}^{EXP} M] {\begin{matrix} ψ_{n_{modes}}^{EXP, ℳ} \\ ψ_{n_{modes}}^{U} \end{matrix}} \cdot w_{n_{modes}} \end{matrix}]

(23)

where

r_{i} (x) : ℝ^{n_{x}} \to ℝ^{N}

equals

[K (α) - λ_{i}^{EXP} M] {\begin{matrix} ψ_{i}^{EXP, ℳ} \\ ψ_{i}^{U} \end{matrix}} \cdot w_{i}, i = 1 \dots n_{modes}

. Again, the Jacobian for

f (x) = r {(x)}^{T} r (x)

D_{x} f \in ℝ^{1 \times n_{x}}

equals

D_{r} f \cdot D_{x} r

from the chain rule. The first Jacobian term is

D_{r} f = 2 r^{T} \in ℝ^{1 \times (N \cdot n_{modes})}

while the second term is

D_{x} r = [\begin{matrix} D_{x} r_{1}; & D_{x} r_{2}; & \dots & D_{x} r_{n_{modes}} \end{matrix}] \in ℝ^{(N \cdot n_{modes}) \times n_{x}}

, where each

D_{x} r_{i} \in ℝ^{N \times n_{x}}

can be expressed as follows

D_{x} r_{i} = [\begin{matrix} D_{α} r_{i} D_{ψ_{1}^{U}} r_{i} \dots D_{ψ_{n_{modes}}^{U}} r_{i} \end{matrix}], i = 1 \dots n_{modes}

(24)

In equation (24), $D_{α} r_{i} \in ℝ^{N \times n_{α}}$ can be derived considering $K (α) = K_{0} + \sum_{j = 1}^{n_{α}} α_{j} K_{j}$ from equation (1)

D_{α} r_{i} = [\begin{matrix} K_{1} {\begin{matrix} ψ_{i}^{EXP, ℳ} \\ ψ_{i}^{U} \end{matrix}} \cdot w_{i} & \begin{matrix} K_{2} {\begin{matrix} ψ_{i}^{EXP, ℳ} \\ ψ_{i}^{U} \end{matrix}} \cdot w_{i} & \dots \end{matrix} & K_{n_{α}} {\begin{matrix} ψ_{i}^{EXP, ℳ} \\ ψ_{i}^{U} \end{matrix}} \cdot w_{i} \end{matrix}], i = 1 \dots n_{modes}

(25)

and

D_{ψ_{j}^{U}} r_{i} \in ℝ^{N \times n_{U}}, j = 1 \dots n_{modes}

in equation (24) can be expressed as

D_{ψ_{j}^{U}} r_{i} = {\begin{matrix} [K (α) - λ_{i}^{EXP} M] {\begin{matrix} 0_{n_{ℳ} \times n_{U}} \\ I_{n_{U}} \end{matrix}} \cdot w_{i} for j = i \\ 0_{N \times n_{U}} for j \neq i \end{matrix}

(26)

where

0_{n_{ℳ} \times n_{U}}

and

0_{N \times n_{U}}

denote zero matrices and

I_{n_{U}}

denotes the identity matrix, respectively.

Introduction of SMU package

Figure 1 shows the overall flowchart of the SMU package for structural model updating (Otsuki et al., 2021b). The GitHub package contains MATLAB subroutines for model updating with example codes for several structures. The main function StructModelUpdating performs a single run of structural model updating using one of the three formulations selected by users. Analytical gradients for all three formulations are implemented in SMU. In order to increase the chance of finding a better local minimum for these nonconvex optimization problems, randomly generated starting points can be assigned by MultiRunModelUpdating. This script calls StructModelUpdating to conduct multiple runs of model updating, that is, performs gradient searches from all the random starting points. To solve each corresponding optimization problem from one assigned starting point in the feasible region of variables, SMU uses the optimization solvers lsqnonlin and fmincon available in the MATLAB optimization toolbox (MathWorks Inc., 2019). The GitHub package provides a number of structural examples, currently including a 4-DOF shear model, an 18-DOF shear model, a steel pedestrian bridge, and a concrete building frame. Model updating is performed for structures using both simulation and experiment data. Additionally, when users want to perform model updating of their structure using the SMU package, users only need to provide the mass and stiffness information of the initial FE model of the structure along with experimentally obtained eigenvalues and eigenvectors for a certain number of modes.

Figure 1.

Overall flowchart of SMU.

Figure 2 shows the recommended optimization solver and algorithm for the model updating formulations. Note that the computational time of the FE model updating depends not only on the selection of optimization solvers and algorithms but also the size of the optimization problem, the use of numerical or analytical Jacobian, and the selection of formulations, which will be illustrated in two example structures in the Example 1: 18-story steel frame section and the Example 2: concrete building frame section. When minimizing least-squares objective functions, lsqnonlin is computationally much more efficient than fmincon and thus recommended for a first try. The lsqnonlin solver minimizes the objective function $f (x) = r {(x)}^{T} r (x) = {‖ r (x) ‖}_{2}^{2}$ in equation (7), that is, the square of the $ℒ_{2}$ -norm of a residual vector function $r (x)$ . The gradient, $\nabla f (x) \in ℝ^{n_{x}}$ and Hessian, $\nabla^{2} f (x) \in ℝ^{n_{x} \times n_{x}}$ of $f (x) : ℝ^{n_{x}} \to ℝ$ can be expressed as follows

\nabla f (x) = 2 {[D_{x} r]}^{T} \cdot r (x)

(27)

\nabla^{2} f (x) = 2 {[D_{x} r]}^{T} \cdot D_{x} r + 2 \sum_{i = 1}^{n_{r}} r_{i} (x) \nabla^{2} r_{i}

(28)

where

D_{x} r \in ℝ^{n_{r} \times n_{x}}

is defined as the Jacobian matrix of the scalar residuals

(r_{i}, i = 1 \dots n_{r})

with respect to the optimization variables

(x_{j}, j = 1 \dots n_{x})

. Neglecting the higher-order second term in

\nabla^{2} f (x)

, the optimization algorithms adopted by lsqnonlin in MATLAB uses

2 {[D_{x} r]}^{T} \cdot D_{x} r

to approximate the Hessian matrix. At each step of the optimization process, by default, lsqnonlin numerically calculates the search gradient,

\nabla f (x)

, of the objective function using the finite forward difference method (LeVeque, 2007). The numerically calculated gradient results are affected by the difference

\nabla x

, that is, the step size of x, and are more prone to inaccuracies. Meanwhile, instead of using the numerically calculated gradient, lsqnonlin also accepts a user-provided analytical formulation of the gradient. Given that the gradient simply equals the transpose of the Jacobian, that is,

\nabla f (x) = {(D_{x} f)}^{T}

, the analytical gradient for each formulation has been derived in the Jacobian derivatives of the formulations section and implemented in SMU.

Figure 2.

Recommended MATLAB solvers and optimization algorithms for each formulation.

At each iteration, the optimization solver evaluates the objective function value once or several times at points near the current value of $x$ , which is called the number of function evaluations or F-count in the MATLAB optimization toolbox (MathWorks Inc., 2019). When using the default numerical gradient, that is, finite forward difference method, each iteration needs to evaluate the objective function $n_{x} + 1$ times in the neighbor of the current $x$ . On the other hand, when the formulation of the analytical gradient is provided, the function evaluation happens only once at each iteration. Therefore, when the computation of the analytical gradient itself is not expensive, the analytical gradient can not only provide better accuracy but also require less computational time due to the smaller number of function evaluations. The advantage of using the analytical gradient is more significant when $n_{x}$ is large, that is, a large-scale optimization problem.

Two optimization algorithms available in lsqnonlin are the Levenberg–Marquardt (L–M) algorithm (Moré, 1978) and the trust-region-reflective (TRR) algorithm (Coleman and Li, 1996). The L–M algorithm is a combination of the steepest decent and Gauss-Newton algorithms and is specially designed to solve nonlinear least-squares problems. At every iteration, the algorithm first linearizes the objective function with respect to the corresponding optimization variables. When the current solution is far from a local optimum, the L–M algorithm approaches the steepest descent algorithm. On the other hand, when the current solution is close to a local optimum, the L–M algorithm approaches the Gauss–Newton algorithm. The TRR algorithm approximates the original problem with a quadratic subproblem within a small region around the current solution point, that is, a trusted region. The quadratic subproblem is formulated using the same gradient and approximated Hessian of the original problem. By solving the quadratic subproblem using the two-dimensional subspace approach, a solution of the current subproblem can be obtained (Branch et al., 1999; Byrd et al., 1988). If the decrease of the objective function evaluated at the current step is within the prescribed upper and lower bounds, the solution will be accepted, and the algorithm will continue with the next iteration. Otherwise, the trusted region at the current iteration will be adjusted, and the quadratic subproblem is solved again with the new region. Iteratively, the optimization converges to a local minimum of the objective function.

The L–M algorithm in the current version of SMU does not allow to set upper/lower bounds of optimization variables. This can be a drawback in model updating using experimental data because appropriate upper/lower bounds ensure a physically meaningful optimal solution. The TRR algorithm in MATLAB does allow setting the upper/lower bounds of optimization variables; however, the TRR algorithm is not applicable for underdetermined problems, that is, when $n_{r} < n_{x}$ . The MAC value approach for model updating (equation (3)) oftentimes provides an underdetermined problem. For these reasons, we recommend to first try the L–M algorithm for the MAC value approach. Whereas for the eigenvector difference approach (equation (5)) and the modal dynamic residual formulation (equation (6)), we recommend to first try the TRR algorithm.

All the model updating formulations (equations (3), (5), and (6)) in general provide nonconvex optimization problems. If the optimization problem is nonconvex, none of the existing algorithms in the MATLAB optimization toolbox can guarantee the global optimality. Therefore, besides starting the search from random initial points in the feasible region (MultiRunModelUpdating), it is also recommended to use different optimization algorithms and compare the model updating results. Notice that the optimization algorithms available in lsqnonlin and fmincon are different. Users may select fmincon when there is a particular interest in one of the algorithms available in fmincon, such as the interior-point method, the sequential-quadratic-programming algorithm, and the active-set method (MathWorks Inc., 2019).

Example 1: 18-story steel frame

Overview of the 18-story steel frame example

This section demonstrates an application example of SMU, using a numerical model based on a one-third scale 18-story steel frame. The physical structure was constructed and tested at the E-Defense in 2014 (Figure 3(a)) (Suita et al., 2015). The structure’s behavior represents that of typical steel high-rise buildings constructed in the 1980s to 1990s in Japan. The plane dimension of each floor is 5 $\times$ 6 m, and the total height is 25.35 m. In this numerical study, the overall structure is simplified to an 18-DOF shear model as shown in Figure 3(b).

Figure 3.

18-story steel frame. (a) Front view. (b) 18-DOF model.

The structural properties for this numerical study are presented in Table 1. It is assumed that the floor masses are accurate and do not require updating. The initial stiffness values are calculated from nominal material properties. The “actual” structure in this simulation corresponds to the assigned stiffness updating variable

α_{j}^{act}

for each story, which is to be identified through model updating. The value of each updating variable represents the relative change of a stiffness parameter from its initial/nominal value. The number of measured DOFs is assumed to be seven, corresponding to the 1st, 4th, 7th, 9th, 12th, 15th, and 18th DOFs. Therefore,

n_{ℳ}

= 7 and

n_{U}

= 11. The number of “experimental” modes available in this model updating is set as four

(n_{modes} = 4)

. For the MAC value and the eigenvector difference approaches, the number of optimization variables is

n_{x} = n_{α} =

18. For the modal dynamic residual formulation, the number is

n_{x} = n_{α} + n_{U} \times n_{modes} = 18 + 11 \times 4 = 62

Table 1.

Structural properties of the 18-DOF model for this numerical study.

Floor	Weight (kN)	Initial inter-story stiffness (kN/m)	Actual inter-story stiffness (kN/m)	Updating variable	Actual updating variable $α_{j}^{act}$
18*	202	36,300	43,560	$α_{18}$	0.20
17	206	49,100	58,920	$α_{17}$	0.20
16	206	56,200	61,820	$α_{16}$	0.10
15*	206	61,900	52,615	$α_{15}$	−0.15
14	206	66,000	69,300	$α_{14}$	0.05
13	206	71,200	60,520	$α_{13}$	−0.15
12*	206	78,800	98,500	$α_{12}$	0.25
11	208	82,400	107,120	$α_{11}$	0.30
10	208	84,000	100,800	$α_{10}$	0.20
9*	208	87,600	78,840	$α_{9}$	−0.10
8	208	93,800	117,250	$α_{8}$	0.25
7*	208	96,300	110,745	$α_{7}$	0.15
6	208	99,000	84,150	$α_{6}$	−0.15
5	208	102,800	113,080	$α_{5}$	0.10
4*	208	102,800	92,520	$α_{4}$	−0.10
3	208	107,300	101,935	$α_{3}$	−0.05
2	208	109,200	114,660	$α_{2}$	0.05
1*	208	115,500	121,275	$α_{1}$	0.05

Note. (* measured/instrumented DOFs).

Procedure illustration for structural model updating using SMU

This section demonstrates the procedures for updating the 18-DOF model. The SMU MATLAB package can be downloaded from GitHub (Wang et al., 2019b). The folder “Structural-Model-Updating\SMU” contains shared MATLAB subroutines for structural model updating. This folder should be added into the MATLAB path first. The MATLAB code that contains structural properties of this 18-DOF simulation study can be found in “Structural-Model-Updating\Examples\EighteenStoryStructure\Simulation”. Figure 4 shows a part of the main script for model updating of the 18-DOF model. The script at first sets the basic parameters of the structure. As mentioned, this example assumes four “experimental” modes $(n_{modes} = 4)$ are available; the measured DOFs are identified. These parameters can be easily changed by users to see the effect on model updating results. The mass matrix $M$ , initial stiffness matrix $K_{0}$ , and influence matrix $K_{j}$ in equation (1) are constructed based on properties in Table 1. Upon construction or loading from a data file, they are assembled in a MATLAB structure array structModel. Experimental eigenvalues $λ_{i}^{EXP}$ and eigenvectors $ψ_{i}^{EXP, ℳ}$ used in the model updating formulations—see equations (3), (5), and (6)—are then generated using the “actual” values of the stiffness updating variables as listed in Table 1, and they are assigned to variables lambdaExp and psiExp_m, containing $λ_{i}^{EXP}$ and $ψ_{i}^{EXP, ℳ}$ , respectively. The script next calls the SMU subroutines to perform structural model updating. Before calling the subroutines, options for optimization algorithms can be specified by users. The example below uses lsqnonlin and the L–M algorithm with analytical gradient.

Figure 4.

Main script for model updating of the 18-DOF model in simulation.

The formID = 1 means that the modal property difference formulation with the MAC value approach in equation (3) is used. For the matching method between simulated and experimental modes, modeMatch and modeIndex are assigned. When modeMatch is set as 1, at every iteration step, the program goes through remaining simulated modes specified by modeIndex and pairs each experimental mode with a remaining simulated mode with the largest MAC value; the identified pairs of $ψ_{i}^{EXP, ℳ}$ and $ψ_{i}^{ℳ}$ are then used to evaluate the objective function and search gradient. When modeMatch is set as 2, the program strictly matches the first experimental mode with the first simulated mode specified by modeIndex, the second experimental mode with the second one in modeIndex, etc. The modeMatch = 2 setting should be used only when users are confident that all the lowest few modes are captured by experimental data, that is, there is no missing or unmeasured/undetected mode from the experimental data.

In this example, the formID = 1 is selected to use the MAC value formulation in equation (3). The upper and lower bounds (x_lb and x_ub) are set as $\pm$ 1 for the stiffness updating variables $α$ . In case the modal dynamic residual formulation in equation (6) is used (formID = 3), the upper and lower bounds are set as $\pm$ 1 for the stiffness updating variables $α$ and $\pm 2$ for the unmeasured entries of eigenvectors $ψ_{i}^{U}$ . Note that the eigenvectors in simulation and experiment should have been normalized so that the largest entry of $ψ_{i}^{EXP, ℳ}$ equals 1. In this demonstration, the weighting factors are all equal to 1 and are assigned to the corresponding field of the structure array expModes, including $w_{λ_{i}}$ , $w_{ψ_{i}}$ , and $w_{i}$ . In practice, depending on the confidence in experimental data and the knowledge of the structure’s dynamics, appropriate weighting factors should be specified. Next, the options for MultiRunModelUpdating are specified. The 100 randomly generated starting points in the feasible regions are used. The variable randSeed is the random seed value for the MATLAB function “rng”. The variable filename will be used for the result file name (.mat file) that contains model updating results for all the starting points. Upon running the code, the iteration history will be shown in the MATLAB Command Window, and the result file (.mat file) will be automatically generated. For this study, we use a PC with an Intel i7-7700 CPU and 16 GB RAM.

Model updating results: MAC value approach

For this simulation, the MAC value approach in equation (3) provides $n_{r} = 2 \cdot n_{modes} = 2 \cdot 4 = 8$ and $n_{x} = n_{α} = 18$ . Hence, $n_{r} < n_{x}$ , that is, the optimization problem is underdetermined; the TRR algorithm is not applicable, and the L–M algorithm is used to solve the optimization problem. When using MATLAB lsqnonlin with the L–M algorithm in the current version of SMU, the optimal result sets that are out of the assigned bounds are rejected and the corresponding starting point is replaced with the next randomly generated point that achieves valid optimal results. For comparison, both numerical and analytical gradients are used. For each optimum parameter set $α^{*}$ , the relative error of the stiffness updating variables is calculated as

e_{i} = \frac{| α_{i}^{*} - α_{i}^{act} |}{1 + α_{i}^{act}} \times 100 %, i = 1 \dots n_{α}

(29)

where

α_{i}^{*}

is the i-th entry of

α^{*}

and

α_{i}^{act}

is the actual value as listed in Table 1. Then, the average relative error

e_{avg}

is calculated for this optimum parameter set

α^{*}

e_{avg} = \frac{1}{n_{α}} \sum_{i = 1}^{n_{α}} e_{i}

(30)

Figure 5(a) plots the $e_{avg}$ of the optimum $α^{*}$ from all of the 100 starting points. In this plot, many of the optimization searches end with a solution that has a large relative error $e_{avg}$ ; some relative errors are close to or even above 30%. The reason is that in this simulation, the maximum number of search iterations from each starting point is set as the default value of 400 by SMU. If desired, the number can be easily increased for better accuracy while consuming longer runtime. Among 100 sets of the optimum vector variable $α^{*}$ , the best parameter set is chosen as the one corresponding to the minimum objective function value. Figure 5(b) shows the relative error of this best parameter set. Both numerical and analytical gradients achieve a solution close to the actual value of updating variables, but the analytical gradient obtains better accuracy. The total computational time is 475.8 s when using the numerical gradient while it is 386.9 s for the analytical gradient.

Figure 5.

Model updating results of 18-DOF model by MAC value approach using L–M algorithm with numerical and analytical gradients. (a) Average relative errors of updating variables. (b) Relative errors of updated stiffnesses.

Model updating results: Eigenvector difference approach

The eigenvector difference approach in equation (5) provides $n_{r} = n_{ℳ} \cdot n_{modes} = 7 \cdot 4 = 28$ and $n_{x} = n_{α} = 18$ . Hence, $n_{r} > n_{x}$ , that is, the optimization problem is determined. Therefore, the TRR algorithm is applicable and used as recommended in Figure 2. Both numerical and analytical gradients are studied for comparison. Figure 6(a) shows the average relative error for each starting point. Figure 6(b) plots the relative error of the stiffness values for the minimum objective function value among 100 optimal sets. It can be concluded that for this relatively simple example, the eigenvector difference approach, using both numerical and analytical gradients, is able to identify accurate stiffness values of the model. The total computational time is 49.2 s when using the numerical gradient and 23.6 s for the analytical gradient. Compared to the MAC value approach, the eigenvector difference approach achieves higher accuracy and much better computational efficiency.

Figure 6.

Model updating results of 18-DOF model by eigenvector difference approach using TRR algorithm with numerical and analytical gradients. (a) Average relative errors of updating variables. (b) Relative errors of updated stiffnesses.

Model updating results: Modal dynamic residual formulation

The modal dynamic residual formulation in equation (6) provides $n_{r} = n_{N} \cdot n_{modes} = 18 \cdot 4 = 72$ and $n_{x} = n_{α} + n_{U} \cdot n_{modes} = 18 + 11 \cdot 4 = 62$ . Hence, $n_{r} > n_{x}$ , that is, the optimization problem is determined. To compare the different optimization algorithms, the L–M algorithm and the TRR algorithm are both used with the analytical gradient. Figure 7(a) shows the average relative error of the updated result from each starting point. Compared to the L–M algorithm, several optimization searches using the TRR algorithm fail to converge at the correct value. Figure 7(b) plots the best solution among 100 optimal sets. The accuracy of the updated results is better than the MAC value approach and nearly the same as the eigenvector difference approach. The total computational time is 265.0 s for the L–M algorithm and 45.3 s for the TRR algorithm.

Figure 7.

Model updating results of 18-DOF model by modal dynamic residual formulation using L–M and TRR algorithm with analytical gradient. (a) Average relative errors of updating variables. (b) Relative errors of updated stiffnesses.

In general, with eigenvectors $ψ_{i}^{U}$ among the optimization variables, the modal dynamic residual formulation can be computationally expensive when the number of unmeasured DOFs is large. In addition, as aforementioned, because the formulation does not directly minimize the discrepancy of modal properties, the updated model may provide modal properties that are more noticeably different from experimental data. Our overall experience is that the modal dynamic residual formulation usually does not perform well with larger-scale experimental data. Therefore, while for completeness the formulation is investigated in this section, it will not be studied further in the next section.

Example 2: concrete building frame

Overview of the concrete building frame example

The second structure studied in this research is a concrete building frame (Figure 8(a)), which simulates a full-scale test structure in the Structural Engineering and Materials Laboratory on the Georgia Tech campus. The test frame structure is representative of low-rise reinforced concrete office buildings in the central and eastern United States built in the 1950s to 1970s (Dong et al., 2016). The structure consists of a set of four identical two-story two-bay concrete frames. Each frame was constructed with a gap from its neighboring frames, allowing free in-plane longitudinal movement and can thus be modeled independently from the other frames. Frame #1 is used in this study. The columns and beams for frame #1 are modeled with frame elements in SAP2000 (Figure 8(b)).

Figure 8.

Concrete building frame (height in z: 2 × 3.66 m; length in x: 2 × 5.49 m; slab width in y: 2.74 m). (a) Front view of the four testing frames, enveloped by two collapse-prevention frames. (b) Model of the 2-story 2-bay concrete building frame and sensor instrumentation.

Corresponding to dense sensor instrumentation, seven segments are allocated per column on each story, and 12 segments per beam in each bay. In SAP2000 to ensure a stiffness contribution from both concrete and rebar, along every column or beam segment one frame element is assigned for the concrete material and another frame element is assigned for the steel reinforcement. Each floor slab is meshed into 175 shell elements. In total, the FE model of the concrete building frame has 2302 DOFs. The modeling software SAP2000 assigns non-zero concentrated mass only to the translational DOFs. As a result, the mass matrix (M) is a diagonal matrix whose diagonal entries associated with rotational DOFs equal zero.

Figure 8(b) also shows the accelerometer instrumentation for this simulation study, and the corresponding measurement directions. A total of 43 DOFs are measured, that is, the length of $ψ_{i}^{EXP, ℳ}$ and $ψ_{i}^{ℳ} (α)$ is $n_{ℳ} = 43$ . The accelerometers measure longitudinal and vertical vibration, that is, x and z directions. Thus, only in-plane vibration mode shapes, that is, in the x-z plane, can be extracted from measurement data. To avoid the side effect of out-of-plane mode shapes (in the y-z plane) on the FE model updating, the vertical and transverse DOFs (y and z direction) at both ends of the three transverse beams (along the y direction) on each slab are constrained. Lastly, at the bottom of the three columns below the first slab, all six DOFs are constrained to represent an ideal fixed support.

This study updates 12 elastic moduli of concrete members, denoted as E₁ ∼ E₁₂ in Figure 8(b). In the first story, E₁ ∼ E₃ represent the concrete elastic moduli of members in the three columns; E₇ and E₈ represent the concrete elastic moduli of longitudinal beam members (along the x direction); and E₁₁ represents the concrete elastic moduli of the first slab and the associated transverse beam members (along the y direction). Similarly, other moduli for the second story can be found in Figure 8(b). While this study only involves simulation, the selection of moduli corresponds to different concrete pours during the construction, and thus is in preparation for future model updating of the as-built structure with experimental data. Compared to concrete, the estimated elastic modulus of steel reinforcement is considered to be sufficiently accurate; therefore, it is not being updated in this study.

For all the concrete moduli being updated, Table 2 lists the nominal and actual values. In total, there are 12 updating variables for this model updating, that is,

n_{α} = 12

. The column

α_{i}^{act}

in Table 2 lists the actual values of

α

, that is, the ideal solutions to be identified through FE model updating. For the model updating of this concrete structure, this paper presents the results of the modal property difference formulations with the MAC value approach (equation (3)) and the eigenvector difference approach (equation (5)). It is assumed that the first three vibration modes

(n_{modes} = 3)

are available. The upper and lower bounds of

α

are set to be 1 and −1, respectively. The weightings are set as

w_{λ_{i}} =

1 and

w_{ψ_{i}} = 1

, respectively. The optimization process is initiated from 100 random starting points within the bounds of

α

. The maximum number of iterations is set as 1000. Same as the simulation for the 18-DOF model, a PC with an Intel i7-7700 CPU and 16 GB RAM is used.

Table 2.

Structural properties of the concrete building frame model for this numerical study.

Stiffness parameters		Nominal value (N/mm²)	Actual value (N/mm²)	Updating variables	Actual updating variable $α_{j}^{act}$
Elastic moduli of concrete members	E ₁	26,890	24,201	$α_{1}$	−0.10
	E ₂	25,511	30,613	$α_{2}$	0.20
	E ₃	25,511	30,613	$α_{3}$	0.20
	E ₄	22,063	20,960	$α_{4}$	−0.05
	E ₅	22,063	26,476	$α_{5}$	0.20
	E ₆	22,063	25,373	$α_{6}$	0.15
	E ₇	22,063	25,373	$α_{7}$	0.15
	E ₈	22,063	24,270	$α_{8}$	0.10
	E ₉	23,442	21,098	$α_{9}$	−0.10
	E ₁₀	23,442	19,926	$α_{10}$	−0.15
	E ₁₁	22,063	26,476	$α_{11}$	0.20
	E ₁₂	23,442	26,958	$α_{12}$	0.15

Model updating results: MAC value approach

The MAC value approach (equation (3)) in this simulation provides $n_{r} = 2 \cdot n_{modes} = 2 \cdot 3 = 6$ while $n_{x} = n_{α} = 12$ . As a result, $n_{r} < n_{x}$ and the problem is underdetermined; the TRR algorithm is not applicable and the L–M algorithm is used. Recall that the L–M algorithm in the current version of SMU does not allow upper/lower bounds. Therefore, when the converged optimal solution is out of the range [−1,1], the solution was discarded, and a new starting point was assigned. For comparison, both numerical and analytical gradients are used during optimization searches. For each of the 100 successful runs, Figure 9(a) displays the average relative error e_avg (equation (30)) for all $n_{α} = 12$ stiffness updating variables after discarding 10 optimal result sets that are out of the [−1,1] bounds during the optimization process. Similar to the case for the MAC value approach of the 18-DOF model, a large number of optimization searches end with a large relative error $e_{avg}$ . For the solution set that achieves the minimum objective function value among the 100 starting points, Figure 9(b) plots the relative error $e_{i}$ of each stiffness parameter. The figure shows that the obtained stiffness parameter values are not reasonable, with the maximum relative error larger than 13% using the numerical gradient and 6.5% using the analytical gradient. It took 18 h and 58 min to obtain the 100 optimal result sets that are within the bounds when using the numerical gradient, while it took 22 h and 58 min when using the analytical gradient. While the analytical gradient provides more accurate solution, due to the computation cost of the analytical gradient for the MAC value approach, the numerical gradient required less computational time than the analytical gradient in this example. In conclusion, for this concrete frame structure, the MAC value approach cannot provide reasonable FE model updating results using either the numerical gradient or analytical gradient during the optimization process.

Figure 9.

Model updating results for the concrete building frame by MAC value approach using the L–M algorithm with numerical and analytical gradients. (a) Average relative errors of updating variables. (b) Relative errors of updated stiffness parameters.

Further study on the nonconvexity of the MAC value approach is conducted, using solution A in Figure 9(a). The average relative error of solution A is 4.71% and the objective function equals 4.71×10⁻¹⁰. While it is impossible to visualize $f (α)$ in 12-dimensional space, Figure 10 displays a walk from solution A, denoted as $α_{A}$ , to the actual value $α^{act}$ shown in Table 2, along a hyperline in $ℝ^{12}$ space. The hyperline is defined with a scalar $θ$ ; along the hyperline, the values of updating variables are $α (θ) = α_{A} + θ \cdot (α^{act} - α_{A})$ . When $θ = 0$ , $α = α_{A}$ ; when $θ = 1$ , $α = α^{act}$ . Accordingly, the y-axis represents the MAC value objective function (equation (3)) evaluated at different values of $α$ along the hyperline

f (α (θ)) = f (α_{A} + θ \cdot (α^{act} - α_{A}))

(31)

Figure 10.

Hyperline walk between solution A to the global optimum.

Had $f (α)$ been convex on $α$ , $f (α (θ))$ should also be convex on $θ$ (see page 68 in the 7th Ed. of reference (Boyd and Vandenberghe, 2004)). Figure 10 clearly shows that there are two valleys along the hyperline, located at $θ = 0$ and $1$ , respectively. Therefore, the optimization problem in equation (3) for the MAC value approach is confirmed to be nonconvex, and off-the-shelf local optimization algorithms cannot ensure global optimality.

Model updating results: eigenvector difference approach

For the eigenvector difference formulation in equation (5), the residual vector length $n_{r}$ equals $n_{ℳ} \cdot n_{modes} = 43 \cdot 3 = 129$ . As a result, with the number of optimization variables $n_{x}$ still equal to $n_{α} = 12$ , the problem is determined and TRR algorithm is used. The numerical and analytical gradients are utilized for comparison. Figure 11(a) plots the average relative error $e_{avg}$ of the stiffness updating variables for 100 random searches. Many starting points end with a large error using the numerical gradient while only three starting points converged with large errors using the analytical gradient. For each gradient case, the best solution is again selected as the one with the minimum objective function value among 100 result sets. The relative errors of the updated stiffness parameters for the best solution are plotted in Figure 11(b). Both the numerical and analytical gradient provide accurate updating results, and the analytical gradient provides no errors. To obtain the 100 optimal solutions, it took 22 h 38 min using the numerical gradient while it took only 52 min using the analytical gradient.

Figure 11.

Model updating results for the concrete building frame by eigenvector difference approach using the TRR algorithm with numerical and analytical gradients. (a) Average relative errors of updating variables. (b) Relative errors of updated stiffness parameters.

Similar to the MAC value approach, the nonconvexity of the eigenvector difference approach is studied using solution B in Figure 11(a). The average relative error of solution B is 20.81% and the objective function equals 0.0054. Figure 12 displays a walk from solution B to the global optimum along a hyperline in $ℝ^{12}$ space. Figure 12 clearly shows two valleys along the hyperline, indicating the nonconvexity of the eigenvector difference approach in equation (5). Hence, global optimality cannot be ensured using off-the-shelf local optimization algorithms and multiple random starting points are necessary to obtain a better solution.

Figure 12.

Hyperline walk between solution B to the global optimum.

In summary, for model updating of the concrete building frame, the MAC value approach could not provide a set of optimal results with acceptable accuracy using either the numerical or analytical gradient during the optimization process. As for the eigenvector difference formulation, using the analytical gradient not only provides more accurate model updating results in general but also can find the correct updating parameter values more efficiently. The study demonstrates the advantage of using the analytical gradient versus the numerical gradient.

Conclusions

An open-source MATLAB package for structural model updating called “SMU” is developed and shared with the research community. Current version 1.1 supports three model updating formulations in frequency domain, namely, (i) the MAC (modal assurance criterion) value approach (equation (3)), (ii) the eigenvector difference approach (equation (5)), and (iii) the modal dynamic residual formulation (equation (6)). The first two belong to the family of modal property difference formulations. To facilitate the optimization process, analytical gradients of the formulations are derived and implemented in SMU. Considering that the objective functions of the model updating formulations are nonconvex, randomly generated starting points are adopted to increase the chance of finding the global minimum.

Among several structural examples provided in SMU, this paper demonstrates the SMU functionality using an 18-DOF model and a concrete building frame model in simulation. The MAC value approach and the modal dynamic residual formulation can correctly identify stiffness values for a relatively simpler structural model (the 18-DOF model), but may fail to provide reasonable results and increase computational cost when the complexity of a structure increases (as in the concrete building frame). The eigenvector difference approach can correctly update stiffness values for both simpler and more complex structures with better computational efficiency than the MAC value approach and the modal dynamic residual formulation. Furthermore, this study demonstrates the advantage of using the analytical gradient versus the numerical gradient. Using the analytical gradient during the optimization process in general not only provides more accurate model updating results but also saves computational time, especially for a large-scale structure. For the concrete building frame, the computational time using the analytical gradient for the eigenvector difference approach is 52 min while it is 22 h 38 min using the numerical gradient. For future work, different formulations and more structural examples are to be added. Moreover, the latest MATLAB 2021b version changed the Levenberg–Marquardt algorithm to support the bounds of optimization variables, and this change will be reflected in the SMU package in the future.

Footnotes

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This research was partially funded by the National Science Foundation (CMMI-1150700 and CMMI-1634483). The first author received scholarship support from the Nakajima Foundation. Any opinions, findings, and conclusions expressed in this publication are those of the authors and do not necessarily reflect the view of the sponsors.

ORCID iDs

Peter Lander

Yang Wang

References

Bakir

Reynders

De Roeck

(2008) An improved finite element model updating method by the global optimization technique ‘Coupled Local Minimizers’. Computers & Structures 86(11–12): 1339–1352.

Boyd

Vandenberghe

(2004) Convex Optimization. Cambridge, UK: Cambridge University Press.

Branch

Coleman

(1999) A subspace, interior, and conjugate gradient method for large-scale bound-constrained minimization problems. SIAM Journal on Scientific Computing 21(1): 1–23.

Byrd

Schnabel

Shultz

(1988) Approximate solution of the trust region problem by minimization over two-dimensional subspaces. Mathematical Programming 40(1–3): 247–263.

Chatzi

Smyth

(2009) The unscented Kalman filter and particle filter methods for nonlinear structural system identification with non‐collocated heterogeneous sensing. Structural Control and Health Monitoring 16(1): 99–123.

Coleman

(1996) An interior trust region approach for nonlinear minimization subject to bounds. SIAM Journal on Optimization 6(2): 418–445.

Dong

Liu

Wright

, et al. (2016) Validation of wireless sensing technology densely instrumented on a full-scale concrete frame structure. In: Proceedings of international conference on smart infrastructure and construction (ICSIC), Cambridge, UK, June 27‐29, pp. 143–148.

Dong

Wang

(2018) Modal property difference formulations and optimization algorithm comparison towards FE model updating. In: Proceedings of SPIE 2018, smart structures and materials + nondestructive evaluation and health monitoring. Denver, CO, March 4‐8, pp. 1059828.

Dong

Wang

(2019) Finite element model updating of a steel pedestrian bridge model. In: ASCE international conference on computing in civil engineering, Atlanta, GA, USA, June 17‐19, pp. 397–404.

10.

Doucet

Tadić

(2003) Parameter estimation in general state-space models using particle methods. Annals of the Institute of Statistical Mathematics 55(2): 409–422.

11.

Farhat

Hemez

(1993) Updating finite element dynamic models using an element-by-element sensitivity methodology. AIAA Journal 31(9): 1702–1711.

12.

Fox

Kapoor

(1968) Rates of change of eigenvalues and eigenvectors. AIAA Journal 6(12): 2426–2429.

13.

Friswell

Mottershead

(1995) Finite Element Model Updating in Structural Dynamics. Dordrecht, the Netherlands: Kluwer Academic Publishers.

14.

Hofmeister

Bruns

Rolfes

(2019) Finite element model updating using deterministic optimisation: a global pattern search approach. Engineering Structures 195: 373–381.

15.

Kane

Zhu

Hirose

, et al. (2014) Development of an extensible dual-core wireless sensing node for cyber-physical systems. In: Proceedings of SPIE 2014, nondestructive characterization for composite materials, aerospace engineering, civil infrastructure, and homeland security, San Diego, CA, March 10‐14, pp. 90611U.

16.

Kosmatka

Ricles

(1999) Damage detection in structures by modal vibration characterization. ASCE Journal of Structural Engineering 125(12): 1384–1392.

17.

LeVeque

(2007) Finite Difference Methods for Ordinary and Partial Differential Equations: Steady-State and Time-dependent Problems. Philadelphia, PA: SIAM.

18.

Dong

Wang

(2018) Model updating using sum of squares (SOS) optimization to minimize modal dynamic residuals. Structural Control and Health Monitoring 25(12): e2263.

19.

Lynch

Loh

(2006) A summary review of wireless sensors and sensor networks for structural health monitoring. The Shock and Vibration Digest 38(2): 91–128.

20.

MathWorks Inc (2019) Optimization Toolbox™ User’s Guide. Natick, MA, USA: MathWorks Inc.

21.

Moller

Friberg

(1998) Updating large finite element models in structural dynamics. AIAA Journal 36(10): 1861–1868.

22.

Moré

(1978) The Levenberg-Marquardt algorithm: implementation and theory. In: Numerical Analysis. Berlin, Germany: Springer pp. 105–116.

23.

Nelson

(1976) Simplified calculation of eigenvector derivatives. AIAA Journal 14(9): 1201–1205.

24.

Otsuki

Dey

, et al. (2021a) Finite element model updating of an 18-story structure using branch-and-bound algorithm with epsilon-constraint. Journal of Civil Structural Health Monitoring 11: 575–592. DOI: 10.1007/s13349-020-00468-3.

25.

Otsuki

Xinjun

, et al. (2021b) SMU – an open-source MATLAB package for structural model updating. In: Proceedings of the 10th international conference on bridge maintenance, safety and management (IABMAS), Sapporo, Japan, April 11‐18, pp. 1621–1628.

26.

Salawu

(1997) Detection of structural damage through changes in frequency: a review. Engineering Structures 19(9): 718–723.

27.

Sato

(1998) Adaptive H∞ filter: its application to structural identification. Journal of Engineering Mechanics 124(11): 1233–1240.

28.

Suita

Suzuki

Takahashi

(2015) Collapse behavior of an 18-story steel moment frame during a shaking table test. International Journal of High-Rise Buildings 4(3): 171–180.

29.

Teughels

Maeck

De Roeck

(2002) Damage assessment by FE model updating using damage functions. Computers & Structures 80(25): 1869–1879.

30.

Wang

Dong

(2019a) A non-convexity study in finite element model updating. In: Proceedings of the 9th international conference on structural health monitoring of intelligent infrastructure (SHMII-9), St. Louis, MO, August 4‐7, pp. 24–29.

31.

Wang

Dong

, et al. (2019b) SMU: MATLAB Package for Structural Model Updating, version 1.1. Available at: https://github.com/ywang-structures/Structural-Model-Updating

(accessed 6 October 2021)

32.

Weng

Xia

Y-L

, et al. (2011) Substructure based approach to finite element model updating. Computers & Structures 89(9–10): 772–782.

33.

Yang

Lin

Huangl

, et al. (2006) An adaptive extended Kalman filter for structural damage identification. Structural Control & Health Monitoring 13(4): 849–867.

34.

Yuen

(2012) Updating large models for mechanical systems using incomplete modal measurement. Mechanical Systems and Signal Processing 28: 297–308.

35.

Zhang

Chang

TYP

(2000) Finite element model updating for structures with parametric constraints. Earthquake Engineering & Structural Dynamics 29(7): 927–944.

36.

Zhang

Sun

(2021) Structural damage identification via physics-guided machine learning: a methodology integrating pattern recognition with finite element model updating. Structural Health Monitoring 20(4): 1675–1688. DOI: 10.1177/1475921720927488.

37.

Zhu

Dong

Wang

(2016) Substructure stiffness and mass updating through minimization of modal dynamic residuals. ASCE Journal of Engineering Mechanics 142(5): 04016013.

38.

Zhu

Tian

, et al. (2021) An enhanced substructure-based response sensitivity method for finite element model updating of large-scale structures. Mechanical Systems and Signal Processing 154: 107359.