# Construction of Stiffness and Flexibility for Substructure-Based Model Updating

Authors: Shun Weng; Hong-Ping Zhu; Yong Xia; Ling Ye; Xian-Yan Hu
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

In substructuring methods, the substructures are independently analyzed under free-free conditions. For a free-free substructure, its stiffness matrix is singular and rank deficient due to rigid body motion. The variables associated with the inverse of the stiffness matrix are not easy to be accurately determined in the usual manner. This study expands on the previous research on the substructuring methods by taking a deeper look at the analysis of a free-free substructure. A well-conditioned stiffness matrix is constructed for the analysis of a free-free structure. Some difficulties associated with the analysis of the free-free substructures can be solved in a simple and effective way. The substructural eigensolutions and eigensensitivity are solved from the well-conditioned stiffness matrix, other than the singular stiffness matrix. The proposed well-conditioned eigenequation is accurate and efficient to calculate the substructural eigensolutions and eigensensitivity. The properties addressed in this paper are not limited to be used for the analysis of a free-free substructure in many substructuring methods, and they are promising to be generalized to a range of analysis relevant to a free-free structure.


---

## Body

## 1. Introduction

In the past several decades, a large number of long-term structural health monitoring (SHM) systems have been designed and implemented worldwide on civil engineering structures such as large-scale bridges and high-rise buildings [1–4]. The accurate and efficient model updating and damage detection are significant for the long-term SHM systems. The substructuring methods have proved to be accurate and efficient for the analysis of large-scale structure, and they have been extensively utilized in a large number of applications, such as the model updating, system identification, and structural control [5]. The substructuring methods possess more advantages than the traditional global methods which analyze a structure as a whole. First, as the global structure is replaced by smaller and more manageable substructures, it is much easier and quicker to analyze the small system matrices. Second, the substructuring methods allow for the analysis of local parts. When the substructuring method is applied in model updating or damage identification, only one or more substructures are involved in an optimization procedure. The size of the model and the number of the uncertain parameters are much smaller than those of the global structure. Finally, in practical testing, the experimental instruments can be saved if it is necessary to measure the whole structure only for one or more substructures [5].In general, the substructuring approach can be utilized in the forward and inverse manners, respectively. The forward substructuring approach is frequently found in the eigenanalysis of a structure [5–14]. The partitioned substructures are analyzed independently to obtain their designated solutions, such as the substructural eigensolutions and eigensensitivity. The substructural solutions are then assembled to recover the solutions of the global structure by imposing constraints on the interfaces. On the other hand, the substructuring approach can be used in an inverse manner to disassemble the properties of the global structure to the substructure level by satisfying the constraints at the interfaces [15–21]. After eliminating the rigid body components, the independent substructures can be singled out to be used for the static analysis, dynamic analysis, nonlinear analysis, fatigue analysis, and so forth.The substructuring methods require dividing the global structure into independent free or fixed substructures. After partition, the substructures are usually analyzed independently under the free-free constraints. Since a free-free structure includes the rigid body motion, its stiffness matrixK is singular and rank deficient, and hence the inverse of the stiffness (K-1) does not exist [22, 23]. In consequence, the variables associated with the inverse of the singular stiffness matrix, such as the modal flexibility, residual flexibility, eigensolutions, and eigensensitivity, are not easy to be determined. Some researchers avoided the rigid body modes (zero-frequency modes) by introducing a small shift in the singular and rank-deficient stiffness matrix [23, 24]. This inevitably introduces some errors. For example, the rigid body modes computed by a shift eigensolver are not perfect 0.0 Hz frequencies, and the computed mode shapes are not “clean” rigid body modes. In consequence, the obtained deformational modes which are orthogonal to the “unclean” rigid body modes are inaccurate. The variables relating to the zero-frequency modes, such as the eigensolutions and modal flexibility, are thereby not accurate [24].When the modal flexibility of a free-free structure is required, the modal flexibility was sometimes computed by extracting the Moore-Penrose pseudoinverse of the stiffness matrix [6, 24]. Due to the numerical roundoffs, the frequencies and mode shapes of the rigid body modes are not accurate. This kind of analysis is not only computationally expensive, but also significantly sensitive to the rank condition when carried out in floating-point arithmetic [25].This paper addresses some frequently encountered difficulties associated with the analysis of the free-free substructures when the authors studied the substructuring methods in the previous research [5, 9–16]. A new full-rank stiffness matrix is proposed, which leads to a well-conditioned eigenequation. Based on the well-conditioned eigenequation, the substructural flexibility, residual flexibility, eigensolutions, and eigensensitivity of a free-free structure are solved in an effective and efficient way. The formulae proposed in this paper are not only useful for the analysis of a free-free substructure in many substructuring methods but also generally applicable in the analysis of a free-free structure.


## 2. Construction of Free-Free Stiffness and Flexibility Matrices

### 2.1. Basic Theory for Eigenanalysis

A structure withN degrees of freedom (DOF) has the eigenequation of
(1)K{ϕi}=λiM{ϕi},
where K and M are the stiffness and mass matrices. λi is the ith eigenvalue of the structure, and {ϕi} is the corresponding eigenvector. They are determined by the physical property of a structure, such as Young’s modulus, density, Poisson ratio, and geometric dimension. The eigensolutions of (1) consist of the eigenvalues Λ=Diag(λ1λ2⋯λN) and the corresponding eigenvectors Φ=[ϕ1ϕ2⋯ϕN]. The eigenvectors are mass-normalized such that they satisfy the following relation:
(2)ΦTKΦ=Λ,ΦTMΦ=I.
The stiffness matrix can be written by the mass-normalized eigenmodes as
(3)K=∑i=1Nλi(Mϕi)(Mϕi)T=MΦΛΦTM.
A flexibility matrix has a very straightforward physical interpretation: the displacement response caused by an applied unit load [26]. The flexibility matrix can also be written by the mass-normalized eigenmodes as
(4)F=∑i=1N1λiϕiϕiT=ΦΛ-1ΦT.
For a fixed structure, the stiffness matrix and flexibility matrix normally formed a dual inverse of each other as
(5)KF=I,F=K-1,K=F-1.
The displacement {x} of a free-free structure can be written as a superposition of the deformational and rigid body motions,
(6){x}={xd}+{xr}=Φd{q}+R{α},
where {xd} is the displacement due to the deformational motion and {xr} is the displacement due to the rigid body motion. Φd is the linear orthogonal deformational modes, R is the orthogonal rigid body modes, and {q} and {α} are the participation factors of the orthogonal modes. The subscripts “d” and “r” hereinafter, respectively, represent the variables associated with the deformational motion and rigid body motion.Accordingly, eigenequation (1) for a free-free structure has two kinds of eigenpairs.(1)
Nr zero eigenvalues pertaining to the rigid body motions: the associated eigenvectors span the null space of the stiffness matrix K, which contribute to the columns of R. Nr is equal to the number of the statically determinate constraints required to prevent all rigid body motion [24].(2)
Nd=N-Nr nonzero eigenvalues λi(i=1,2,…,Nd): the associated orthogonal deformational eigenvectors Φd=[ϕ1ϕ2⋯ϕNd] span the range space of K.The rigid body modes, and deformational modes satisfy the orthogonal condition of
(7)RTMR=I,ΦdTMΦd=I,RTMΦd=0.
In this research, the rigid body modes are proposed to be formulated by the geometric node locations of the structure, other than being extracted from a shift eigensolver or determined by the null space of the rank deficient stiffness matrix. For a two-dimensional structure having N nodes, the three independent rigid body modes are the x translation (Rx=1, Ry=0), the y translation (Rx=0, Ry=1), and the z rotation (Rx=-y, Ry=x); that is,
(8)RT=[1001⋯000100⋯10-y1x11-y2⋯xN1].
The columns of R can be orthogonalized and be normalized with respect to mass matrix [24].


### 2.2. The Formulation of a Well-Conditioned Eigenequation

As the eigenvalues associated with the rigid body modes are zeros, the stiffness matrix (3) of a free-free structure is rewritten by the deformational eigenmodes as
(9)K=∑i=1Nλi(Mϕi)(Mϕi)T=∑i=1Ndλi(Mϕi)(Mϕi)T.
Mathematically, the flexibility matrix of a structure is defined as the inverse of the stiffness matrix, and it is expressed as F=∑i=1N(1/λi)ϕiϕiT=ΦΛ-1ΦT (as (4)). Since the rigid body eigenvalues are zeros, the flexibility matrix is positive infinity for a free-free structure. Physically, the flexibility is defined as the displacement response of a structure when a unit force is applied to it. A unit force applied to a free-free structure will make the structure move freely. In this viewpoint, the flexibility matrix of a free-free structure does not exist and cannot be determined like the usually fixed structure in civil engineering [26]. Here, a modal flexibility matrix is employed which is contributed by the deformational modes solely as
(10)F=∑i=1Nd1λiϕiϕiT=ΦdΛd-1ΦdT.
The stiffness matrix and modal flexibility matrix are orthogonal to the rigid body modes
(11)KR=RTK=0,FMR=RTMF=0.
The stiffness matrix and modal flexibility matrix satisfy the relation of
(12)KF=(MΦdΛdΦdTM)(ΦdΛd-1ΦdT)=MΦdΦdT=I-MRRT.
Mathematically, the stiffness matrix and modal flexibility matrix, which are formed from the deformational eigenmodes, are singular and rank deficient for a free-free structure. The stiffness matrix and modal flexibility matrix are Moore-Penrose pseudoinverse of each other. The extraction of Moore-Penrose pseudoinverse is computationally expensive and not accurate. Herein, a new form of the stiffness and flexibility matrices is defined as
(13)K-=K+α(MR)(MR)T=M(∑i=1NdλiϕiϕiT+∑r=1NrαrϕrϕrT)M,F-=F+βRRT=∑i=1Nd1λiϕiϕiT+∑r=1NrβrϕrϕrT,
where K- and F- are well conditioned and full rank, and they are hereinafter called the generalized stiffness matrix and the generalized flexibility matrix, respectively. αr and βr are the participation factors of the rigid body modes to the generalized stiffness and flexibility matrices, with αr>0 and βr>0 (r=1,2,…,Nr).The generalized stiffness and flexibility matrices are related to the rigid body modes by(14)K-R=RTK-=αMR,F-MR=RTMF-=βR,RTK-R=αI,RTMF-MR=βI,(15)K-F-=(K+α(MR)(MR)T)(F+βRRT)=KF+αβMRRT.In particular, if the participation factors are chosen asαr=1, βr=1, (r=1,2,…,Nr), that is,(16a)K-=K+(MR)(MR)T=M(∑i=1NdλiϕiϕiT+∑r=1NrϕrϕrT)M(16b)F-=F+RRT=∑i=1Nd1λiϕiϕiT+∑r=1NrϕrϕrT, the generalized stiffness and flexibility matrices are related to the rigid body modes as
(17)K-R=RK-=MR,F-MR=RMF-=R(18)RK-R=I,RMF-MR=I(19)K-F-=(K+(MR)(MR)T)(F+RRT)=MΦdΦdT+MRRT=I.
Considering (19), the generalized stiffness matrix K- and generalized flexibility matrix F- are dual inverse of each other as
(20)(K+(MR)(MR)T)-1=F+RRT=∑i=1N1λiϕiϕiT,(F+RRT)-1=K+(MR)(MR)T=M(∑i=1NλiϕiϕiT)M.
The modal flexibility matrix is determined from the stiffness matrix by
(21)F=(K+(MR)(MR)T)-1-RRT,
and the modal flexibility can be transformed into the stiffness matrix in the form of
(22)K=(F+RRT)-1-(MR)(MR)T.
The generalized stiffness and flexibility matrices render a simple and effective transformation between the stiffness and modal flexibility matrices, avoiding the expensive pseudoinverse.


### 2.3. A Well-Conditioned Eigenproblem

The generalized stiffness matrixK- leads to the eigenequation
(23)(K--λM)Φ=0.
Equation (23) has the identical eigenvectors to those of (1), and only the eigenvalues of the rigid body modes are changed from 0 to 1. This well-conditioned eigenequation inherently avoids the zero eigenvalues, which leads to a faster and more accurate extraction of eigensolutions than the traditional strategy with a shift eigensolver [23]. Equation (23) extracts the real eigenvalues and eigenvectors of a free-free structure, whereas the eigenvalues of the rigid body modes are changed to 1.The well-conditioned eigenequation can also be employed for the calculation of the eigenvalue derivatives and eigenvector derivatives using Nelson’s method [27] or the modal method [28]. For the sake of simplicity, the detailed calculation of eigensensitivity will not be demonstrated here. The proposed method for eigensensitivity could be more efficient and accurate than the traditional method which introduces a small shift in the stiffness matrix to calculate the eigensensitivity of a free-free structure [23].


### 2.4. Calculation of Residual Flexibility and Its Derivative

In some forward substructuring methods [3–14], only a few lower modes of a substructure are calculated to assemble the global structure while the residual flexibility matrix is required for the compensation of the higher modes. In consequence, the residual flexibility and its derivative matrices are required for the calculation of the eigensolutions and eigensensitivity [9–12]. The detailed substructuring method will be described in the next section. Based on the proposed well-conditioned stiffness and flexibility matrices, the general formulation of the first-order and high-order residual flexibility matrices and their derivatives is derived in this section for a free-free structure.


#### 2.4.1. Residual Flexibility

The complete eigenmodes of a structure are divided into theNm master modes Φm, which are usually the lower modes of a structure, and the residual Ns slave modes Φs. For a free-free structure, the master modes Φm include the Nr rigid body modes R and the (Nm-Nr) deformational master modes Φm-r. The relation between the master modes, slave modes, rigid body modes and deformational modes is described in Figure 1.Figure 1
Division of master and slave modes in a substructure.In the substructuring methods, only some master modes are calculated, while the slave modes are discarded and compensated with a residual flexibility. According to (21), the first-order residual flexibility matrix can be expressed by the generalized stiffness matrix and master modes as
(24)G(1)=ΦsΛs-1ΦsT=(K+(MR)(MR)T)-1-Φm-rΛm-r-1Φm-rT-RRT=K--1-Φm-rΛm-r-1Φm-rT-RRT.
The second-order residual flexibility matrix can be expressed in terms of the generalized stiffness matrix and master modes as
(25)G(2)=ΦsΛs-2ΦsT=(ΦsΛs-1ΦsT)M(ΦsΛs-1ΦsT)=K--1MK--1-Φm-rΛm-r-2Φm-rT-RRT.
In general, the k-order residual flexibility is given by
(26)G(k)=ΦsΛs-kΦsT=(ΦsΛs-1ΦsT)(M(ΦsΛs-1ΦsT))k-1=K--1(MK--1)k-1-Φm-rΛm-r-kΦm-rT-RRT,
where the subscript k indicates the kth-order residual flexibility.Due to the orthogonal properties of the eigenmodes, thekth-order residual flexibility can also be generally expressed by the lower-order residual flexibility matrices as(27a)G(k)=G(1)MG(k-1)=G(1)(MG(1))k-1
or
(27b)G(k)=G(c)MG(k-c)=G(c)(MG(1))k-c(c=1,2,…,k).


#### 2.4.2. Derivative of the Residual Flexibility

Differentiating (24) with respect to an elemental parameter a, the derivative matrix of the first-order residual flexibility can be expressed by the generalized stiffness matrix and master modes as
(28)∂G(1)∂a=∂(ΦsΛs-1ΦsT)∂a=∂(K+(MR)(MR)T)-1∂a-∂(Φm-rΛm-r-1Φm-rT)∂a=-K--2∂K∂a-(2Φm-rΛm-r-1∂(Φm-rT)∂a00000000000000+Φm-r∂(Λm-r-1)∂aΦm-rT).
In general, the kth-order residual flexibility (26) has the derivative matrix
(29)∂G(k)∂a=∂(ΦsΛs-kΦsT)∂a=∂(K--1(MK--1)k-1)∂a-(2Φm-rΛm-r-k∂(Φm-rT)∂a+Φm-r∂(Λm-r-k)∂aΦm-rT).
According to (27a), the derivative of the kth-order residual flexibility can also be obtained by those of the lower order residual flexibility matrices as
(30)∂G(k)∂a=∂G(1)∂aMG(k-1)+G(1)M∂G(k-1)∂a,∂G(k)∂a=∂G(c)∂aMG(k-c)+G(c)M∂G(k-c)∂a∂G(k)∂a=∂G(c)∂aMG(k-c)+(c=1,2,…,k).


### 2.5. Illustrative Example: Spring-Mass Model

A four-DOF spring-mass model without constraint (Figure2) is employed to illustrate the formulation of the modal flexibility and residual flexibility from the proposed well-conditioned stiffness matrix.Figure 2
The spring-mass model under free-free condition.The stiffness parameters of the three springs are set tok1=k2=k3=20 N/m. The four masses are set to m1=1 kg, m2=2 kg, m3=2 kg, and m4=1 kg. The stiffness matrix of the model is
(31)K=[20-2000-2040-2000-2040-2000-2020].
The stiffness matrix is singular and rank deficient. To form the proposed full-rank stiffness matrix, the mass normalized rigid body mode is constructed according to (8) and is given by
(32)RT=[0.40820.40820.40820.4082].
In consequence, the generalized stiffness matrix is formed as
(33)K-=K+(MR)(MR)T=[20.1667-19.66670.3333-0.1667-19.666740.6667-19.33330.33330.3333-19.333340.6667-19.66670.16670.3333-19.666720.1667].
The eigenequation (K--λM)Φ=0 leads to the eigenvalues and eigenvectors
(34)Λ=[1000010000030000040],Φd=[0.4082-0.57740.57740.40820.4082-0.2887-0.2887-0.40820.40820.28870.28870.40820.40820.5774-0.5774-0.4082].
The eigenequation formulated by the generalized stiffness matrix K- is full rank and well conditioned. The eigenvalue of the rigid body mode is 1 as expected.If the first three modes are chosen as the master modes, the residual flexibility is formulated by the generalized stiffness matrix as(35)G=K--1-∑i=2,31λiϕiϕiT-RRT=[0.0042-0.00420.0042-0.0042-0.00420.0042-0.00420.00420.0042-0.00420.0042-0.0042-0.00420.0042-0.00420.0042].


## 2.1. Basic Theory for Eigenanalysis

A structure withN degrees of freedom (DOF) has the eigenequation of
(1)K{ϕi}=λiM{ϕi},
where K and M are the stiffness and mass matrices. λi is the ith eigenvalue of the structure, and {ϕi} is the corresponding eigenvector. They are determined by the physical property of a structure, such as Young’s modulus, density, Poisson ratio, and geometric dimension. The eigensolutions of (1) consist of the eigenvalues Λ=Diag(λ1λ2⋯λN) and the corresponding eigenvectors Φ=[ϕ1ϕ2⋯ϕN]. The eigenvectors are mass-normalized such that they satisfy the following relation:
(2)ΦTKΦ=Λ,ΦTMΦ=I.
The stiffness matrix can be written by the mass-normalized eigenmodes as
(3)K=∑i=1Nλi(Mϕi)(Mϕi)T=MΦΛΦTM.
A flexibility matrix has a very straightforward physical interpretation: the displacement response caused by an applied unit load [26]. The flexibility matrix can also be written by the mass-normalized eigenmodes as
(4)F=∑i=1N1λiϕiϕiT=ΦΛ-1ΦT.
For a fixed structure, the stiffness matrix and flexibility matrix normally formed a dual inverse of each other as
(5)KF=I,F=K-1,K=F-1.
The displacement {x} of a free-free structure can be written as a superposition of the deformational and rigid body motions,
(6){x}={xd}+{xr}=Φd{q}+R{α},
where {xd} is the displacement due to the deformational motion and {xr} is the displacement due to the rigid body motion. Φd is the linear orthogonal deformational modes, R is the orthogonal rigid body modes, and {q} and {α} are the participation factors of the orthogonal modes. The subscripts “d” and “r” hereinafter, respectively, represent the variables associated with the deformational motion and rigid body motion.Accordingly, eigenequation (1) for a free-free structure has two kinds of eigenpairs.(1)
Nr zero eigenvalues pertaining to the rigid body motions: the associated eigenvectors span the null space of the stiffness matrix K, which contribute to the columns of R. Nr is equal to the number of the statically determinate constraints required to prevent all rigid body motion [24].(2)
Nd=N-Nr nonzero eigenvalues λi(i=1,2,…,Nd): the associated orthogonal deformational eigenvectors Φd=[ϕ1ϕ2⋯ϕNd] span the range space of K.The rigid body modes, and deformational modes satisfy the orthogonal condition of
(7)RTMR=I,ΦdTMΦd=I,RTMΦd=0.
In this research, the rigid body modes are proposed to be formulated by the geometric node locations of the structure, other than being extracted from a shift eigensolver or determined by the null space of the rank deficient stiffness matrix. For a two-dimensional structure having N nodes, the three independent rigid body modes are the x translation (Rx=1, Ry=0), the y translation (Rx=0, Ry=1), and the z rotation (Rx=-y, Ry=x); that is,
(8)RT=[1001⋯000100⋯10-y1x11-y2⋯xN1].
The columns of R can be orthogonalized and be normalized with respect to mass matrix [24].


## 2.2. The Formulation of a Well-Conditioned Eigenequation

As the eigenvalues associated with the rigid body modes are zeros, the stiffness matrix (3) of a free-free structure is rewritten by the deformational eigenmodes as
(9)K=∑i=1Nλi(Mϕi)(Mϕi)T=∑i=1Ndλi(Mϕi)(Mϕi)T.
Mathematically, the flexibility matrix of a structure is defined as the inverse of the stiffness matrix, and it is expressed as F=∑i=1N(1/λi)ϕiϕiT=ΦΛ-1ΦT (as (4)). Since the rigid body eigenvalues are zeros, the flexibility matrix is positive infinity for a free-free structure. Physically, the flexibility is defined as the displacement response of a structure when a unit force is applied to it. A unit force applied to a free-free structure will make the structure move freely. In this viewpoint, the flexibility matrix of a free-free structure does not exist and cannot be determined like the usually fixed structure in civil engineering [26]. Here, a modal flexibility matrix is employed which is contributed by the deformational modes solely as
(10)F=∑i=1Nd1λiϕiϕiT=ΦdΛd-1ΦdT.
The stiffness matrix and modal flexibility matrix are orthogonal to the rigid body modes
(11)KR=RTK=0,FMR=RTMF=0.
The stiffness matrix and modal flexibility matrix satisfy the relation of
(12)KF=(MΦdΛdΦdTM)(ΦdΛd-1ΦdT)=MΦdΦdT=I-MRRT.
Mathematically, the stiffness matrix and modal flexibility matrix, which are formed from the deformational eigenmodes, are singular and rank deficient for a free-free structure. The stiffness matrix and modal flexibility matrix are Moore-Penrose pseudoinverse of each other. The extraction of Moore-Penrose pseudoinverse is computationally expensive and not accurate. Herein, a new form of the stiffness and flexibility matrices is defined as
(13)K-=K+α(MR)(MR)T=M(∑i=1NdλiϕiϕiT+∑r=1NrαrϕrϕrT)M,F-=F+βRRT=∑i=1Nd1λiϕiϕiT+∑r=1NrβrϕrϕrT,
where K- and F- are well conditioned and full rank, and they are hereinafter called the generalized stiffness matrix and the generalized flexibility matrix, respectively. αr and βr are the participation factors of the rigid body modes to the generalized stiffness and flexibility matrices, with αr>0 and βr>0 (r=1,2,…,Nr).The generalized stiffness and flexibility matrices are related to the rigid body modes by(14)K-R=RTK-=αMR,F-MR=RTMF-=βR,RTK-R=αI,RTMF-MR=βI,(15)K-F-=(K+α(MR)(MR)T)(F+βRRT)=KF+αβMRRT.In particular, if the participation factors are chosen asαr=1, βr=1, (r=1,2,…,Nr), that is,(16a)K-=K+(MR)(MR)T=M(∑i=1NdλiϕiϕiT+∑r=1NrϕrϕrT)M(16b)F-=F+RRT=∑i=1Nd1λiϕiϕiT+∑r=1NrϕrϕrT, the generalized stiffness and flexibility matrices are related to the rigid body modes as
(17)K-R=RK-=MR,F-MR=RMF-=R(18)RK-R=I,RMF-MR=I(19)K-F-=(K+(MR)(MR)T)(F+RRT)=MΦdΦdT+MRRT=I.
Considering (19), the generalized stiffness matrix K- and generalized flexibility matrix F- are dual inverse of each other as
(20)(K+(MR)(MR)T)-1=F+RRT=∑i=1N1λiϕiϕiT,(F+RRT)-1=K+(MR)(MR)T=M(∑i=1NλiϕiϕiT)M.
The modal flexibility matrix is determined from the stiffness matrix by
(21)F=(K+(MR)(MR)T)-1-RRT,
and the modal flexibility can be transformed into the stiffness matrix in the form of
(22)K=(F+RRT)-1-(MR)(MR)T.
The generalized stiffness and flexibility matrices render a simple and effective transformation between the stiffness and modal flexibility matrices, avoiding the expensive pseudoinverse.


## 2.3. A Well-Conditioned Eigenproblem

The generalized stiffness matrixK- leads to the eigenequation
(23)(K--λM)Φ=0.
Equation (23) has the identical eigenvectors to those of (1), and only the eigenvalues of the rigid body modes are changed from 0 to 1. This well-conditioned eigenequation inherently avoids the zero eigenvalues, which leads to a faster and more accurate extraction of eigensolutions than the traditional strategy with a shift eigensolver [23]. Equation (23) extracts the real eigenvalues and eigenvectors of a free-free structure, whereas the eigenvalues of the rigid body modes are changed to 1.The well-conditioned eigenequation can also be employed for the calculation of the eigenvalue derivatives and eigenvector derivatives using Nelson’s method [27] or the modal method [28]. For the sake of simplicity, the detailed calculation of eigensensitivity will not be demonstrated here. The proposed method for eigensensitivity could be more efficient and accurate than the traditional method which introduces a small shift in the stiffness matrix to calculate the eigensensitivity of a free-free structure [23].


## 2.4. Calculation of Residual Flexibility and Its Derivative

In some forward substructuring methods [3–14], only a few lower modes of a substructure are calculated to assemble the global structure while the residual flexibility matrix is required for the compensation of the higher modes. In consequence, the residual flexibility and its derivative matrices are required for the calculation of the eigensolutions and eigensensitivity [9–12]. The detailed substructuring method will be described in the next section. Based on the proposed well-conditioned stiffness and flexibility matrices, the general formulation of the first-order and high-order residual flexibility matrices and their derivatives is derived in this section for a free-free structure.


### 2.4.1. Residual Flexibility

The complete eigenmodes of a structure are divided into theNm master modes Φm, which are usually the lower modes of a structure, and the residual Ns slave modes Φs. For a free-free structure, the master modes Φm include the Nr rigid body modes R and the (Nm-Nr) deformational master modes Φm-r. The relation between the master modes, slave modes, rigid body modes and deformational modes is described in Figure 1.Figure 1
Division of master and slave modes in a substructure.In the substructuring methods, only some master modes are calculated, while the slave modes are discarded and compensated with a residual flexibility. According to (21), the first-order residual flexibility matrix can be expressed by the generalized stiffness matrix and master modes as
(24)G(1)=ΦsΛs-1ΦsT=(K+(MR)(MR)T)-1-Φm-rΛm-r-1Φm-rT-RRT=K--1-Φm-rΛm-r-1Φm-rT-RRT.
The second-order residual flexibility matrix can be expressed in terms of the generalized stiffness matrix and master modes as
(25)G(2)=ΦsΛs-2ΦsT=(ΦsΛs-1ΦsT)M(ΦsΛs-1ΦsT)=K--1MK--1-Φm-rΛm-r-2Φm-rT-RRT.
In general, the k-order residual flexibility is given by
(26)G(k)=ΦsΛs-kΦsT=(ΦsΛs-1ΦsT)(M(ΦsΛs-1ΦsT))k-1=K--1(MK--1)k-1-Φm-rΛm-r-kΦm-rT-RRT,
where the subscript k indicates the kth-order residual flexibility.Due to the orthogonal properties of the eigenmodes, thekth-order residual flexibility can also be generally expressed by the lower-order residual flexibility matrices as(27a)G(k)=G(1)MG(k-1)=G(1)(MG(1))k-1
or
(27b)G(k)=G(c)MG(k-c)=G(c)(MG(1))k-c(c=1,2,…,k).


### 2.4.2. Derivative of the Residual Flexibility

Differentiating (24) with respect to an elemental parameter a, the derivative matrix of the first-order residual flexibility can be expressed by the generalized stiffness matrix and master modes as
(28)∂G(1)∂a=∂(ΦsΛs-1ΦsT)∂a=∂(K+(MR)(MR)T)-1∂a-∂(Φm-rΛm-r-1Φm-rT)∂a=-K--2∂K∂a-(2Φm-rΛm-r-1∂(Φm-rT)∂a00000000000000+Φm-r∂(Λm-r-1)∂aΦm-rT).
In general, the kth-order residual flexibility (26) has the derivative matrix
(29)∂G(k)∂a=∂(ΦsΛs-kΦsT)∂a=∂(K--1(MK--1)k-1)∂a-(2Φm-rΛm-r-k∂(Φm-rT)∂a+Φm-r∂(Λm-r-k)∂aΦm-rT).
According to (27a), the derivative of the kth-order residual flexibility can also be obtained by those of the lower order residual flexibility matrices as
(30)∂G(k)∂a=∂G(1)∂aMG(k-1)+G(1)M∂G(k-1)∂a,∂G(k)∂a=∂G(c)∂aMG(k-c)+G(c)M∂G(k-c)∂a∂G(k)∂a=∂G(c)∂aMG(k-c)+(c=1,2,…,k).


## 2.4.1. Residual Flexibility

The complete eigenmodes of a structure are divided into theNm master modes Φm, which are usually the lower modes of a structure, and the residual Ns slave modes Φs. For a free-free structure, the master modes Φm include the Nr rigid body modes R and the (Nm-Nr) deformational master modes Φm-r. The relation between the master modes, slave modes, rigid body modes and deformational modes is described in Figure 1.Figure 1
Division of master and slave modes in a substructure.In the substructuring methods, only some master modes are calculated, while the slave modes are discarded and compensated with a residual flexibility. According to (21), the first-order residual flexibility matrix can be expressed by the generalized stiffness matrix and master modes as
(24)G(1)=ΦsΛs-1ΦsT=(K+(MR)(MR)T)-1-Φm-rΛm-r-1Φm-rT-RRT=K--1-Φm-rΛm-r-1Φm-rT-RRT.
The second-order residual flexibility matrix can be expressed in terms of the generalized stiffness matrix and master modes as
(25)G(2)=ΦsΛs-2ΦsT=(ΦsΛs-1ΦsT)M(ΦsΛs-1ΦsT)=K--1MK--1-Φm-rΛm-r-2Φm-rT-RRT.
In general, the k-order residual flexibility is given by
(26)G(k)=ΦsΛs-kΦsT=(ΦsΛs-1ΦsT)(M(ΦsΛs-1ΦsT))k-1=K--1(MK--1)k-1-Φm-rΛm-r-kΦm-rT-RRT,
where the subscript k indicates the kth-order residual flexibility.Due to the orthogonal properties of the eigenmodes, thekth-order residual flexibility can also be generally expressed by the lower-order residual flexibility matrices as(27a)G(k)=G(1)MG(k-1)=G(1)(MG(1))k-1
or
(27b)G(k)=G(c)MG(k-c)=G(c)(MG(1))k-c(c=1,2,…,k).


## 2.4.2. Derivative of the Residual Flexibility

Differentiating (24) with respect to an elemental parameter a, the derivative matrix of the first-order residual flexibility can be expressed by the generalized stiffness matrix and master modes as
(28)∂G(1)∂a=∂(ΦsΛs-1ΦsT)∂a=∂(K+(MR)(MR)T)-1∂a-∂(Φm-rΛm-r-1Φm-rT)∂a=-K--2∂K∂a-(2Φm-rΛm-r-1∂(Φm-rT)∂a00000000000000+Φm-r∂(Λm-r-1)∂aΦm-rT).
In general, the kth-order residual flexibility (26) has the derivative matrix
(29)∂G(k)∂a=∂(ΦsΛs-kΦsT)∂a=∂(K--1(MK--1)k-1)∂a-(2Φm-rΛm-r-k∂(Φm-rT)∂a+Φm-r∂(Λm-r-k)∂aΦm-rT).
According to (27a), the derivative of the kth-order residual flexibility can also be obtained by those of the lower order residual flexibility matrices as
(30)∂G(k)∂a=∂G(1)∂aMG(k-1)+G(1)M∂G(k-1)∂a,∂G(k)∂a=∂G(c)∂aMG(k-c)+G(c)M∂G(k-c)∂a∂G(k)∂a=∂G(c)∂aMG(k-c)+(c=1,2,…,k).


## 2.5. Illustrative Example: Spring-Mass Model

A four-DOF spring-mass model without constraint (Figure2) is employed to illustrate the formulation of the modal flexibility and residual flexibility from the proposed well-conditioned stiffness matrix.Figure 2
The spring-mass model under free-free condition.The stiffness parameters of the three springs are set tok1=k2=k3=20 N/m. The four masses are set to m1=1 kg, m2=2 kg, m3=2 kg, and m4=1 kg. The stiffness matrix of the model is
(31)K=[20-2000-2040-2000-2040-2000-2020].
The stiffness matrix is singular and rank deficient. To form the proposed full-rank stiffness matrix, the mass normalized rigid body mode is constructed according to (8) and is given by
(32)RT=[0.40820.40820.40820.4082].
In consequence, the generalized stiffness matrix is formed as
(33)K-=K+(MR)(MR)T=[20.1667-19.66670.3333-0.1667-19.666740.6667-19.33330.33330.3333-19.333340.6667-19.66670.16670.3333-19.666720.1667].
The eigenequation (K--λM)Φ=0 leads to the eigenvalues and eigenvectors
(34)Λ=[1000010000030000040],Φd=[0.4082-0.57740.57740.40820.4082-0.2887-0.2887-0.40820.40820.28870.28870.40820.40820.5774-0.5774-0.4082].
The eigenequation formulated by the generalized stiffness matrix K- is full rank and well conditioned. The eigenvalue of the rigid body mode is 1 as expected.If the first three modes are chosen as the master modes, the residual flexibility is formulated by the generalized stiffness matrix as(35)G=K--1-∑i=2,31λiϕiϕiT-RRT=[0.0042-0.00420.0042-0.0042-0.00420.0042-0.00420.00420.0042-0.00420.0042-0.0042-0.00420.0042-0.00420.0042].


## 3. Substructure-Based Model Updating

In the sensitivity-based model updating procedure, the general objective function combining the modal properties of the frequencies and the mode shapes is usually denoted by [29, 30]
(36)J(r)=∑iWλi2[λi({a})M-λiE]2+∑iWϕi2∑j[ϕji({a})M-ϕjiE]2,
where λiE represents the ith experimental frequencies and ϕjiE is the ith experimental mode shape at the jth point. λiM and ϕjiM denote the corresponding frequencies and mode shapes from the analytical FE model, which are expressed as the function of the uncertain physical parameters {a}. Wλi2 and Wϕi2 are the weight coefficients due to the different measurement accuracy of the frequencies and mode shapes.The eigensolutions are used to form the objective function. The objective function, formed from the residuals between the eigensolutions of the FE model and the modal properties of the practical structures, is minimized by adjusting continuously the parameters{a} of the initial FE model through the optimization searching techniques. The eigensensitivity calculates the changes in the eigensolutions caused by the perturbations of the design parameters of a structural model. It serves for indicating the searching direction of an optimization algorithm, which endows the more sensitive parameter (with respect to the objective function) a higher priority.


### 3.1. Substructure Method to Eigensolutions

The global structure withN DOF is divided into NS substructures. Treating the jth substructure of n(j) DOF (j=1,2,…,NS) as an independent structure, it has the stiffness matrix K(j) and mass matrix M(j). If the jth substructure is a free-free structure after division, the stiffness matrix K(j) is singular and rank deficient. Then the generalized stiffness matrix K-(j) is used for the free-free substructure, replacing the stiffness matrix K(j). The eigenequation for the jth substructure is written as
(37)(K(j)-λ(j)M(j))Φ(j)=0.
It is noted that K-(j) is used in (37) for a free-free substructure, replacing the stiffness matrix K(j). Hereinafter, only K(j) is used in the formula for clearance, and it means the proposed generalized stiffness matrix (K-(j)) for a free-free substructure. The substructural eigenequation has n(j) pairs of eigenvalues and eigenvectors, which are orthogonal to the stiffness and mass matrices as [5]
(38)Λ(j)=Diag[λ1(j),λ2(j),…,λn(j)(j)],Φ(j)=[ϕ1(j),ϕ2(j),…,ϕn(j)(j)],[Φ(j)]TK(j)Φ(j)=Λ(j),[Φ(j)]TM(j)Φ(j)=I(j),(j=1,2,…,NS).Based on the principle of virtual work and geometric compatibility, the substructuring method [5, 8] reconstructs the eigensolutions of the global structure by imposing the constraints at the interfaces as(39a)[Λp-λ-I-Γ-ΓT0]{zτ}={00}.
In this equation,
(39b)Γ=[CΦp]T,Λp=Diag[Λ(1),Λ(2),…,Λ(NS)],Φp=Diag[Φ(1),Φ(2),…,Φ(NS)],where C is a rectangular connection matrix constraining the interface DOF of the adjacent substructures to move jointly [5]. τ is the internal connection forces of the adjacent substructures. λ- is the eigenvalue of the global structure. z acts as the participation factor of the substructural eigenmodes, and the eigenvectors of the global structure can be recovered by Φ-=Φp{z}. Superscript “p” denotes the diagonal assembly of the independent substructural matrices before constraining the independent substructures at the interfaces.From the viewpoint of energy conservation, all modes of the substructures contribute to the eigenmodes of the global structure; that is, the complete eigensolutions of all substructures are required to assemble the primitive form ofΛp and Φp. It is inefficient and not worthwhile to calculate all modes of the substructures, as only a few eigenmodes are generally of interest for a large-scale structure. To overcome this difficulty, only the master modes, corresponding to the lower vibration modes, are calculated to assemble the eigenequation of the global structure, while the slave modes (residual higher modes) are discarded and compensated by the residual flexibility in the later calculations. From here on, subscripts “m” and “s” will denote the “master” and “slave” modes, respectively.Eigenequation (39a) and (39b) is rewritten according to the master modes and slave modes as
(40)[Λmp-λ-I0-Γm0Λsp-λ-I-Γs-ΓmT-ΓsT0]{zmzsτ}={000},
where Γm=[CΦmp]T,Γs=[CΦsp]T, Λmp and Φmp include the master eigenvalues and eigenvectors of the independent substructures, Λsp and Φsp include the slave eigenvalues and eigenvectors of the independent substructures, and zm and zs are the mode participation factors of the master modes and slave modes.With the second line of (40), the slave part of the mode participation factor can be expressed as
(41)zs=(Λsp-λ-I)-1Γsτ.
Substituting (41) into (40) gives
(42)[Λmp-λ-I-Γm-ΓmT-ΓsT(Λsp-λ-I)-1Γs]{zmτ}={00}.
In (42), Taylor expansion of the nonlinear item (Λsp-λ-I)-1 has
(43)(Λsp-λ-I)-1=(Λsp)-1+λ-(Λsp)-2+λ-2(Λsp)-3+⋯.In general, the required eigenvaluesλ- correspond to the lowest modes of the global structure and are far less than the values in Λsp when the master modes are appropriately chosen. In that case, retaining only the first item of the Taylor expansion gives ΓsT(Λsp-λ-I)-1Γs≈ΓsT(Λsp)-1Γs. In consequence, (42) is reduced into [5](44a)[(Λmp-λ-Im)+Γmς-1ΓmT]zm=0,(44b)ς=ΓsT(Λsp)-1Γs.The size of the reduced eigenequation (44a) and (44b) is equal to the number of the retained master modes NPm×NPm, which is much smaller than the original one of NP×NP in (39a) and (39b). λ- and zm can be solved from this reduced eigenequation using the common eigensolvers. As before, the eigenvalues of the global structure are λ-, and the eigenvectors of the global structure are recovered by Φ-=Φmpzm. ς=ΓsT(Λsp)-1Γs is associated with the first-order residual flexibility that can be calculated using the master modes of the substructures as(45)ΓsT(Λsp)-1Γs=CΦsp(Λsp)-1[Φsp]TCT,Φsp(Λsp)-1[Φsp]T=[(K(1))-1-Φm(1)(Λm(1))-1[Φm(1)]T000⋱000(K(NS))-1-Φm(NS)(Λm(NS))-1[Φm(NS)]T].


### 3.2. Eigensensitivity with Substructuring Method

The eigensensitivity of theith mode (i=1,2,…,N) with respect to an elemental parameter will be derived in this section. The elemental parameter is chosen to be the stiffness parameter, such as the bending rigidity of an element, and denoted by parameter a in the Ath substructure. The reduced eigenequation (44a) and (44b) is rewritten for the ith mode as
(46)[(Λmp-λ-iIm)+Γmς-1ΓmT]{zi}={0}.
Equation (46) is differentiated with respect to parameter a as
(47)[(Λmp-λ-iIm)+Γmς-1ΓmT]∂{zi}∂a+∂[(Λmp-λ-iIm)+Γmς-1ΓmT]∂a{zi}={0}.
Since [(Λmp-λ-iIm)+Γmς-1ΓmT] is symmetric, premultiplying {zi}T on both sides of (47) gives the eigenvalue derivative of the ith mode as
(48)∂λ-i∂a={zi}T[∂Λmp∂a+∂(Γmς-1ΓmT)∂a]{zi},
where
(49)∂(Γmς-1ΓmT)∂a=∂Γm∂aς-1ΓmT-Γmς-1∂ς∂aς-1ΓmT+Γmς-1∂ΓmT∂a.
In this equation, the derivative matrices ∂Λmp/∂a, ∂Γm/∂a, and ∂ς/∂a are formed using the eigenvalue derivatives, eigenvector derivatives, and residual flexibility derivatives of the substructures.Since the substructures are independent, these derivative matrices are calculated within the Ath substructure solely, while those in other substructures are zeros; that is,
(50)∂Λmp∂a=[0000∂Λm(A)∂a0000],∂ΓmT∂a=C∂Φmp∂a=C[0000∂Φm(A)∂a0000]∂ς∂r=C×Diag[0000∂((K(A))-1-Φm(A)(Λm(A))-1[Φm(A)]T)∂a0000]×CT.{zi}, Γm, and ς-1 have been computed in the previous section for eigensolutions and can be reused here directly. ∂Λm(A)/∂a and ∂Φm(A)/∂a are the eigensolution derivatives of the master modes in the Ath substructures. They can be calculated with common methods, such as Nelson’s method [27], by treating the Ath substructure as one independent structure. Subsequently, the eigenvalue derivative of the global structure can be obtained from (48), and it solely relies on a particular substructure (the Ath substructure).The eigenvectors of the global structure are recovered byΦ-=Φmpzm. Hence, the ith eigenvector of the global structure can be expressed as
(51)Φ-i=Φmp{zi}.
Differentiating (51) with respect to the elemental parameter a, one can obtain the eigenvector derivative of the ith mode as
(52)∂Φ-i∂a=∂Φmp∂a{zi}+Φmp{∂zi∂a}.
In (52), Φmp and {zi} have been obtained when calculating the eigensolutions. ∂Φmp/∂a is associated with the eigenvector derivatives of the Ath substructure as (50). {∂zi/∂a} can be obtained from the reduced eigenequation (47), as described in the following.{
∂
z
i
/
∂
a
} is separated into the sum of a particular part and a homogeneous part as
(53){∂zi∂a}={vi}+ci{zi},
where ci is a participation factor. Substituting (53) into (47) gives(54a)Ψ{vi}={Yi},
where
(54b)Ψ=[(Λmp-λ-iIm)+Γmς-1ΓmT],{Yi}=-∂[(Λmp-λ-iIm)+Γmς-1ΓmT]∂a{zi}.Since the items inΨ and {Yi} are available when calculating the eigenvalue derivatives, the vector {vi} can be solved from (54a) and (54b) effortlessly.The eigenvector{zi} of the reduced eigenequation (44a) and (44b) satisfies the orthogonal condition of
(55){zi}T{zi}=1.
Differentiating (55) with respect to a gives
(56)∂{zi}T∂a{zi}+{zi}T∂{zi}∂a=0.
Substituting (53) into (56), the participation factor ci is therefore obtained as
(57)ci=-12({vi}T{zi}+{zi}T{vi}).
After the vector {vi} and the factor ci have been achieved, the eigenvector derivative of the global structure can be calculated from (52). Since the reduced eigenequation (44a) and (44b) is smaller in size compared to that of the global structure, calculation of {∂zi/∂a} can be processed much faster than that in the conventional Nelson’s method [27]. As calculation of the eigenvector derivatives dominates the whole model updating process, the substructuring method will improve the computational efficiency significantly [9, 11, 12].With the proposed substructuring method, the eigenvalue and eigenvector derivatives with respect to an elemental parameter are computed solely within the substructure that contains the element, whereas the derivative matrices of all other substructures with respect to the parameter are zero, thus allowing a significant reduction in computational cost.Based on the proposed full-rank well-conditioned substructural eigenequation, the substructure-based model updating is proceeded by the following procedure.(1)
Divide the global structure into several manageable substructures.(2)
Calculate the rigid body modes (R) for the free-free substructures according to (8).(3)
Construct the generalized stiffness matrix for the free-free substructures byK-=K+(MR)(MR)T.(4)
Construct the full-rank well-conditioned substructural eigenequation for the free-free substructures as (23). Based on the full-rank well-conditioned substructural eigenequation, the substructural eigensolutions and eigensensitivity of the master mode are calculated for the free-free substructures(5)
Calculate the generalized flexibility for the free-free substructures byF-=F+RRT. Based on the generalized flexibility matrix, the residual flexibility and its derivatives are calculated for the free-free substructures.(6)
Based on the substructural eigensolutions, eigensensitivity, and residual flexibility, the eigensolutions of the global structure is calculated by (44a) and (44b), and eigensensitivity of the global structure are calculated by (48) and (52). The eigensolutions of the global structure are used to construct the objective function in the model updating process, while the eigensensitivity is used for indicating the searching direction of the optimization process.The accuracy and efficiency of the full-rank well-conditioned substructural eigenequation in substructure-based model updating will be investigated by two examples in the following section.


## 3.1. Substructure Method to Eigensolutions

The global structure withN DOF is divided into NS substructures. Treating the jth substructure of n(j) DOF (j=1,2,…,NS) as an independent structure, it has the stiffness matrix K(j) and mass matrix M(j). If the jth substructure is a free-free structure after division, the stiffness matrix K(j) is singular and rank deficient. Then the generalized stiffness matrix K-(j) is used for the free-free substructure, replacing the stiffness matrix K(j). The eigenequation for the jth substructure is written as
(37)(K(j)-λ(j)M(j))Φ(j)=0.
It is noted that K-(j) is used in (37) for a free-free substructure, replacing the stiffness matrix K(j). Hereinafter, only K(j) is used in the formula for clearance, and it means the proposed generalized stiffness matrix (K-(j)) for a free-free substructure. The substructural eigenequation has n(j) pairs of eigenvalues and eigenvectors, which are orthogonal to the stiffness and mass matrices as [5]
(38)Λ(j)=Diag[λ1(j),λ2(j),…,λn(j)(j)],Φ(j)=[ϕ1(j),ϕ2(j),…,ϕn(j)(j)],[Φ(j)]TK(j)Φ(j)=Λ(j),[Φ(j)]TM(j)Φ(j)=I(j),(j=1,2,…,NS).Based on the principle of virtual work and geometric compatibility, the substructuring method [5, 8] reconstructs the eigensolutions of the global structure by imposing the constraints at the interfaces as(39a)[Λp-λ-I-Γ-ΓT0]{zτ}={00}.
In this equation,
(39b)Γ=[CΦp]T,Λp=Diag[Λ(1),Λ(2),…,Λ(NS)],Φp=Diag[Φ(1),Φ(2),…,Φ(NS)],where C is a rectangular connection matrix constraining the interface DOF of the adjacent substructures to move jointly [5]. τ is the internal connection forces of the adjacent substructures. λ- is the eigenvalue of the global structure. z acts as the participation factor of the substructural eigenmodes, and the eigenvectors of the global structure can be recovered by Φ-=Φp{z}. Superscript “p” denotes the diagonal assembly of the independent substructural matrices before constraining the independent substructures at the interfaces.From the viewpoint of energy conservation, all modes of the substructures contribute to the eigenmodes of the global structure; that is, the complete eigensolutions of all substructures are required to assemble the primitive form ofΛp and Φp. It is inefficient and not worthwhile to calculate all modes of the substructures, as only a few eigenmodes are generally of interest for a large-scale structure. To overcome this difficulty, only the master modes, corresponding to the lower vibration modes, are calculated to assemble the eigenequation of the global structure, while the slave modes (residual higher modes) are discarded and compensated by the residual flexibility in the later calculations. From here on, subscripts “m” and “s” will denote the “master” and “slave” modes, respectively.Eigenequation (39a) and (39b) is rewritten according to the master modes and slave modes as
(40)[Λmp-λ-I0-Γm0Λsp-λ-I-Γs-ΓmT-ΓsT0]{zmzsτ}={000},
where Γm=[CΦmp]T,Γs=[CΦsp]T, Λmp and Φmp include the master eigenvalues and eigenvectors of the independent substructures, Λsp and Φsp include the slave eigenvalues and eigenvectors of the independent substructures, and zm and zs are the mode participation factors of the master modes and slave modes.With the second line of (40), the slave part of the mode participation factor can be expressed as
(41)zs=(Λsp-λ-I)-1Γsτ.
Substituting (41) into (40) gives
(42)[Λmp-λ-I-Γm-ΓmT-ΓsT(Λsp-λ-I)-1Γs]{zmτ}={00}.
In (42), Taylor expansion of the nonlinear item (Λsp-λ-I)-1 has
(43)(Λsp-λ-I)-1=(Λsp)-1+λ-(Λsp)-2+λ-2(Λsp)-3+⋯.In general, the required eigenvaluesλ- correspond to the lowest modes of the global structure and are far less than the values in Λsp when the master modes are appropriately chosen. In that case, retaining only the first item of the Taylor expansion gives ΓsT(Λsp-λ-I)-1Γs≈ΓsT(Λsp)-1Γs. In consequence, (42) is reduced into [5](44a)[(Λmp-λ-Im)+Γmς-1ΓmT]zm=0,(44b)ς=ΓsT(Λsp)-1Γs.The size of the reduced eigenequation (44a) and (44b) is equal to the number of the retained master modes NPm×NPm, which is much smaller than the original one of NP×NP in (39a) and (39b). λ- and zm can be solved from this reduced eigenequation using the common eigensolvers. As before, the eigenvalues of the global structure are λ-, and the eigenvectors of the global structure are recovered by Φ-=Φmpzm. ς=ΓsT(Λsp)-1Γs is associated with the first-order residual flexibility that can be calculated using the master modes of the substructures as(45)ΓsT(Λsp)-1Γs=CΦsp(Λsp)-1[Φsp]TCT,Φsp(Λsp)-1[Φsp]T=[(K(1))-1-Φm(1)(Λm(1))-1[Φm(1)]T000⋱000(K(NS))-1-Φm(NS)(Λm(NS))-1[Φm(NS)]T].


## 3.2. Eigensensitivity with Substructuring Method

The eigensensitivity of theith mode (i=1,2,…,N) with respect to an elemental parameter will be derived in this section. The elemental parameter is chosen to be the stiffness parameter, such as the bending rigidity of an element, and denoted by parameter a in the Ath substructure. The reduced eigenequation (44a) and (44b) is rewritten for the ith mode as
(46)[(Λmp-λ-iIm)+Γmς-1ΓmT]{zi}={0}.
Equation (46) is differentiated with respect to parameter a as
(47)[(Λmp-λ-iIm)+Γmς-1ΓmT]∂{zi}∂a+∂[(Λmp-λ-iIm)+Γmς-1ΓmT]∂a{zi}={0}.
Since [(Λmp-λ-iIm)+Γmς-1ΓmT] is symmetric, premultiplying {zi}T on both sides of (47) gives the eigenvalue derivative of the ith mode as
(48)∂λ-i∂a={zi}T[∂Λmp∂a+∂(Γmς-1ΓmT)∂a]{zi},
where
(49)∂(Γmς-1ΓmT)∂a=∂Γm∂aς-1ΓmT-Γmς-1∂ς∂aς-1ΓmT+Γmς-1∂ΓmT∂a.
In this equation, the derivative matrices ∂Λmp/∂a, ∂Γm/∂a, and ∂ς/∂a are formed using the eigenvalue derivatives, eigenvector derivatives, and residual flexibility derivatives of the substructures.Since the substructures are independent, these derivative matrices are calculated within the Ath substructure solely, while those in other substructures are zeros; that is,
(50)∂Λmp∂a=[0000∂Λm(A)∂a0000],∂ΓmT∂a=C∂Φmp∂a=C[0000∂Φm(A)∂a0000]∂ς∂r=C×Diag[0000∂((K(A))-1-Φm(A)(Λm(A))-1[Φm(A)]T)∂a0000]×CT.{zi}, Γm, and ς-1 have been computed in the previous section for eigensolutions and can be reused here directly. ∂Λm(A)/∂a and ∂Φm(A)/∂a are the eigensolution derivatives of the master modes in the Ath substructures. They can be calculated with common methods, such as Nelson’s method [27], by treating the Ath substructure as one independent structure. Subsequently, the eigenvalue derivative of the global structure can be obtained from (48), and it solely relies on a particular substructure (the Ath substructure).The eigenvectors of the global structure are recovered byΦ-=Φmpzm. Hence, the ith eigenvector of the global structure can be expressed as
(51)Φ-i=Φmp{zi}.
Differentiating (51) with respect to the elemental parameter a, one can obtain the eigenvector derivative of the ith mode as
(52)∂Φ-i∂a=∂Φmp∂a{zi}+Φmp{∂zi∂a}.
In (52), Φmp and {zi} have been obtained when calculating the eigensolutions. ∂Φmp/∂a is associated with the eigenvector derivatives of the Ath substructure as (50). {∂zi/∂a} can be obtained from the reduced eigenequation (47), as described in the following.{
∂
z
i
/
∂
a
} is separated into the sum of a particular part and a homogeneous part as
(53){∂zi∂a}={vi}+ci{zi},
where ci is a participation factor. Substituting (53) into (47) gives(54a)Ψ{vi}={Yi},
where
(54b)Ψ=[(Λmp-λ-iIm)+Γmς-1ΓmT],{Yi}=-∂[(Λmp-λ-iIm)+Γmς-1ΓmT]∂a{zi}.Since the items inΨ and {Yi} are available when calculating the eigenvalue derivatives, the vector {vi} can be solved from (54a) and (54b) effortlessly.The eigenvector{zi} of the reduced eigenequation (44a) and (44b) satisfies the orthogonal condition of
(55){zi}T{zi}=1.
Differentiating (55) with respect to a gives
(56)∂{zi}T∂a{zi}+{zi}T∂{zi}∂a=0.
Substituting (53) into (56), the participation factor ci is therefore obtained as
(57)ci=-12({vi}T{zi}+{zi}T{vi}).
After the vector {vi} and the factor ci have been achieved, the eigenvector derivative of the global structure can be calculated from (52). Since the reduced eigenequation (44a) and (44b) is smaller in size compared to that of the global structure, calculation of {∂zi/∂a} can be processed much faster than that in the conventional Nelson’s method [27]. As calculation of the eigenvector derivatives dominates the whole model updating process, the substructuring method will improve the computational efficiency significantly [9, 11, 12].With the proposed substructuring method, the eigenvalue and eigenvector derivatives with respect to an elemental parameter are computed solely within the substructure that contains the element, whereas the derivative matrices of all other substructures with respect to the parameter are zero, thus allowing a significant reduction in computational cost.Based on the proposed full-rank well-conditioned substructural eigenequation, the substructure-based model updating is proceeded by the following procedure.(1)
Divide the global structure into several manageable substructures.(2)
Calculate the rigid body modes (R) for the free-free substructures according to (8).(3)
Construct the generalized stiffness matrix for the free-free substructures byK-=K+(MR)(MR)T.(4)
Construct the full-rank well-conditioned substructural eigenequation for the free-free substructures as (23). Based on the full-rank well-conditioned substructural eigenequation, the substructural eigensolutions and eigensensitivity of the master mode are calculated for the free-free substructures(5)
Calculate the generalized flexibility for the free-free substructures byF-=F+RRT. Based on the generalized flexibility matrix, the residual flexibility and its derivatives are calculated for the free-free substructures.(6)
Based on the substructural eigensolutions, eigensensitivity, and residual flexibility, the eigensolutions of the global structure is calculated by (44a) and (44b), and eigensensitivity of the global structure are calculated by (48) and (52). The eigensolutions of the global structure are used to construct the objective function in the model updating process, while the eigensensitivity is used for indicating the searching direction of the optimization process.The accuracy and efficiency of the full-rank well-conditioned substructural eigenequation in substructure-based model updating will be investigated by two examples in the following section.


## 4. Case Studies

### 4.1. Three-Span Frame Structure

The accuracy of the proposed well-conditioned eigenequation for calculation of substructural residual flexibility, eigensolutions, and eigensensitivity will be illustrated by a frame structure. The global frame is shown in Figure3. The material constants are chosen as bending rigidity (EI)=170×106 Nm2, axial rigidity (EA)=2500×106 N, mass per unit length (ρA)=110kg/m, and Poisson’s ratio = 0.3. The frame is discretized into 160 two-dimensional beam elements each 2.5 m long, which results in 140 nodes and 408 DOF (N=408). The frame is disassembled into three substructures (NS=3) when it is torn at 8 nodes as shown in Figure 3. After division, there are 51, 55, and 42 nodes in the three substructures with the DOF of n1=153, n2=165, and n3=114, respectively.The frame structure. →: the points and directions of experimental measurement.
(a)
The global structure
(b)
The divided substructuresIn this example, the first substructure is analyzed to investigate the accuracy of the proposed well-conditioned eigensolver in calculation of substructural eigensolutions and eigensensitivity for a free-free substructure. The first substructure has 153 DOF. The first 30 modes are calculated as the master modes to assemble the eigenequation of the global structure, and the other slave modes are compensated by the residual flexibility. As the first substructure is free after partition, the system matrices of the first substructure are singular and rankdeficient. Traditionally, the Moore-Penrose pseudoinverse is usually used for the analysis of rank-deficient matrix, or a small shift is introduced for the rank-deficient eigenequation to avoid the ill-conditioned eigenproblem. In the following, the results of the proposed well-conditioned eigenequation will be compared with the two traditional methods to investigate its accuracy in analysis of the free-free substructures.First, the substructural eigensolutions and eigensensitivity are calculated from the proposed well-conditioned eigenequation. The mode shapes of the rigid body modes are calculated from the geometric node locations of the first substructure according to (8). The well-conditioned eigenequation is formulated from (23). Based on the well-conditioned eigenequation, the eigensolutions of the first substructure are calculated. Since the zero-frequency eigensolutions are usually difficult to be accurately calculated, only the zero-frequency eigensolutions are listed in Table 1. The two-dimension substructure contains three rigid body modes. The eigenvalues of all the three rigid body modes are exactly 1 as expected. The proposed stiffness matrix K-(1) is well conditioned and full rank, based on which the residual flexibility is calculated from G=(K-(1))-1-Φm(1)(Λm(1))-1[Φm(1)]T. For comparison, all the eigensolutions of the first substructure are calculated, and the residual flexibility directly from the slave modes G~=Φs(1)(Λs(1))-1[Φs(1)]T is regarded as exact. The accuracy of the proposed substructure method in calculation of the residual flexibility is evaluated by the difference of the residual flexibility in terms of
(58)diff(G)=norm(G-G~)norm(G~),
where G is the residual flexibility calculated from the substructural stiffness matrix and master modes and G~ is the actual residual flexibility calculated from the slave modes. The difference of the residual flexibility between the proposed method and the exact one is 0.0% as shown in Table 1. The proposed method is exact in calculation of the substructural residual flexibility and eigensolutions. Based on the accurate eigensolutions and residual flexibility, the substructural eigensensitivity can be accurately calculated by common methods [10, 27].Table 1
Accuracy of the rigid body eigenvalues and the residual flexibility of the first substructure.
With the proposed well-conditioned eigenequation
With Moore-Penrose pseudo-inverse for rank-deficient matrix
With a small shift for rank-deficient matrix
Eigenvalue
Inverse of eigenvalue
Eigenvalue
Inverse of eigenvalue
Eigenvalue
Inverse of eigenvalue
Rigid body modes
1.00
1.00
1.48
E
-
09
6.74
E
+
08
0.1
10
1.00
1.00
4.62
E
-
09
2.17
E
+
08
0.1
10
1.00
1.00
1.64
E
-
08
6.07
E
+
07
0.1
10
Error of residual flexibility
0.00%
0.30%
1.60%Afterwards, the substructural solutions are calculated from the rank-deficient eigenequation directly by MATLAB eigensolver, in which the Moore-Penrose pseudoinverse is used for the singular and rank-deficient matrix. The eigensolutions are obtained and the zero-frequency eigensolutions are listed in Table1. Due to the numerical roundoffs, the rigid body modes are not perfect 0.0 Hz frequencies, and the rigid body eigenvalue is about 10−8. The inverse of the rigid body eigenvalues (Λm(1))-1 is a large value with order of 108. The residual flexibility is calculated by G=(K(1))+-Φm(1)(Λm(1))-1[Φm(1)]T, in which the inverse of the stiffness matrix is calculated from Moore-Penrose pseudoinverse. In consequence, the residual flexibility is calculated based on the Moore-Penrose pseudoinverse of the stiffness matrix and the inaccurate inverse of the eigenvalues. The difference in the residual flexibility between the proposed method and the exact ones is 0.3% in terms of (58). The accuracy of the eigensolutions and residual flexibility is significantly sensitive to the rank condition when carried out in floating-point arithmetic.Finally, the substructural eigensolutions and eigensensitivity are calculated from the rank-deficient eigenequation with a small shift introduced. The small shift is set to be 0.1, and the eigenvalues for the rigid body modes are 0.1 as shown in Table1. The residual flexibility is calculated by
(59)ΦsΛs-1ΦsT≅Φs(Λs+ε)-1ΦsT=(K+εM)-1-Φm(Λm+ε)-1ΦmT.
Since the small shift is far less than the slave eigenvalues (ε≪Λs), it is acceptable to be used to calculate the substructural residual flexibility. The difference between the calculated residual flexibility and the exact residual flexibility in terms of (58) is 1.6%. The small shift introduces a small error in the calculation of residual flexibility, and it will inevitably influence the accuracy of the global eigensolutions and eigensensitivity.The substructural eigensolutions and eigensensitivity are assembled to the objective function and sensitivity matrix of the global structure for model updating. The model updating process is performed based on the eigensolutions and eigensensitivity by the above three methods. In model updating, the simulated “experimental” modal data are obtained by intentionally introducing damages on some elements, and then the analytical model is updated to identify these damages [10, 30]. In the present paper, the simulated frequencies and mode shapes, which are treated as the “experimental” data, are calculated from the FE model by intentionally reducing the bending rigidity of Element 139 and Element 140 by 25% (denoted in Figure 3(a)). The first 10 “experimental” modes are available, and the measurements are obtained at the points and directions denoted in Figure 3(a). Both the “experimental” frequencies and mode shapes are utilized to update the analytical model. The mode shapes have been normalized with respect to the mass matrix.The first 30 modes are retained as master modes in each substructure to calculate the first 10 eigensolutions and eigensensitivities of the global structure. It is noted that using the proposed substructuring method, the eigensolutions and eigensensitivities are calculated based on the reduced equation (15) with size of 90 × 90, rather than on the original global system matrices with size of 408 × 408. The bending rigidities of all column elements are chosen as the updating parameters. Accordingly, there are 64 updating parameters in total. The optimization is processed with the trust-region method provided by the MATLAB Optimization Toolbox [29–33]. The algorithm can automatically select the steps and searching directions according to the objective function (discrepancy of eigensolutions) and the provided eigensensitivity matrices. To compare the accuracy of the above three methods, the model updating process stops when 12 iterations are performed for all the three methods.The identified changes of the elemental stiffness parameters are displayed in Figure4. The stiffness parameters of Element 139 and Element 140 are reduced by 25%, which agree with the simulated reduction in the elemental parameters. The three methods have different accuracy in calculation of the residual flexibility, eigensolutions, and eigensensitivity of the free-free substructure. Some small values observed in other elements are due to the errors in calculation of the substructural eigensolutions and eigensensitivity by the three methods.Location and severity of identified stiffness changes by the substructure-based model updating.
(a)
With proposed well-conditioned eigenequation
(b)
With Moore-Penrose pseudoinverse
(c)
With a small shift in eigenequationFigure4(a) shows that the proposed well-conditioned eigenequation is accurate to be used in substructure-based model updating, as it calculates the substructural residual flexibility, eigensolutions, and eigensensitivity with a high accuracy. The stiffness of Element 139 and Element 140 is identified as −25%. The identified changes of the elemental parameters are exactly consistent with the simulated damage. The proposed method has invisible values in other undamaged elements. On the other hand, the model updating results from traditional methods, which employ the Moore-Penrose pseudoinverse or a small shift for the rank-deficient eigenequation, introduce some small changes of undamaged elements as shown in Figures 4(b) and 4(c). And the introduction of a small shift for the rank-deficient eigenequation is less accurate than the use of Moore-Penrose pseudoinverse in calculation of substructural residual flexibility and eigensolutions. It is noted that a small shift 0.1 is introduced in this example. The selection of 0.1 in this example is not necessarily the best case, and a different selection of the shift value might contribute to better computational accuracy.Table1 and Figure 4 show that the proposed well-conditioned eigenequation has high accuracy in calculation of substructural eigensolutions and eigensensitivity for the free-free substructure, and it is accurate to be used in substructure-based model updating. The accurate calculation of eigensolutions and eigensensitivity is significant and helpful for the convergence of the model updating process. The efficiency of the proposed method in substructure-based model updating will be illustrated in the following case study.


### 4.2. Canton Tower

To illustrate the computational efficiency of the proposed substructuring method in large-scale structures, the FE model of the Canton Tower is employed here. The Canton Tower is a supertall structure of 610 m height. It consists of a main tower (454 m) and an antennary mast (156 m), as shown in Figure5(a). The main tower comprises a reinforced concrete inner tube and a steel outer tube of concrete-filled-tube (CFT) columns [34]. The outer tube consists of 24 CFT columns, uniformly spaced in an oval while being inclined in the vertical direction. The columns are interconnected transversely by steel ring beams and bracings. The analytical finite element (FE) model of the structure (Figure 5(b)) includes 8,738 three-dimensional elements, 3,671 nodes (each of which has six DOF), and 21,690 DOF.Canton Tower and the FE model.
(a)
Landscape view
(b)
Global model
(c)
Divided substructures
(d)
Concerned substructureThe global structure is divided into 10 substructures along the vertical direction as in Figure5(c). The “experimental” frequencies and mode shapes are simulated on the global structure by intentionally reducing the bending rigidity of 48 column elements of the outer tube in the local area (denoted in Figure 5) by 30%. The first 10 “experimental” modes are generated from the structure. The mode shapes are normalized to the mass matrix.The analytical model is updated by employing the substructuring method. The first 20 modes of the independent substructures are selected as the master modes to calculate the eigensolutions and eigensensitivity of the global structure. The bending rigidities of all the column elements of the outer tube in the concerned local area (the second substructure) are chosen as the updating parameters. Accordingly, there are a total of 144 updating parameters. The second substructure to the tenth substructure are free-free substructures after partition, and the system matrices of the nine substructures are singular and rank deficient. Three methods are utilized to handle the rank-deficient eigenproblem, namely, with the construction of the proposed well-conditioned eigenequation, with the Moore-Penrose pseudoinverse of rank-deficient matrix, and with a small shift for rank-deficient matrix.In each iteration, the substructural residual flexibility, eigensolutions, and its derivative matrix are calculated from the independent substructural model. The second substructure is taken as an example to illustrate the computation efficiency of the above three methods for the analysis of the rank-deficient eigenproblem. The system matrices of the second substructure take the size 2,736 × 2,736. Table2 shows the computation time of the substructure-based model updating process by the three methods for the analysis of rank-deficient eigenproblem. With the construction of the proposed well-conditioned eigenequation, it costs about 4.1590 seconds to calculate the first 20 eigensolutions of the second substructure via an ordinary personal computer with a 3.40 GHz CPU and 16 GB memory. The calculation of the substructural eigensensitivity with respect to the 144 updating parameters costs 655.59 seconds. The substructural eigensolutions and eigensensitivity are assembled to form the eigensolutions and eigensensitivity of the global structure for model updating process. One iteration of model updating consumes about 0.31 hours. Since the proposed full-rank well-conditioned eigenequation is accurate for the calculation of the substructural residual flexibility, eigensolutions, and eigensensitivity, the norm of objective function reaches 10−7 after 13 iterations as in Figure 6. The whole model updating process takes about 3.8 hours.Table 2
Computation time for substructure-based model updating.
With the proposed well-conditioned eigenequation
With Moore-Penrose pseudo-inverse for rank-deficient matrix
With a small shift for rank-deficient matrix
Calculation of substructural eigensolution (seconds)
4.1590
16.6342
4.1433
Calculation of substructural eigensensitivity (seconds)
655.59
1683.21
642.77
One iteration of model updating (hours)
0.3103
0.7328
0.3087
Entire time of substructure-based model updating (hours)
3.8
11.5
5.9Figure 6
Convergence of substructure-based model updating by three methods.Calculation of Moore-Penrose pseudoinverse for rank-deficient matrix is usually computationally expensive. If the Moore-Penrose pseudoinverse is used for the rank-deficient system matrix, the calculation of substructural eigensolution costs about 16.6342 seconds, and the calculation of substructural eigensensitivity with respect to 144 updating parameters costs about 16.6342 seconds as in Table2. One iteration of model updating consumes 0.7328 hours. The whole model updating process takes about 11.5 hours until the norm of objective function reaches 10−6 after 16 iterations as in Figure 6. The computation time is about 2.5 times of that consumed by the proposed method. The calculation of Moore-Penrose pseudoinverse of the rank-deficient system matrix is computationally expensive, while the results are not as accurate as those of the construction of the proposed well-conditional eigenequation.Finally, the substructural residual flexibility, eigensolutions, and eigensensitivity of the free substructures are calculated by introducing a small shift for rank-deficient eigenequation. The small shift is selected to be 0.1. Calculation of the substructural eigensolutions of the second substructure costs about 4.1433 seconds, and the calculation of substructural eigensensitivity costs about 642.77 seconds. The computation time is a little shorter than the method of the proposed well-conditioned eigenequation. However, the small shift in the original eigenequation inevitably introduces some errors. These errors hinder the convergence of model updating process. In Figure6, the norm of objective function reaches 10−5 after 20 iterations. It is noted that the convergence is difficult to achieve the precision of 10−7, as the small shift introduces errors in calculation of substructural solutions and hinders the convergence of model updating. To reach the insufficient accurate results, the model updating process costs 5.9 hours, which is longer than the proposed well-conditioned eigenequation.In consequence, the construction of the proposed well-conditioned eigenequation for the free-free substructure is accurate and efficient for the calculation of substructural residual flexibility, eigensolutions, and eigensensitivity and is thus efficient for the substructure-based model updating process. The proposed well-conditioned eigenequation is more accurate and efficient than the traditional Moore-Penrose pseudoinverse method for the analysis of the free-free substructures. As compared to the introduction of a small shift, the proposed well-conditioned eigenequation achieves a similar efficiency in calculation of substructural eigensolutions and eigensensitivity, but it has much higher accuracy. The proposed method is more efficient than that of the introduction of a small shift when they are utilized in the substructure-base model updating.


## 4.1. Three-Span Frame Structure

The accuracy of the proposed well-conditioned eigenequation for calculation of substructural residual flexibility, eigensolutions, and eigensensitivity will be illustrated by a frame structure. The global frame is shown in Figure3. The material constants are chosen as bending rigidity (EI)=170×106 Nm2, axial rigidity (EA)=2500×106 N, mass per unit length (ρA)=110kg/m, and Poisson’s ratio = 0.3. The frame is discretized into 160 two-dimensional beam elements each 2.5 m long, which results in 140 nodes and 408 DOF (N=408). The frame is disassembled into three substructures (NS=3) when it is torn at 8 nodes as shown in Figure 3. After division, there are 51, 55, and 42 nodes in the three substructures with the DOF of n1=153, n2=165, and n3=114, respectively.The frame structure. →: the points and directions of experimental measurement.
(a)
The global structure
(b)
The divided substructuresIn this example, the first substructure is analyzed to investigate the accuracy of the proposed well-conditioned eigensolver in calculation of substructural eigensolutions and eigensensitivity for a free-free substructure. The first substructure has 153 DOF. The first 30 modes are calculated as the master modes to assemble the eigenequation of the global structure, and the other slave modes are compensated by the residual flexibility. As the first substructure is free after partition, the system matrices of the first substructure are singular and rankdeficient. Traditionally, the Moore-Penrose pseudoinverse is usually used for the analysis of rank-deficient matrix, or a small shift is introduced for the rank-deficient eigenequation to avoid the ill-conditioned eigenproblem. In the following, the results of the proposed well-conditioned eigenequation will be compared with the two traditional methods to investigate its accuracy in analysis of the free-free substructures.First, the substructural eigensolutions and eigensensitivity are calculated from the proposed well-conditioned eigenequation. The mode shapes of the rigid body modes are calculated from the geometric node locations of the first substructure according to (8). The well-conditioned eigenequation is formulated from (23). Based on the well-conditioned eigenequation, the eigensolutions of the first substructure are calculated. Since the zero-frequency eigensolutions are usually difficult to be accurately calculated, only the zero-frequency eigensolutions are listed in Table 1. The two-dimension substructure contains three rigid body modes. The eigenvalues of all the three rigid body modes are exactly 1 as expected. The proposed stiffness matrix K-(1) is well conditioned and full rank, based on which the residual flexibility is calculated from G=(K-(1))-1-Φm(1)(Λm(1))-1[Φm(1)]T. For comparison, all the eigensolutions of the first substructure are calculated, and the residual flexibility directly from the slave modes G~=Φs(1)(Λs(1))-1[Φs(1)]T is regarded as exact. The accuracy of the proposed substructure method in calculation of the residual flexibility is evaluated by the difference of the residual flexibility in terms of
(58)diff(G)=norm(G-G~)norm(G~),
where G is the residual flexibility calculated from the substructural stiffness matrix and master modes and G~ is the actual residual flexibility calculated from the slave modes. The difference of the residual flexibility between the proposed method and the exact one is 0.0% as shown in Table 1. The proposed method is exact in calculation of the substructural residual flexibility and eigensolutions. Based on the accurate eigensolutions and residual flexibility, the substructural eigensensitivity can be accurately calculated by common methods [10, 27].Table 1
Accuracy of the rigid body eigenvalues and the residual flexibility of the first substructure.
With the proposed well-conditioned eigenequation
With Moore-Penrose pseudo-inverse for rank-deficient matrix
With a small shift for rank-deficient matrix
Eigenvalue
Inverse of eigenvalue
Eigenvalue
Inverse of eigenvalue
Eigenvalue
Inverse of eigenvalue
Rigid body modes
1.00
1.00
1.48
E
-
09
6.74
E
+
08
0.1
10
1.00
1.00
4.62
E
-
09
2.17
E
+
08
0.1
10
1.00
1.00
1.64
E
-
08
6.07
E
+
07
0.1
10
Error of residual flexibility
0.00%
0.30%
1.60%Afterwards, the substructural solutions are calculated from the rank-deficient eigenequation directly by MATLAB eigensolver, in which the Moore-Penrose pseudoinverse is used for the singular and rank-deficient matrix. The eigensolutions are obtained and the zero-frequency eigensolutions are listed in Table1. Due to the numerical roundoffs, the rigid body modes are not perfect 0.0 Hz frequencies, and the rigid body eigenvalue is about 10−8. The inverse of the rigid body eigenvalues (Λm(1))-1 is a large value with order of 108. The residual flexibility is calculated by G=(K(1))+-Φm(1)(Λm(1))-1[Φm(1)]T, in which the inverse of the stiffness matrix is calculated from Moore-Penrose pseudoinverse. In consequence, the residual flexibility is calculated based on the Moore-Penrose pseudoinverse of the stiffness matrix and the inaccurate inverse of the eigenvalues. The difference in the residual flexibility between the proposed method and the exact ones is 0.3% in terms of (58). The accuracy of the eigensolutions and residual flexibility is significantly sensitive to the rank condition when carried out in floating-point arithmetic.Finally, the substructural eigensolutions and eigensensitivity are calculated from the rank-deficient eigenequation with a small shift introduced. The small shift is set to be 0.1, and the eigenvalues for the rigid body modes are 0.1 as shown in Table1. The residual flexibility is calculated by
(59)ΦsΛs-1ΦsT≅Φs(Λs+ε)-1ΦsT=(K+εM)-1-Φm(Λm+ε)-1ΦmT.
Since the small shift is far less than the slave eigenvalues (ε≪Λs), it is acceptable to be used to calculate the substructural residual flexibility. The difference between the calculated residual flexibility and the exact residual flexibility in terms of (58) is 1.6%. The small shift introduces a small error in the calculation of residual flexibility, and it will inevitably influence the accuracy of the global eigensolutions and eigensensitivity.The substructural eigensolutions and eigensensitivity are assembled to the objective function and sensitivity matrix of the global structure for model updating. The model updating process is performed based on the eigensolutions and eigensensitivity by the above three methods. In model updating, the simulated “experimental” modal data are obtained by intentionally introducing damages on some elements, and then the analytical model is updated to identify these damages [10, 30]. In the present paper, the simulated frequencies and mode shapes, which are treated as the “experimental” data, are calculated from the FE model by intentionally reducing the bending rigidity of Element 139 and Element 140 by 25% (denoted in Figure 3(a)). The first 10 “experimental” modes are available, and the measurements are obtained at the points and directions denoted in Figure 3(a). Both the “experimental” frequencies and mode shapes are utilized to update the analytical model. The mode shapes have been normalized with respect to the mass matrix.The first 30 modes are retained as master modes in each substructure to calculate the first 10 eigensolutions and eigensensitivities of the global structure. It is noted that using the proposed substructuring method, the eigensolutions and eigensensitivities are calculated based on the reduced equation (15) with size of 90 × 90, rather than on the original global system matrices with size of 408 × 408. The bending rigidities of all column elements are chosen as the updating parameters. Accordingly, there are 64 updating parameters in total. The optimization is processed with the trust-region method provided by the MATLAB Optimization Toolbox [29–33]. The algorithm can automatically select the steps and searching directions according to the objective function (discrepancy of eigensolutions) and the provided eigensensitivity matrices. To compare the accuracy of the above three methods, the model updating process stops when 12 iterations are performed for all the three methods.The identified changes of the elemental stiffness parameters are displayed in Figure4. The stiffness parameters of Element 139 and Element 140 are reduced by 25%, which agree with the simulated reduction in the elemental parameters. The three methods have different accuracy in calculation of the residual flexibility, eigensolutions, and eigensensitivity of the free-free substructure. Some small values observed in other elements are due to the errors in calculation of the substructural eigensolutions and eigensensitivity by the three methods.Location and severity of identified stiffness changes by the substructure-based model updating.
(a)
With proposed well-conditioned eigenequation
(b)
With Moore-Penrose pseudoinverse
(c)
With a small shift in eigenequationFigure4(a) shows that the proposed well-conditioned eigenequation is accurate to be used in substructure-based model updating, as it calculates the substructural residual flexibility, eigensolutions, and eigensensitivity with a high accuracy. The stiffness of Element 139 and Element 140 is identified as −25%. The identified changes of the elemental parameters are exactly consistent with the simulated damage. The proposed method has invisible values in other undamaged elements. On the other hand, the model updating results from traditional methods, which employ the Moore-Penrose pseudoinverse or a small shift for the rank-deficient eigenequation, introduce some small changes of undamaged elements as shown in Figures 4(b) and 4(c). And the introduction of a small shift for the rank-deficient eigenequation is less accurate than the use of Moore-Penrose pseudoinverse in calculation of substructural residual flexibility and eigensolutions. It is noted that a small shift 0.1 is introduced in this example. The selection of 0.1 in this example is not necessarily the best case, and a different selection of the shift value might contribute to better computational accuracy.Table1 and Figure 4 show that the proposed well-conditioned eigenequation has high accuracy in calculation of substructural eigensolutions and eigensensitivity for the free-free substructure, and it is accurate to be used in substructure-based model updating. The accurate calculation of eigensolutions and eigensensitivity is significant and helpful for the convergence of the model updating process. The efficiency of the proposed method in substructure-based model updating will be illustrated in the following case study.


## 4.2. Canton Tower

To illustrate the computational efficiency of the proposed substructuring method in large-scale structures, the FE model of the Canton Tower is employed here. The Canton Tower is a supertall structure of 610 m height. It consists of a main tower (454 m) and an antennary mast (156 m), as shown in Figure5(a). The main tower comprises a reinforced concrete inner tube and a steel outer tube of concrete-filled-tube (CFT) columns [34]. The outer tube consists of 24 CFT columns, uniformly spaced in an oval while being inclined in the vertical direction. The columns are interconnected transversely by steel ring beams and bracings. The analytical finite element (FE) model of the structure (Figure 5(b)) includes 8,738 three-dimensional elements, 3,671 nodes (each of which has six DOF), and 21,690 DOF.Canton Tower and the FE model.
(a)
Landscape view
(b)
Global model
(c)
Divided substructures
(d)
Concerned substructureThe global structure is divided into 10 substructures along the vertical direction as in Figure5(c). The “experimental” frequencies and mode shapes are simulated on the global structure by intentionally reducing the bending rigidity of 48 column elements of the outer tube in the local area (denoted in Figure 5) by 30%. The first 10 “experimental” modes are generated from the structure. The mode shapes are normalized to the mass matrix.The analytical model is updated by employing the substructuring method. The first 20 modes of the independent substructures are selected as the master modes to calculate the eigensolutions and eigensensitivity of the global structure. The bending rigidities of all the column elements of the outer tube in the concerned local area (the second substructure) are chosen as the updating parameters. Accordingly, there are a total of 144 updating parameters. The second substructure to the tenth substructure are free-free substructures after partition, and the system matrices of the nine substructures are singular and rank deficient. Three methods are utilized to handle the rank-deficient eigenproblem, namely, with the construction of the proposed well-conditioned eigenequation, with the Moore-Penrose pseudoinverse of rank-deficient matrix, and with a small shift for rank-deficient matrix.In each iteration, the substructural residual flexibility, eigensolutions, and its derivative matrix are calculated from the independent substructural model. The second substructure is taken as an example to illustrate the computation efficiency of the above three methods for the analysis of the rank-deficient eigenproblem. The system matrices of the second substructure take the size 2,736 × 2,736. Table2 shows the computation time of the substructure-based model updating process by the three methods for the analysis of rank-deficient eigenproblem. With the construction of the proposed well-conditioned eigenequation, it costs about 4.1590 seconds to calculate the first 20 eigensolutions of the second substructure via an ordinary personal computer with a 3.40 GHz CPU and 16 GB memory. The calculation of the substructural eigensensitivity with respect to the 144 updating parameters costs 655.59 seconds. The substructural eigensolutions and eigensensitivity are assembled to form the eigensolutions and eigensensitivity of the global structure for model updating process. One iteration of model updating consumes about 0.31 hours. Since the proposed full-rank well-conditioned eigenequation is accurate for the calculation of the substructural residual flexibility, eigensolutions, and eigensensitivity, the norm of objective function reaches 10−7 after 13 iterations as in Figure 6. The whole model updating process takes about 3.8 hours.Table 2
Computation time for substructure-based model updating.
With the proposed well-conditioned eigenequation
With Moore-Penrose pseudo-inverse for rank-deficient matrix
With a small shift for rank-deficient matrix
Calculation of substructural eigensolution (seconds)
4.1590
16.6342
4.1433
Calculation of substructural eigensensitivity (seconds)
655.59
1683.21
642.77
One iteration of model updating (hours)
0.3103
0.7328
0.3087
Entire time of substructure-based model updating (hours)
3.8
11.5
5.9Figure 6
Convergence of substructure-based model updating by three methods.Calculation of Moore-Penrose pseudoinverse for rank-deficient matrix is usually computationally expensive. If the Moore-Penrose pseudoinverse is used for the rank-deficient system matrix, the calculation of substructural eigensolution costs about 16.6342 seconds, and the calculation of substructural eigensensitivity with respect to 144 updating parameters costs about 16.6342 seconds as in Table2. One iteration of model updating consumes 0.7328 hours. The whole model updating process takes about 11.5 hours until the norm of objective function reaches 10−6 after 16 iterations as in Figure 6. The computation time is about 2.5 times of that consumed by the proposed method. The calculation of Moore-Penrose pseudoinverse of the rank-deficient system matrix is computationally expensive, while the results are not as accurate as those of the construction of the proposed well-conditional eigenequation.Finally, the substructural residual flexibility, eigensolutions, and eigensensitivity of the free substructures are calculated by introducing a small shift for rank-deficient eigenequation. The small shift is selected to be 0.1. Calculation of the substructural eigensolutions of the second substructure costs about 4.1433 seconds, and the calculation of substructural eigensensitivity costs about 642.77 seconds. The computation time is a little shorter than the method of the proposed well-conditioned eigenequation. However, the small shift in the original eigenequation inevitably introduces some errors. These errors hinder the convergence of model updating process. In Figure6, the norm of objective function reaches 10−5 after 20 iterations. It is noted that the convergence is difficult to achieve the precision of 10−7, as the small shift introduces errors in calculation of substructural solutions and hinders the convergence of model updating. To reach the insufficient accurate results, the model updating process costs 5.9 hours, which is longer than the proposed well-conditioned eigenequation.In consequence, the construction of the proposed well-conditioned eigenequation for the free-free substructure is accurate and efficient for the calculation of substructural residual flexibility, eigensolutions, and eigensensitivity and is thus efficient for the substructure-based model updating process. The proposed well-conditioned eigenequation is more accurate and efficient than the traditional Moore-Penrose pseudoinverse method for the analysis of the free-free substructures. As compared to the introduction of a small shift, the proposed well-conditioned eigenequation achieves a similar efficiency in calculation of substructural eigensolutions and eigensensitivity, but it has much higher accuracy. The proposed method is more efficient than that of the introduction of a small shift when they are utilized in the substructure-base model updating.


## 5. Conclusions

This paper provides a deep look at the properties of a free-free structure. It addresses the difficulties associated with the analysis of a free-free substructure that were frequently encountered in the substructuring methods.The well-conditioned stiffness and flexibility matrices are formulated to overcome the difficulty in analyzing singular stiffness and flexibility matrices. The generalized stiffness and flexibility matrices are constructed to be the dual inversion of each other, such that the stiffness matrix and modal flexibility matrix are transformed into each other efficiently and effectively avoiding the expensive pseudoinverse. The proposed generalized stiffness matrix is full rank, which is helpful for the analysis of a free-free structure in many aspects, such as the analysis of eigenequation with singular stiffness matrix and the calculation of the modal flexibility, residual flexibility, and their derivatives.The accuracy of the proposed method for the analysis of the free-free substructures and its accuracy in the substructure-based model updating are verified through application to a frame structure. The construction of the proposed well-conditioned eigenequation is accurate in calculation of the substructural residual flexibility, eigensolutions, and eigensensitivity, and the substructure-based model updating results are exactly in agreement with the predefined damaged cases. The efficiency of the proposed method is illustrated through a supertall structure. The proposed full-rank well-conditioned eigenequation is more accurate and efficient than the Moore-Penrose pseudoinverse and the introduction of a small shift for the analysis of the free-free substructural model. The proposed full-rank well-conditioned eigenequation is accurate and efficient to be used in substructure-based model updating.Although the present research intends to assist the analysis of the free-free substructures in substructuring methods, the proposed well-conditioned eigenequation can be generalized to the analysis of a general free-free structure.



---

*Source: 706798-2013-12-04.xml*