Structural damage assessment and finite element model refinement using measured modal data

MISSING IMAGE

Material Information

Title:
Structural damage assessment and finite element model refinement using measured modal data
Physical Description:
xiv, 167 leaves : ill. ; 29 cm.
Language:
English
Creator:
Kaouk, Mohamed
Publication Date:

Subjects

Genre:
bibliography   ( marcgt )
theses   ( marcgt )
non-fiction   ( marcgt )

Notes

Thesis:
Thesis (Ph. D.)--University of Florida, 1993.
Bibliography:
Includes bibliographical references (leaves 162-166).
Statement of Responsibility:
by Mohamed Kaouk.
General Note:
Typescript.
General Note:
Vita.

Record Information

Source Institution:
University of Florida
Rights Management:
All applicable rights reserved by the source institution and holding location.
Resource Identifier:
aleph - 001933243
notis - AKA9310
oclc - 30784235
System ID:
AA00003663:00001

Full Text










STRUCTURAL DAMAGE ASSESSMENT AND FINITE ELEMENT
MODEL REFINEMENT USING MEASURED MODAL DATA














By

MOHAMED KAOUK


A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY

UNIVERSITY OF FLORIDA


1993

































A Mes Parents

Tout le m6rite de ce travail, s'il en est, vous revient.













ACKNOWLEDGEMENTS


I would like to express my sincere gratitude toward my advisor, Dr. David Zimmerman,

for his expert guidance, friendship, endless encouragement and support. I will forever be

indebted to him for inspiring me in my research and for the priceless education I acquired

from him. I am also grateful for the financial support he provided me during the course of my

graduate studies.

Words are not enough to express my deepest gratitude to my parents for their love, moral

and financial support; it is these things that have made this work possible. I also wish to

acknowledge my sisters and brothers for their continuous encouragement.

I would like to thank the members of my supervisory committee, Drs. Norman Fitz-Coy,

Lawrence Malvern, Bhavani Sankar, and Kermit Sigmon, for their helpful advice.

I wish to thank my colleagues of the Dynamic Systems and Control Laboratory for their

consideration and entertaining discussions. In particular, William Leath and Cinnamon

Larson for their friendship and encouragement.

I would like to thank the staff of the Aerospace Engineering, Mechanics, and

Engineering Science department for their assistance throughout the years, especially Shirley

Robinson for making my life easier during registration, John Young for his prompt responses

in fixing my computer problems, and Ronald Brown for his assistance in the machine shop

and for many stimulating discussions.

I would like to acknowledge the financial support received from Harris Corporation,

NASA/Florida Space Grant Consortium and Florida High Technology Council. I would like

to thank Dr. T. A. Kashangaki of the NASA Langley Research Center and Dr. S. W. Smith of

the University of Kentucky for providing the data of the Eight-Bay Truss used in this study.

Last, but not least, I am grateful to my good friends Joel Payabyab and Fadel Abdallah

for their continuous support.
















TABLE OF CONTENTS


paga

ACKNOWLEDGEMENTS ................................................................................ iii
LIST OF TABLES ................................................................................................. ix
LIST OF FIGURES ............................................................................................... x

A B ST R A C T ............................................................................................................ xiii
CHAPTERS

1 INTRODUCTION ................................................................................. 1

1.1 Finite Element Model Refinement ................................................. 1
1.1.1 O verview ........................................................................... 1
1.1.2 Literature Survey ............................................................. 3
1.2 Structural Damage Assessment ....................................................... 6
1.3 Objective of the Present Study ....................................................... 8

2 MATHEMATICAL PRELIMINARIES AND PRACTICAL
ISSUES TO THE PROBLEMS OF MODEL
REFINEMENT AND DAMAGE DETECTION ............................... 11
2.1 Introduction ..................................................................................... 11
2.2 The Eigenvalue Problem of Discrete Systems ................................ 11
2.2.1 Undamped Models .............................................................. 15
2.2.2 Proportional Damped Models ............................................ 16
2.3 Experimental Modal Analysis ....................................................... 16
2.4 Analytical/Experimental Model Dimensions Correlation ............. 19
2.4.1 Model Reduction Methods .................................................. 19
2.4.1.1 Static Reduction ................................................. 21
2.4.1.2 IRS Reduction ....................................................... 22
2.4.1.3 Dynamic Reduction ............................................... 22








2.4.2 Eigenvector Expansion Methods ........................................ 23
2.4.2.1 Dynamic Expansion ............................................ 23
2.4.2.2 Orthogonal Procrustes Expansion ........................ 24
2.5 Eigenvector Orthogonalization ..................................................... 26
2.5.1 Optimal Weighted Orthogonalization .................................. 27
2.5.2 Selective Optimal Orthogonalization .................................. 27
2.6 Load Path Preservation .................................................................. 28

3 INVERSE / HYBRID PROBLEM APPROACH FOR FINITE
ELEMENT MODEL REFINEMENT ................................................ 31

3.1 Introduction ................................................................................... 31
3.2 Theoretical Formulation ................................................................. 31
3.3 Numerical Illustrations .................................................................. 35
3.3 Sum m ary ......................................................................................... 37

4 SYMMETRIC EIGENSTRUCTURES ASSIGNMENT
MODEL REFINEMENT ALGORITHM ....................... 39
4.1 Introduction .................................................................................... 39
4.2 Problem Formulation ..................................................................... 39
4.2.1 Standard Eigenstructure Assignment Formulation ............. 40
4.2.2 Symmetric Eigenstructure Assignment Formulation ......... 42
4.2.3 Best Achievable Eigenvectors ............................................. 44
4.2.4 Selection of Bo : The Subspace Rotation Method .............. 46
4.3 Numerical Illustrations ................................................................. 47
4.3.1 Damage Detection: Kabe's Problem ................................... 48
4.3.1.1 Local to Global Mode Change ............................ 49
4.3.1.2 Consistent Modes ................................................. 52
4.3.2 Model Refinement of a Cantilever Beam:
Experimental Study .......................................................... 55
4.3.2.1 Modal Test Description ......................................... 55
4.3.2.2 Finite Element Model Description ........................ 56
4.3.2.3 Application of SEAMRA ...................................... 57
4.4 Discussion of the SEAMRA's Formulation ................................... 57
4.5 Sum m ary ......................................................................................... 59








5 DAMAGE LOCATION: THE SUBSPACE ROTATION
ALG O RITH M ................................................................................. 60

5.1 Introduction .................................................................................. 60
5.2 The Subspace Rotation Algorithm: The Direct Method ................ 60
5.3 The Subspace Rotation Algorithm: The Angle
Perturbation Method ................................................................. 63
5.4 Practical Issues .............................................................................. 64
5.4.1 Cumulative Damage Vectors ............................................... 64
5.4.2 Eigenvector Filtering Algorithm ......................................... 65
5.5 Sum m ary ......................................................................................... 66

6 THE MINIMUM RANK PERTURBATION THEORY............... 67
6.1 B background .................................................................................. 67
6.2 The Minimum Rank Perturbation Theory:
Theoretical Background ............................................................ 68
6.3 Damage Extent: Undamped Structures ......................................... 72
6.3.1 Damage Extent: Mass Properties ........................................ 73
6.3.2 Damage Extent: Stiffness Properties .................................. 74
6.3.3 Damage Extent: Mass and Stiffness Properties ................... 76
6.3.3.1 Application of the MRPT ...................................... 76
6.3.3.2 Decomposition of Matrix B ................................. 77
6.4 Damage Extent: Proportionally Damped Structures ..................... 79
6.4.1 Damage Extent: Stiffness and Damping Properties ............ 79
6.4.2 Damage Extent: Mass and Damping Properties .................. 82
6.4.3 Damage Extent: Mass and Stiffness Properties ................... 84
6.4.4 Damage Extent: Mass, Damping and Stiffness Properties .... 85
6.5 Damage Extent: Nonproportionally Damped Structures .............. 87
6.5.1 Damage Extent: Damping and Stiffness Properties ........... 89
6.5.2 Damage Extent: Mass and Damping Properties .................. 90
6.5.3 Damage Extent: Mass and Stiffness Properties .................. 92
6.6 Practical Issues .............................................................................. 94
6.6.1 The Concept of "Best" Modes ............................................. 95
6.6.2 Application of the Eigenvector Filtering Algorithm ............ 97
6.7 Sum m ary ....................................................................................... 97










7 VALIDATION AND ASSESSMENT OF THE SUBSPACE
ROTATION ALGORITHM AND THE MINIMUM
RANK PERTURBATION THEORY ............................................. 98

7.1 Introduction .................................................................................. 98
7.2 Kabe's Problem .......................................................................... 98
7.2.1 Problem Description ........................................................... 98
7.2.2 Damage Location ............................................................... 100
7.2.3 Damage Extent ............................................................... 103
7.3 Damage Detection: Fifty-Bay Two-Dimensional
Truss: Undamped FEM ........................................................ 105
7.3.1 Problem Description .......................................................... 105
7.3.2 Damage Location ............................................................... 106
7.3.3 Damage Extent .................................................................. 107
7.4 Experimental Study: The NASA 8-Bay Truss ............................ 112
7.4.1 Problem Description .......................................................... 112
7.4.2 Refinement of the Original FEM ....................................... 115
7.4.3 Damage Location ............................................................... 118
7.4.4 Damage Extent .................................................................. 120
7.4.4.1 The Brute Force Method ..................................... 121
7.4.4.2 The Damage Consistent Method ........................ 121
7.4.4.3 Application of the Eigenvector Filtering
Algorithm ........................................................ 122
7.5 Experimental Study: Mass Loaded Cantilevered Beam .............. 140
7.5.1 Problem Description .......................................................... 140
7.5.2 Analytical and Experimental Models Dimension
C orrelation ................................................................... 141
7.5.3 Refinement of the Original FEM ...................................... 141
7.5.4 Damage Location .............................................................. 142
7.5.5 Damage Extent .................................................................. 143
7.6 Fifty-Bay Two-Dimensional Truss: Proportionally
D am ped FEM ............................................................................. 144
7.6.1 Problem Description .......................................................... 144
7.6.2 Damage Location ............................................................... 145
7.6.3 Damage Extent ................................................................... 146








7.7 Eight-Bay Two-Dimensional Mass-Loaded
Cantilevered Truss ................................................................... 148
7.7.1 Problem Description .......................................................... 148
7.7.2 Proportionally Damped Configuration: Damage of
Small Order of Magnitude ............................................... 149
7.7.2.1 Damage Location ................................................. 150
7.7.2.2 Decomposition of Matrix B ................................. 150
7.7.2.3 Damage Extent .................................................... 152
7.7.3 Undamped Configuration: Damage of Large
Order of Magnitude ......................................................... 153
7.7.3.1 Noise Free Eigendata ........................................... 153
7.7.3.2 Noisy Eigendata .................................................. 154
7.8 Sum m ary ........................................................................................ 156

8 CONCLUSION AND SUGGESTIONS FOR FUTURE WORK ......... 159

REFEREN CES .................................................................................................... 162

BIOGRAPHICAL SKETCH ............................................................................... 167















LIST OF TABLES


Table Pa


3.1 Kabe's Problem: Elemental Stiffness Components ................................. 38

4.1 Properties of the Cantilever Beam ............................................................ 56

4.2 Measured Natural Frequencies and Damping Ratios
of the Cantilever Beam ............................................................................ 56

4.3 Measured Mode Shapes of the Cantilevered Beam ................................... 56

7.1 Fifty-Bay Truss: Summary of Damage Detection
Results using the M RPT ........................................................................ Ill

7.2 Mass Properties of the Eight Bay Truss .................................................... 113

7.3 Strut Properties of the Eight Bay Truss ..................................................... 114

7.4 NASA 8-Bay: Truss Damage Case Definitions ..................................... 115

7.5 Comparison of Analytical and Experimental Frequencies ....................... 116

7.6 Summary of the Filtering Process for Single Member Damage Cases ......... 123

7.7 NASA 8-Bay Truss: Summary of the Damage Assessment Results ........... 139

7.8 Mass Loaded Cantilevered Beam Properties ......................................... 140

7.9 Analytical and Experimental Frequencies of the "Healthy" Structure ...... 142

7.10 50-Bay 2-Dimensional Truss: Summary of the Percentage
Error with Respect to the Exact Damage .......................................... 147

7.11 Problem 7.7: Percentage Error of Damage Estimate with
Respect to Exact Dam age ...................................................................... 156














LIST OF FIGURES


Figure Page


1.1 Overview of Finite Element Model Refinement ...................................... 2

1.2 Overview of FEM Model Refinement Process Used for
D am age A ssessm ent ............................................................................ 7
2.1 Components of a Vibration Measurement System for
M odal A analysis .................................................................................. 17
2.2 A simple Experimental Modal Analysis Setup ...................................... 18

2.3 Flow Chart of the Iterative Load Preservation Path Algorithm ............. 29

3.1 Kabe's Problem: Analytical Test Structure .............................................. 36

4.1 Best Achievable Eigenvector Projection ................................................ 45

4.2 Rotation of the Achievable Subspace ..................................................... 47

4.3 K abe's Problem ..................................................................................... 48

4.4 Results for Kabe's Problem using the 1st Mode,
Full Eigenvector ............................................................................... 50

4.5 Results for Kabe's Problem Modes 1, 2, 3 and
Eigenvectors Components 1, 2, 3 ......................................................... 51

4.6 Results for Kabe's Problem using Load Path Preservation,
Modes 1, 2, 3 and Eigenvectors Components 1, 2, 3 ............................. 52

4.7 Results for Kabe's Problem using the 1 st Mode,
Full Eigenvector ................................................................................. 53

4.8 Results for Kabe's Problem Modes 1,2,3 and
Eigenvectors Components 1, 2, 6 ........................................................ 54

4.9 Experimental Cantilever Beam ............................................................. 55

4.10 Experimental and Analytical Frequency Response Function
of the Cantilever Beam ......................................................................... 58








7.1 K abe's Problem ...................................................................................... Y'

7.2 Kabe's Problem: Damage Location Results using the Subspace
Rotation Direct Method with the Eigendata of the 1st Mode ............... 101

7.3 Kabe's Problem: Damage Location Results using the Angle
Perturbation Method with the Eigendata of the 1st Mode ................... 101

7.4 Kabe's Problem: Damage Location Results using Lin's
Algorithm with the Eigendata of the 1st Mode .................................... 102

7.5 Kabe's Problem: Damage Location Results using the Angle
Perturbation Method with the Eigendata of the 1st and 2nd Modes ....... 102

7.6 Kabe's Problem: Damage Location Results using Lin's
Algorithm with the Eigendata of the 1st and 2nd Modes .................... 103

7.7 Kabe's Problem: Damage Extent Results using the MRPT
with the Eigendata of M ode 2 ............................................................... 104

7.8 Kabe's Problem: Damage Extent Results using Baruch's Method ........... 105

7.9 Fifty-Bay Two-Dimensional Truss ...................................................... 106

7.10 Fifty-Bay Truss: Damage Location Results using the Subspace
Rotation Algorithm with the Eigendata of the First Ten Modes .......... 107

7.11 Fifty-Bay Truss: Damage Extent Results using the
MRPT with the Eigendata of Modes 8 and 9 ..................................... 109

7.12 Fifty-Bay Truss: Damage Extent Results using the
MRPT with the Eigendata of the First Ten Modes ............................. 110

7.13 Fifty-Bay Truss: Damage Extent Results using Baruch's Algorithm ....... 110

7.14 The NASA Eight-Bay Hybrid-Scaled Truss: Damage Cases ................ 112

7.15 The NASA 8-Bay Truss: Lacing Pattern ................................................ 113

7.16 NASA 8-Bay Truss: Typical Frequency Response Comparison .............. 117

7.17 NASA 8-Bay Truss: Perturbation to the Original Stiffness
Matrix that Resulted From the Refinement Process .............................. 118

7.18 NASA 8-Bay Truss: Cumulative Damage Vector Associated
w ith C ase F ......................................................................................... 124

7.19 NASA 8-Bay Truss: Damage Assessment of Case A .............................. 125

7.20 NASA 8-Bay Truss: Damage Assessment of Case C .............................. 126










7.21 NASA 8-Bay Truss: Damage Assessment of Case D .................................. z/

7.22 NASA 8-Bay Truss: Damage Assessment of Case E ............................... 128

7.23 NASA 8-Bay Truss: Damage Assessment of Case G ............................... 129

7.24 NASA 8-Bay Truss: Damage Assessment of Case H ............................... 130

7.25 NASA 8-Bay Truss: Damage Assessment of Case I ............................ 131

7.26 NASA 8-Bay Truss: Damage Assessment of Case J ............................ 132

7.27 NASA 8-Bay Truss: Damage Assessment of Case K ............................ 133

7.28 NASA 8-Bay Truss: Damage Assessment of Case L ............................ 134

7.29 NASA 8-Bay Truss: Damage Assessment of Case M ............................ 135

7.30 NASA 8-Bay Truss: Damage Assessment of Case N ............................ 136

7.31 NASA 8-Bay Truss: Damage Assessment of Case 0 ............................ 137

7.32 NASA 8-Bay Truss: Damage Assessment of Case P ............................ 138

7.33 The Mass Loaded Cantilevered Beam .................................................... 140

7.34 Mass Loaded Cantilevered Beam: Damage Assessment ....................... 143

7.35 50-Bay 2-Dimensional Truss .................................................................. 145

7.36 50-Bay 2-Dimensional Truss: Damage Location .................................... 146

7.37 50-Bay 2-Dimensional Truss: Damage Extent ....................................... 148

7.38 The Eight-Bay Two-Dimensional Mass-Loaded Cantilevered Truss ....... 149

7.39 Problem 7.7: Cumulative Damage Location Vector: First Four Modes ...... 150

7.40 Problem 7.7: Cumulative Vectors Associated with the Exact
and Computed Bin, Bd, Bk: First Three Modes .................................... 151

7.41 Problem 7.7: Exact and Computed AMd, ADd, AKd ............................... 152

7.42 Problem 7.7: Cumulative Vector Associated with B, Bm, Bd, Bk
B Computed using Modes 1-4
Bm, Bd, Bk Computed using Modes 3,4 & 5 ..................................... 155

7.43 Problem 7.7: Exact and Computed AMd, ADd, AKd ............................... 157













Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirement for the Degree of Doctor of Philosophy

STRUCTURAL DAMAGE ASSESSMENT AND FINITE ELEMENT
MODEL REFINEMENT USING MEASURED MODAL DATA

By

Mohamed Kaouk

August 1993


Chairperson: Dr. David C. Zimmerman
Major Department: Aerospace Engineering, Mechanics and Engineering Science

This study investigates the problems of model refinement and structural damage

assessment. The essence of the model refinement problem is to adjust finite element models

(FEMs) of structures with the intent of producing a correlation between experimental and

analytical modal properties. Recently, the framework of model refinement has been adopted

to determine the location and extent of structural damage. Damage will result in changes to

the modal properties of the healthy structure. A further refinement of an accurate FEM of the

structure using damaged modal parameters is expected to generate adjustments to the FEM at

locations associated with damage. Analysis of these adjustments can then be used to assess

damage. In this investigation, four algorithms relevant to the subjects of model refinement

and damage detection are presented. The development of a model refinement algorithm with

its basis in inverse theory is first presented. The algorithm, termed the inverse/hybrid

approach, is illustrated in a comparative computer simulated study. Next, an existing

eigenstructure assignment model refinement algorithm is improved to better approach the

damage assessment problem. The enhanced algorithm is evaluated and compared to other

techniques using simulated and experimental data. The algorithm is shown to perform well








in assessing damage and refining FEMs. A damage location algorithm that bypasses the

general framework of model refinement is discussed. The damage location algorithm,

termed the subspace rotation, is similar to the modal force error criteria proposed by several

researchers. Greater insight to the modal force error criteria, along with a new viewpoint that

reduces the effects of measurement noise, is discussed. Furthermore, an efficient damage

extent algorithm based on a minimum rank perturbation theory (MRPT) is developed. The

formulation of the MRPT is consistent with the effect of most structural damage on FEMs.

The characteristics of the subspace rotation algorithm and the minimum rank perturbation

theory are illustrated using simulated and experimental testbeds. The decomposition of the

damage assessment problem into location and extent subproblems is shown to be

advantageous from both for computational efficiency and for engineering insight.














CHAPTER 1
INTRODUCTION




1.1 Finite Element Model Refinement


1.1.1 Overview

An important aspect in the design process of a structure is the evaluation of its

performance under expected dynamic loading conditions. Dynamic performance can be

analyzed by using either analytical or experimental techniques. Experimental analyses are

generally very costly, time consuming and can encounter technical difficulties. One example

of a commonly encountered technical barrier is the ground base reproduction of a weightless

environment for the testing of space structures. This shortcoming and the advent of

computers have sparked a growing interest in the uses of analytical techniques. This type of

analysis utilizes a mathematical model of the actual structure for computer simulated

evaluation of structural performance. Presently, most analytical models used in industry are

finite element models generated by using Finite Element Methods. The accuracy of a finite

element model (FEM) can be improved to some extent by increasing the number of degrees

of freedom (DOFs) included in the model. However, accuracy of the FEM may still be

lacking due to uncertainties in material properties, geometric layout and fabrication induced

errors. Thus, it is essential to "validate" these FEMs prior to their acceptance as a basis for

analysis. One way to validate them is to compare their modal properties (eigenvalues/

eigenvectors) with the measured modal parameters of the actual structure as obtained from

experimental modal analysis (EMA). A FEM is acceptable when these two sets of modal

properties are in agreement. Unfortunately, this agreement rarely occurs. As a result, the








FEM must be adjusted to produce a correlation of analytical and experimental modal

properties. An overview of this procedure is schematically illustrated in Figure 1.1

(Zimmerman and Smith, 1992).


Experimental
Analytical Frequency
Response
Mass Functions
Damping or
Stiffness Hankel
Matrices / Matrices







Modal Parameters j--- Modal Parameters



Yes

Done
(Confidence in
FEM established)



Figure 1.1 Overview of Finite Element Model Refinement
(Zimmerman and Smith, 1992)



In the past, the process of adjusting FEMs was performed on an ad-hoc basis aided by

engineering experience. This practice was naturally time consuming and in most cases

inadequate for large-order, complex FEMs that are typically needed for accurate dynamic

modelling. In the past years, growing interest has been focused on developing systematic

procedures to produce correlated FEMs. These efforts have resulted in the development of a

large number of algorithms. Finite element model adjustment procedures have been








commonly referred to in the literature as, FEM refinement, FEM adjustment, FEM

correction, FEM correlation, and FEM identification.

1.1.2 Literature Survey

Comprehensive literature surveys covering a large portion of the work that addresses the

model refinement problem can be found in the book chapter by Zimmerman and Smith

(1992), and in the papers by Ibrahim and Saafan (1987) and Heylen and Saas (1987). For

completeness and to properly underline the objective of the current study, it is appropriate to

include a brief survey of the development made in this area.

The concept of using experimental modal data in analytical studies was initiated by an

early work presented by Rodden (1967). In his work, Rodden explored the possibilities of

generating mass and stiffness property matrices by using experimentally measured modal

data. The experimentally generated mass and stiffness matrices were nonsymmetric. Brock

(1968) improved the work of Rodden by proposing a strategy to insure symmetry of property

matrices.

The essence of the model refinement concept, as adopted by most researchers, is to

modify the finite element model (FEM) of systems with the intention of producing a

correlation of test and analytical modal parameters. In general, this task has been approached

through two different philosophies. The first amounts to modifying globally the finite

element property matrices of the system. This approach is commonly referred to as the

property matrix update. The alternative approach is to individually correct parameters of

each element of the finite element model. This will be addressed here as the physical

parameter update.

Property matrix update algorithms can be classified into two subclasses of

methodologies: optimal-matrix update and control-based eigenstructure assignment

techniques.








The basic philosophy of the optimal-matrix update is to minimize the correction to the

FEM property matrices to accomplish the analytical/test modal correlation. The pioneering

work in this area can be credited to Baruch and Bar Itzhack (1978). In their formulation, by

assuming that the mass matrix is correct, the refinement of the stiffness matrix of an

undamped FEM is cast as a constrained minimization problem. The objective of their

formulation is to determine the minimal Frobenius-norm symmetric stiffness adjustment

that satisfies the eigenvalue problem in terms of the experimentally measured eigenvalues

and eigenvectors. A computationally efficient closed form solution was developed for the

updated stiffness matrix. Berman and Nagy (1983) extended the Baruch and Bar Itzhack

approach to the refinement of both the mass and stiffness matrices. This same methodology

was further extended by Fuh et al. (1984) to update the mass, damping and stiffness matrices

of damped FEMs. Fuh and his colleagues used cross-orthogonality relationships to correct

the mass and damping matrices and a constrained minimization problem similar to the one

proposed by Baruch and Bar Itzhack to determine the refined stiffness matrix. The problem

of adjusting the mass, damping and stiffness matrices was also attempted by Hanagud et al.

(1984). In their approach, all three property matrices of a nonproportionally damped FEM

are incorporated in the constrained minimization problem.

The previously stated algorithms do not preserve the load path (sparsity pattern) of the

original analytical property matrices. Kabe (1985) proposed a reformulation of the Baruch

and Bar Itzhack algorithm that constrains the updated stiffness matrix to preserve the load

path of the original stiffness matrix. In addition, in his reformulation, he adopted a percent

weighting on the stiffness changes instead of mass matrix weighting as used in the Baruch

and Bar Itzhack algorithm. Kammer (1987) proposed an alternative solution approach to the

problem defined by Kabe that uses projection matrix theory and the Moore-Penrose

generalized inverse. The algorithms formulated by Kabe and Kammer involve an excessive

amount of computational effort. Another alternative and more efficient formulation of the








Kabe problem that utilized a generalization of Marwill-Toint updates was developed by

Smith and Beattie (1991).

The other subclass of property matrix updates is based on the general framework of

control-based eigenstructure assignment algorithms. The essence of this approach is to

determine pseudo-controllers that would assign the experimentally measured modal

parameters to the original analytical FEM. The pseudo-controllers are then translated into

matrix adjustment to the original property matrices. This approach was initially formulated

by Inman and Minas (1990) to adjust the damping and stiffness matrices of the original FEM.

In their formulation, the symmetry of the adjusted property matrices was enforced through

an iterative process that involves a numerical nonlinear optimization process. Zimmerman

and Widengren (1989, 1990) proposed a formulation that replaces the iterative process of

Inman and Minas with a generalized algebraic Riccati equation. Further details about this

approach are treated in Chapter 4.

The alternative philosophy to the model refinement problem is the physical parameter

update. The basic idea of this alternative approach consists of adjusting some or all structural

physical parameters to produce test/analysis modal correlation. Structural parameters are

design variables such as material densities, modulus of elasticity, cross-sectional areas,

element lengths and others. This type of procedure results in corrected FEMs that are

consistent with the framework of the original FEM.

Most methodologies that have adopted this philosophy have used the sensitivity

derivatives of the system eigenvalues and/or eigenvectors with respect to the structural

parameters. Generally, the refinement process amounts to solving for the corrected

structural parameters through an iterative optimization problem directed by the sensitivity

derivatives. Some notable work in this area was accomplished by Collins et al. (1974), Chen

and Garba (1980), Adelman and Haftka (1986), Creamer and Hendricks (1987), Flanigan

(1991), Martinez et al. (1991), to name only a few. In these studies, different sensitivity

formulations and iterative optimization strategies are proposed. A direct approach that








bypasses the use of closed-form sensitivity derivatives was presented in the papers of Hajela

and Soeiro (1990) and Soeiro (1990).

White and Maytum (1976) set forth another alternative methodology that can be

classified in the physical parameter update group. In their approach, adjustments to the

stiffness matrix are viewed as a collection of known submatrices, elements or group of

elements. Each submatrix is multiplied by an unknown scaling factor. Correction of the

original stiffness matrix amounts to determining the scaling factors that would improve the

correlation of analytical and experimental modal properties. The scaling factors can be

physically viewed as functions of the stiffness structural parameters. A large

underdetermined system of equations in terms of the scaling factor is generated. The scaling

factors are then computed by using a pseudo-inverse. An improvement to the White and

Maytum approach is presented in the paper by Lim (1990).



1.2 Structural Damage Assessment


Structures, in general, are prone to structural damage during their service lives, that

could affect their overall performance and could result in catastrophic failures. This is of

critical concern when expensive systems and/or human lives are at stake. On this basis, it is

highly desirable to monitor the structural health of systems such as airplanes, space systems,

bridges, buildings, oil platforms, among others, to prevent such catastrophic events.

Currently, the most common structural health monitoring methods involve visual

inspections supplemented with X-rays, acoustic emission, magnetic resonance and

ultrasonic testing. These approaches can be time consuming, costly and difficult to perform

on inaccessible structural components. Additionally, none of these approaches provide a

quantitative assessment of the magnitude of the damage.

In response to these shortcomings, structural health monitoring approaches based on

the framework of existing model refinement techniques have been recently proposed by

several researchers (Smith and Hendricks, 1987; Chen and Garba, 1988; Ricles and








Kosmatka, 1992; Zimmerman and Kaouk, 1992a,b; Kaouk and Zimmerman, 1993a,b).

These techniques utilize the vibration signature of the pre- and post- damaged structure, in

conjunction with an analytical model of the original structure, to assess both the location and

extent of the structural damage. The pre-damage modal parameters are used to correct

(refine) the original finite element model (FEM) to determine an "accurate" reference

baseline. Once damaged, the post-damage modal properties of the structure are used to

further "refine" the refined analytical model. This results in perturbations to the refined

analytical model. Analysis of the perturbations could indicate the damage location and

extent. An overview of the application of model refinement algorithms in assessing

structural damage is summarized in Figure 1.2.



Undamaged Model Experimental
(Correlated) Frequency
Response
Mass Functions
Damping or
Stiffness Hankel
Matrices Matrices



No



Modal Parameters Modal Parameters



Yes

Done
Structure Healthy



Figure 1.2 Overview of Finite Element Model Refinement
Process Used for Damage Assessment








Notable exceptions to the direct use of FEM refinement algorithms to the damage

detection problem are the work of Lin (1990), Ojalvo and Pilon (1988), and Gysin (1990). In

the work of Lin, a flexibility matrix is determined using experimental data. This matrix is

then multiplied by the original stiffness matrix, with those rows and/or columns that differ

significantly from a row and/or column of the identity matrix indicating which degrees of

freedom have been most affected by the damage. It is then assumed that damage has occurred

in structural elements connecting those degrees of freedom. An overview of the work by

Ojalvo and Pilon (1988), and Gysin (1990) is given in Chapter 5.

Although the problems of damage assessment and model refinement are related to one

another, they have quite different characteristics. In essence, the model refinement concept

is based on the philosophy of the minimum change to the original FEM. Naturally, this

minimum change constraint has a tendency of smearing the changes throughout the entire

FEM. This aspect is inconsistent with the effects of structural damage on FEMs which are

usually localized perturbations of possibly large magnitude. In fact, structural damage often

occurs at discrete locations and only affects a few elements of the FEM.



1.3 Objective of the Present Study


The present study investigates the development of new and promising model refinement

and damage detection methodologies. Although considerable research has been done in

these areas, no methodology has been fully successful in dealing with the refinement or the

damage detection problem of "real life" systems. The main objective of this investigation is

to formulate efficient model refinement algorithms that are consistent with the effect of

structural damage.

In chapter 2, practical concepts and issues related to the general areas of model

refinement are presented. First, the concept of the eigenvalue problem of discrete structural

models is reviewed with emphasis on the associated cross-orthogonality conditions. A brief

discussion of modal analysis follows. The problem associated with incomplete eigenvector








measurements is then investigated and already existing techniques to deal with this problem

are presented. Furthermore, eigenvector orthogonalization techniques useful to a large

number of model refinement algorithms are reported. Finally, an iterative approach to the

problem of load path preservation encountered in a large number of property matrix update

algorithms is discussed.

In chapter 3, the author proposes and formulates a new model refinement algorithm.

The newly developed algorithm, named the inverse/hybrid method, is naturally based on the

inverse problem. The model type under consideration is proportionally damped and the

inaccuracies of the analytical model is assumed to be due to errors in the modelling of the

damping and stiffness properties. A complete hybrid set of modal data is generated by

approximating the unmeasured modal information by the corresponding analytical modes.

Symmetry of the corrected stiffness and damping matrices is enforced by mass

orthogonalization the complete hybrid eigenvector matrix. An orthogonalization procedure

that assigns more confidence on the measured eigenvectors is proposed. A comparative

study of the inverse/hybrid method and the algorithm proposed by Baruch and Bar Itzhack

(1978) shows that both algorithms give similar results. However, it is shown that the

inverse/hybrid approach is less computationally efficient.

Chapter 4 proposes improvements to the symmetric eigenstructure assignment model

refinement algorithm (SEAMRA) formulated by Zimmerman and Widengren (1989). The

author develops a technique to enhance eigenvector assignability. The technique, termed the

subspace rotation method, is based on rotating the achievable eigenvector subspaces into the

experimental eigenvectors. The subspace rotation method results in both a decrease in the

computational burden and an increase in the accuracy of the assigned eigenvector. The

enhanced SEAMRA is then evaluated and compared to other algorithms using both

computer simulated and experimental testbeds. It is shown that the enhanced algorithm is

suitable for damage detection applications.








In Chapter 5, an efficient damage location algorithm that utilizes modal data

information but bypasses the general framework of the model refinement problem is

presented. This location algorithm is an outgrowth of the subspace rotation method used to

enhance eigenvector assignability of the SEAMRA. The proposed location algorithm is

similar to the modal force error criteria presented by several researchers (Ojalvo and Pilon,

1988; Gysin, 1990). Further interpretation of the algorithm operation is given. Additionally,

the author proposes and formulates a new viewpoint that reduces the effect of measurement

noise for certain types of structures. Based on this formulation, an eigenvector filtering

algorithm is also developed.

Chapter 6 presents the formulation of computationally attractive damage extent

algorithms that are based on a minimum rank perturbation theory (MRPT) developed by the

author. The formulation of the MRPT is consistent with the effect of many classes of

structural damage on FEMs. Several MRPT based algorithms are formulated to

accommodate structures with undamped, proportionally damped, and nonproportionally

damped FEMs. For each type of FEM, several damage scenarios are considered.

Discussions of the characteristics and properties of the developed algorithms are presented

along with practical issues that can be used to improve their performance.

In Chapter 7, the algorithms developed in Chapters 5 and 6 are demonstrated and

evaluated using both computer simulated and actual experimental data. The main objective

of these example problems is to illustrate the potential of these algorithms in assessing

structural damage. In all examples, the location of damage is first determined by using the

algorithm presented in Chapter 5. An MRPT based algorithm is then utilized to assess the

extent of the damage. Several key points made throughout the formulation in Chapters 5 and

6 are emphasized. In particular, it is shown that the damage extent calculations can be greatly

enhanced by making use of the damage location algorithm.

In Chapter 8, summaries and conclusions of the issues discussed in this study are

presented along with some suggestions for future work.














CHAPTER 2
MATHEMATICAL PRELIMINARIES AND PRACTICAL ISSUES RELATED TO
THE PROBLEMS OF MODEL REFINEMENT AND DAMAGE DETECTION




2.1 Introduction

In this chapter, general mathematical preliminaries and practical issues relevant to the

areas of model refinement and damage detection are discussed. In Section 2.2, an overview

of the general eigenvalue problem for discrete models is given. Further, the concept of modal

analysis is introduced in Section 2.3. Section 2.4 deals with the concept of incomplete

degrees of freedom measurement. Two alternative approaches are discussed as possible

solutions to the incomplete measurement problem. Two eigenvector orthogonalization

techniques are discussed in Section 2.5. Finally, in Section 2.6, an algorithm to preserve

sparsity in updated property matrices is presented. The concepts discussed in this chapter

will be frequently referred to in the course of the forthcoming chapters.


2.2 The Eigenvalue Problem of Discrete Systems

In practice, most engineering structures are continuous systems with spatially

distributed material properties. The vibration motion of these structures in terms of their

distributed properties is usually governed by one or several partial differential equations.

With complex systems, both the development and the analysis of partial differential

equations of motion are tedious and in many cases impossible. These shortcomings and the

advent of digital computers have motivated the development of approximate modelling of

continuous systems in order to simplify the equations of motion. The general idea behind

these approximations is to represent the exact distributed model of a system by a discrete one.








This concept is known as spatial discretization which eliminates the continuous spatial

dependence of the distributed properties. The discrete model is of finite order and is

described by a finite number of variables known as degrees of freedom (DOFs). The number

of DOFs used in the discrete model depends on the desired accuracy in representing the

continuous model. Commonly, the vibration motion of systems in terms of their discrete

models is described by a set of simultaneous ordinary differential equations that are usually

simpler to develop and analyze than the partial differential equations of the continuous

system. One of the most commonly used approximate discrete modelling techniques is the

finite element method. The model generated by a finite element method is called a finite

element model (FEM). More elaborate discussions of the concepts of continuous and

discrete models as well as finite element methods are covered in detail in the books of

Meirovitch (1986, 1980), Inman (1989), and Hughes (1987).

Commonly, the free vibration motion of a structure in terms of an nth order discrete

model is represented by the following set of simultaneous ordinary differential equations

Mw(t) + Dw(t) + Kw(t) = 0 (2.2.1)

where M, D and K are termed, respectively, the mass, damping and stiffness matrices. They

are models of the mass, damping and stiffness properties of the structure. Since the system

model is order n, these matrices are of dimension nxn and are generally real. The variable

w(t) represents the n displacements of the n-DOF model of the structure. The overdots

represent differentiation with respect to time. The mass matrix, M, is always symmetric

positive definite. The stiffness (K) and damping (D) matrices of nongyroscopic and

noncirculatory systems are symmetric. In general, the modelling of the mass and stiffness

properties of the structure is simpler and more accurate than the modelling of the damping

properties. In the forthcoming discussion, it is assumed that the system under consideration

is nongyroscopic and noncirculatory. The standard solution to Eq. (2.2.1) is of the form


w(t) = vext


(2.2.2)








where v is a constant nxl vector and k is a constant scalar. Substituting Eq. (2.2.2) into Eq.

(2.2.1) and dividing the resultant equation by ext yields the condition

My X2 + Dv X + Kv = 0 (2.2.3)

There are n sets of nontrivial complex conjugate solutions (Xi, y,) to Eq. (2.2.3). Note that

since the property matrices (M, D, K) are real, if (iQ, yi) is a solution set to Eq. (2.2.3), the

complex conjugate of that set is also a solution. The problem of solving for these solutions

is commonly known as the eigenvalue problem and is sometimes referred to as the

characteristic value problem. The scalar Xi and the vector vi are known, respectively, as the

eigenvalue and eigenvector of the ith mode of vibration of the structure. For a general

damped system, eigenvalues and eigenvectors are both complex. Note that Eq. (2.2.3) can

be easily rearranged in the more general mathematical eigenvalue problem format as

[ 0 Inxn i [ ]
M-'K M-'D [Y = ki yi (2.2.4)

where Inxn is the nxn identity matrix. Equation (2.2.4) is called a state space representation

of Eq. (2.2.3).

The eigenvalue and eigenvector can be related to some physical characteristic properties

of structures. For that matter, the ith eigenvalue is written as

Xi = -- ii + j (Ii -- (2.2.5)


where j = 1 Note that in writing this equation it is assumed that the system is

underdamped. The real scalar oi and ti are the natural frequency and damping ratio (or

damping factor), respectively, of the ith mode of the structure. The eigenvector yi indicates

the "shape" of vibration of the ith mode of the structure. The sets of frequencies, damping

ratios and mode shapes are sometimes referred to as modal parameters.

The symmetric nature of the property matrices (M, D, K) constrains the eigenvectors to

satisfy some cross-orthogonality relationships. For the purpose of discussing these







cross-orthogonality relationships, consider the following alternative state space
representation of the eigenvalue problem in Eq. (2.2.3)
[M O][VA VA ][A 0]+ [D K VA VA] = [ (2.2.6)
0 -K V V 0 ?A K 0 V V 0

where V = [ v n ]
A = diag(X1, X2, ..., Xn)
where the overbar denotes the matrix complex conjugate operator. Based on this particular
representation and the fact that the state matrices of Eq. (2.2.6) are symmetric, by proper
normalization of the eigenvectors, the cross-orthogonality relationship associated with the
system are given by

[VA VA M 0][VA = 2nx2n (2.2.7)
V V 0 K V V 2nx2n

VA VAN D K "VA VA] [A 0]
V V K 0 V V 0 X

Equations (2.2.7) and (2.2.8), respectively, clearly imply the following relations

VA] [M _][VA] =[0] (2.2.9)
V 0 K V 0

[VA]'[D K VA (2.2.10)
V K O V 0
T
where [ ] denotes the complex conjugate transpose operator ( [ ] = ). Note that,
contrary to Eqs. (2.2.7-8), no particular normalization of the eigenvector is needed in order
for Eqs. (2.2.9-10) to be satisfied.
Another state space representation of Eq. (2.2.3) is given by







0 M VA VA A 0 M 0] VA VA [0 (2
[M DN V V 0 ] 0 K[ V V 0jj
Based on the same argument discussed earlier, the cross-orthogonality conditions that arise
from this representation are

M ]T[ ][ ] = 2nx2n (2.2.12)
V v M Dr V Vx]

VA VA T[ M ][ VA VA] [A 0]
V V 0 K V V 0O

Again, the following two relationships follow,

VA][ [0 M][VA] = [0] (2.2.14)
V M D V 0

VAI -M ][VA] = [] (2.2.15)
V [ 0 K V 0

2.2.1 Undamped Models
In the modelling of structures, it is often assumed that the damping is negligible and
hence is set to zero. For these type of system models, the eigenvalue problem in matrix form
is given by
MVA2 + KV = [0] (2.2.16)
The matrices V and A are as defined earlier. For undamped systems, the eigenvalues X, are
pure imaginary and the eigenvectors yi are real. Note that the eigenvalues are related to the
system natural frequency by
X2 = co? (2.2.17)
Furthermore, by proper normalization of the eigenvectors, the cross-orthogonality relations
associated with this type of systems are








VTMV = Inxn (2.2.18)

VTKV = diag(wO, o, .2) (2.2.19)

2.2.2 Proportionally Damped Models

When the damping of the structure is accounted for, it is sometimes modelled to be

proportional,

D = aM + p3K (2.2.20)

where a and 03 are real scalars. The eigenvectors of a proportionally damped system are the

same as the eigenvectors associated with the corresponding undamped system. The

cross-orthogonality relationships associated with proportionally damped models are

VTMV = Inxn (2.2.21)

VTDV = diag(2(-o1, 22W2, 2ae0n) (2.2.22)

VTKV = diag(ao, o Wn) (2.2.23)

where all variables have the same definitions as in the previous discussion.

A more detailed development of the eigenvalue problem can be found in the books of

Meirovitch (1986, 1980) and Inman (1989).


2.3 Experimental Modal Analysis

The vibration characteristics of structures can also be measured experimentally. Recall

that the structural vibration characteristics are given by natural frequencies, damping ratios

and mode shapes (eigenvectors). The process of measuring the modal parameters is known

as modal testing or experimental modal analysis. An introductory treatment of the concept

of experimental modal analysis can be found in the book of Inman (1989); a more rigorous

coverage is treated in the book of Ewins (1986).

The hardware components needed in a modal analysis experiment are identified in

Figure 2.1 (Inman, 1989). A schematic of a simple modal vibration measurement test setup








is shown in Figure 2.2. Brief descriptions of some of the components and their functions are

given below.


T := Transducer
SC := Signal Conditioner


Figure 2.1 Components of a Vibration Measurement System
For Modal Analysis (Inman, 1989).



The exciter is used to deliver the driving force that puts the structure in motion. The

two most commonly used exciters are the shaker and the impulse hammer.


Transducers are devices that measure the driving force as well as the response of the

structure. They sense mechanical forces or motions and, then, convert them into

electrical signals. Mechanical forces are usually recorded by a force transducers.

Commonly, the vibration response of structures is measured by accelerometers.

Accelerometers are mounted directly on the structure and, naturally, record the

local accelerations.


Signal conditioners are used to match the signal, received from the transducers, to

the specification of the analyzer. These devices can also be used to amplify the

signals.















Accelerometer


Force
Transducer


Signal
Conditioner


Accelerometer


Signal
Conditioner


Signal Analyzer


Figure 2.2 A Simple Experimental Modal Analysis Setup.

The signal analyzer processes the electrical signal received from the signal

conditioners. The standard type of analyzer allows time domain signals to be

viewed in the frequency domain via a Fast Fourier Transform (FFT) algorithm. In

an FFT, the signals are first filtered, digitized and then transformed into discrete

frequency spectra. The frequency spectra can then be manipulated to compute the

modal properties of the structure.


It is important to note that the experimentally identified modal parameters are usually

affected by unpredictable measurement noise. Typically, natural frequencies are identified

to within 1 to 2% repeatability; damping ratios to within 5 to 15% repeatability, and mode

shapes to within 5 to 10% repeatability. In practice, due to testing limitations, the set of


Impulse
Hammer








structural modal parameters identified experimentally is incomplete with respect to the

analytical model. Experimental incompleteness is manifested in two forms, (i) a limited

number of measured modes of vibration, (ii) a limited number of eigenvector DOF

measurements.



2.4 Analytical and Experimental Model Dimensions Correlation


One major prerequisite common to most model refinement algorithms is to match the

number of degrees of freedom (DOFs) in the experimentally measured eigenvector and in

the discrete analytical model. Two different approaches have been commonly used to

resolve this matching problem when an incomplete set of eigenvector DOFs have been

measured. The first approach consists of reducing the dimension of the discrete analytical

model to the number of the measured DOFs. The other approach is to expand the measured

eigenvector to the size of the analytical model. A good evaluative survey of a number of

analytical model reduction techniques has been compiled by McGowan (1991). The subject

of eigenvector expansion is discussed in fair detail in the papers by Gysin (1990), and

Zimmerman and Kaouk (1992a). In the next two sections, some commonly used model

reduction and eigenvector expansion algorithms are reported and described.

2.4.1 Model Reduction Methods

In this section the general framework of the model reduction concept is first presented.

Then formulations of three commonly used model reduction techniques are summarized.

Mostly, the concept of model reduction has only been studied for undamped models. In this

presentation no attempt has been made to generalize these concepts to damped models.

For the purpose of reporting the general concept of model reduction, consider the

eigenvalue problem associated with an undamped model,


MVA2 + KV = [0]


(2.4.1)








where, as defined earlier, M and K are, respectively, the mass and stiffness matrix; V is the

eigenvector matrix and A is the diagonal eigenvalue matrix. Assume that only a subset of

the eigenvector DOFs has been experimentally measured.

Equation (2.4.1) can be reordered such that the DOFs associated with the measured

DOFs are in the upper rows of the equation,


MoVoA2 + KoVo = [0] (2.4.2)


Vm Kmm Kmu Mmm Mmu
where Vo = Vu Ko = Kum Kuu Mo Mum Mu


The matrices Mo, Ko, and Vo are, respectively, the reordered mass, stiffness and

eigenvector matrices. In the above equation, the subscripts "m" and "u" denote, respectively,

the components associated with the measured and unmeasured DOFs. A transformation

matrix P, that relates matrices Vu and Vm can be defined as

Vu = PVm (2.4.3)

A substitution of this relationship in the reordered eigenvector matrix Vo results in


Vo = p Vm = TVm (2.4.4)


Substituting Eq. (2.4.4) into Eq. (2.4.2) and premultiplying by TT yields the eigenvalue

problem of the reduced model,

MrVmA2 + KrVm = [0] (2.4.5)

where Mr = TTMoT

Kr = TTKoT

where Mr and Kr are the reduced mass and stiffness matrices, respectively. In terms of the

partitioned matrices, the reduced matrices are defined as








Mr = Mmm + pTMum + MmuP + pTMuuP (2.4.6)

Kr = Kmm + PTKum + KmuP + PTKuuP (2.4.7)

At this point, the only condition placed on P is the relationship of Eq. (2.4.3). Naturally,

matrix P can be computed directly from Eq. (2.4.3) if the eigenvectors of the system are

available. This approach is known as the exact reduction method and has been discussed in

the papers by Kammer (1987) and O'Callahan et al. (1989). The exact reduction method

requires solving for a large number of eigenvectors, which can be computationally

expensive. The three reduction methods that are presented in the forthcoming discussions

propose alternative techniques to compute matrix P. The first two do not require the

computation of the system eigenvalue problem. The last one requires the knowledge of one

eigenvalue which is computationally admissible.

2.4.1.1 Static Reduction

This reduction method is often referred to as Guyan (1965) reduction. In the static

reduction, the mass properties associated to the unmeasured DOFs are assumed negligible.

With that assumption, Eq. (2.4.2) can be written as

Mmm 0 Vm 2 Kmm Kmu Vm 0
[Mm O][VMA2 + [mm :] : = (2.4.8)
0 0 Vu Kum Kuu Vu 0

The second row of this matrix equation can then be manipulated as

Vu = K-uKumVm (2.4.9)

From comparing Eq. (2.4.9) to Eq. (2.4.3), it can deduced that the transformation

matrix, P, computed using the Guyan approach is given by

Pg = KualKum (2.4.10)

The reduced mass and stiffness matrices can then be computed by substituting matrix Pg

for matrix P in Eqs. (2.4.6) and (2.4.7). Naturally, the Guyan assumption (Eq. (2.4.8))








suggests that if the mass properties of the omitted DOFs are not small, the accuracy of the

Guyan reduced model could be lacking.

2.4.1.2 IRS Reduction

The improved reduction method (IRS) was formulated by O'Callahan (1989). It is an

improvement over the Guyan reduction in that it accounts for the mass properties of the

unmeasured DOFs. In the formulation of the IRS method, the Guyan reduced model is

corrected to include the mass influence of the unmeasured DOFs. This formulation is

somewhat lengthy and the interested reader is referred back to the paper of O'Callahan

(1989) or the thesis of McGowan (1991). The transformation matrix P computed using the

IRS reduction is

PIRS = Pg + KuU Mum + MuuPg MrKr, (2.4.11)

The reduced IRS model is then computed by substituting matrix PIRS in Eqs. (2.4.6-7).

2.4.1.3 Dynamic Reduction

The dynamic reduction was proposed as another improvement to the Guyan reduction

(Kidder, 1973; Miller, 1980; and Paz, 1984). This reduction utilizes the dynamic equation

associated with a single mode of vibration to compute the transformation matrix P. In this

technique, the transformation matrix P is arrived at by considering the reordered dynamic

equation associated with the ith mode,

[2Mmm + Kmm Mmu + Kmu m] (2.4.12)
I -(2.4.12)
?2Mum + Kum XMuu + Kuu JYu 0

where Xj is the ith eigenvalue; vmi and vYu are, respectively, the measured and unmeasured

eigenvector DOFs associated with the ith mode of vibration. Based on this partition, two

equations can be generated. By using the rows that correspond to the unmeasured DOFs

(lower rows), the following relationship is obtained,








Yu = [Muuux + Kuu] Mum + Kum] I, (2.4.13)

Hence, the transformation matrix associated with the ith mode is defined by

Pd, = [MuuX? + Kuu] [Mumk32 + Kum] (2.4.14)

The reduced mass and stiffness matrices are then computed by using Eqs. (2.4.6) and (2.4.7).
Note that different modes will result in different transformation matrices Pd and, hence,

different reduced mass and stiffness matrices.

2.4.2 Eigenvector Expansion Methods

Alternatively, the dimension of the measured eigenvectors can be correlated to the
dimension of the analytical model by using eigenvector expansion algorithms. The common
basis of these algorithms is the interpolation of the unmeasured eigenvector components. In
the forthcoming sections, two eigenvector expansion algorithms are reviewed.

2.4.2.1 Dynamic Expansion

The dynamic expansion technique (Berman and Nagy, 1983) is one of the most
commonly used eigenvector expansion algorithms. A slight modification of the Berman and
Nagy formulation is presented here to accommodate damped systems (Fuh et al., 1984). In
the formulation of the dynamic expansion, it is assumed that the measured modes satisfy the
eigenvalue problem involving the property matrices of the original model (M, D, K). For the
ith measured mode, this assumption takes the form

(k2iM + XeD + K) ve, = 0 (2.4.15)

where Xe. and ve. are, respectively, the ith experimental eigenvalue and eigenvector. The

matrices M, D and K have the same definitions as in the earlier sections. Assume that only
a subset of the DOFs of eigenvector ve. has been measured. Equation (2.4.15) can be

reordered, as in Section 2.4.1, such that the measured eigenvector DOFs reside in the upper
half of the equation,








Mmm Mmu] Dmm Dmu Kmm Kmu 1 e[ 0ei n211
Mum MuuJ e Dum Duu Kum Kuu LVXeuJ,


where Yemand Vu, are, respectively, the measured and unmeasured DOFs of ve.. The

subscripts "m" and "u" denote measured and unmeasured components. A rearrangement of
Eq. (2.4.16) yields


[Mmmn + keDmm + Kmm Mmu + ke,Dmu + Kmu em,i
I = (2.4.17)
L Mum + XeDum + Kum 2Muu + ke.Duu + Kuu jYeu,

From Eq. (2.4.17), two matrix equations in function of Vem, and veu, can be generated.

By using the equation associated with the second row of the partition, the unmeasured
components of the ith experimental eigenvector are computed to be


eu = [Muu + DXe + Ku] [MumX2J + Dumke, + Kum vemi (2.4.18)

Note that this expansion works on a single mode at a time. Also, notice that it involves
the original analytical model (M, D, K). This implies that the accuracy of the expansion
depends on the validity of the original analytical model.


2.4.2.2 Orthogonal Procrustes Expansion


Another expansion process that has shown great promise is the orthogonal Procrustes
(OP) expansion method presented in the papers of Smith and Beattie (1990) and Zimmerman
and Kaouk (1992a). The technique uses the general mathematical framework of the
orthogonal Procrustes problem (Golub and Van Loan, 1983).
Let Vem be the experimentally measured eigenvector component matrix and Vain be the
corresponding analytical eigenvector component matrix. The essence of the OP expansion is
to find a unitary transformation matrix Pop that closely rotates Vam into Vem. This is
attempted by solving the following problem,








Minimize 11 Vem VamPop IF
(2.4.19)
subject to pTpPop =

The solution to this general problem is discussed in the book of Golub and Van Loan and

is given by

Pop = YZT (2.4.20)

where Y and Z are, respectively the left and right singular matrices of matrix D defined by

(D = VTVem (2.4.21)

Let Vau be the eigenvector matrix associated with the unmeasured DOFs. In the

orthogonal Procrustes expansion, it is assumed that the transformation matrix computed

above also rotates Vau into the unmeasured "experimental" eigenvector component matrix,

Veu, as

Veu = VauPop (2.4.22)

Two different approaches have been defined to generate the expanded experimental

eigenvector matrix. The first is suggested in the paper by Smith and Beattie (1990),

Vam
Ve = Vau Pop (2.4.23)


In this definition, the expanded "experimental" eigenvector matrix is the rotated analytical

eigenvector matrix. The actual experimental measured eigenvector components Vem are

replaced by matrix VamPop. The advantage of the approach is that the resulting
"experimental" eigenvector satisfies the cross-orthogonality conditions (see Section 2.4).

The second viewpoint was proposed by Zimmerman and Kaouk (1992a),

[Vemn
Ve = [VauPopj (2.4.24)


Here, the unaltered eigenvector components measured experimentally are inserted in Ve. In

this viewpoint, if cross-orthogonalization of the expanded experimental eigenvectors is








required, a separated orthogonality algorithm can be used (Section 2.5). In the paper by

Smith at al. (1993), it is shown that for actual model refinement problems, both viewpoints

give equivalent eigenvector expansion results. However, for damage detection problems,

a preliminary study indicates that eigenvectors expanded by using the viewpoint defined in

Eq. (2.4.24) give better assessment of the damage.


2.5 Eigenvector Orthogonalization


Most matrix update algorithms require the measured eigenvectors to satisfy a

cross-orthogonality condition. This is especially true in matrix update algorithms in which

(i) the model of the structure is assumed undamped and the modelling errors are assumed to

be in only one of the two property matrices (M or K is assumed correct) or (ii) the system is

modeled by a proportionally damped model with errors in only two of the three property

matrices (M, D or K). In these situations, in order to insure symmetry of the updated property

matrices, it is required that the expanded experimental eigenvectors be orthogonal with

respect to the property matrix assumed accurate. This situation is encountered in the optimal

update algorithms developed by Baruch and Bar Itzhack (1978), Kabe (1985), Kammer

(1985), and Smith and Beattie (1991), Zimmerman and Kaouk (1992b), Kaouk and

Zimmerman (1993b) among others. In most of these algorithms, it is assumed that the mass

matrix is correct. This assumption is used in a number of the model refinement algorithms

since the inertial properties of structures are known to a good extent. In these cases, one

would expect the expanded experimental eigenvectors to be mass orthogonal. However,

because of measurements errors, this condition rarely occurs. For this reason, a great deal of

effort was focused on the development of mass orthogonalization techniques. Some of the

most notable work in that area was performed by Targoff (1976), Baruch and Bar Itzhack

(1978), and Baruch (1979). In the next two sections, two orthogonalization techniques

(Baruch and Bar Itzhack, 1978; Baruch, 1979) are discussed. Both techniques are mass

orthogonalization techniques; however, with obvious modifications, these techniques can be








adopted to solve the orthogonalization problem of the eigenvectors with respect to the

stiffness or the damping matrices.

2.5.1 Optimal Weighted Orthogonalization

The essence of the standard mass orthogonalization technique is to modify the measured

eigenvectors such that the mass cross-orthogonality condition is satisfied. Baruch and Bar

Itzhack (1978) proposed an elegant solution to that problem. An overview of their problem

statement and solution is given below.

Assuming that Ve is a matrix of expanded experimental eigenvectors that need to be

mass orthogonalized. The present formulation searches for the optimal mass weighted

change of matrix Ve such that the mass cross-orthogonality condition is satisfied. This

problem is cast as

Minimize N (Veo Ve) 1F (2.5.1)

subject to VT M Veo = I (2.5.2)

where N=M112 and M is the mass matrix. By means of a Lagrange multiplier, Eq. (2.5.2)

can be incorporated into Eq. (2.5.1); then the application of the optimality conditions yields

the following expression for Veo,

Veo= Ve(V M Ve)-1/2 (2.5.3)

Before being incorporated into the orthogonalization process, the measured

eigenvectors have to be unit mass normalized, i.e.,

Vej = ve,(vT M ve 12 (2.5.4)

where ve. is the it expanded experimental eigenvector (ith column of Ve).

2.5.2 Selective Optimal Orthogonalization

Some structures exhibit rigid body modes (modes with zero eigenvalues). It is desirable

to preserve these rigid body modes in the refinement process. However, some matrix update








algorithms require the rigid body modes and the experimental eigenvectors to be mass

orthogonal to insure symmetry of the updated property matrices (see Chapter 6). Naturally,

the rigid body modes will be corrupted if they are incorporated along with the expanded

experimental eigenvectors in the above orthogonalization process. Thus, Baruch (1979)

presented a modification of the procedure in Section 2.5.1 to deal with such a problem. The

resulting problem is a selective orthogonalization and is formulated as follows,

minimize 11 N (Veo Ve) 1F (2.5.5)


Subject to Veo M Veo = I (2.5.6)

and Veo M Vr = [0] (2.5.7)

In the above equations, Ve and Vr are, respectively, the expanded experimental eigenvector

matrix and the rigid body mode matrix. Again, the Lagrange multiplier is used, and the

orthogonalized experimental eigenvector matrix that satisfies the conditions in Eqs. (2.5.6)

and (2.5.7) is found to be


Veo= Q(QT M Q)-/2
(2.5.8)
where Q = Ve VrVTMVe

Note that, as in previous process, the expanded experimental eigenvectors have to be

unit mass normalized.


2.6 Load Path Preservation

Many matrix update algorithms introduce additional load paths in their updated models,

i.e., elements of the mass, damping or stiffness matrices that were originally zero may

become nonzero. Whether or not preserving the original load path is a practical problem is

still a matter of current debate. It seems that for damage assessment of truss structures it is

desired to maintain load paths. In the paper of Zimmerman and Kaouk (1992a), an iterative

approach to preserve the load path of the original property matrices was developed. The








approach was presented in the context of the symmetric eigenstructure assignment model

refinement algorithm (discussed in Chapter 4); however, its application can also be extended

to other model refinement algorithms. In Figure 2.3, a flow chart of the iterative load

preservation algorithm is presented. The procedure is illustrated for a general model

refinement scenario in which all three property matrices (M, D, K) are being updated.

However, it can be easily modified to accommodate other refinement problems.


Figure 2.3 Flow Chart of the Iterative Load Path Preservation Algorithm



In the flow chart, the matrices Ma,m, Da,m, and Ka,m are respectively the adjusted

masked mass, damping and stiffness matrices defined by








Ma,m = Ma 0 Mm

Da,m = Da 0 Dm (2.6.2)

Ka,m = Ka 0 Km

where Ma, Da, and Ka are the adjusted mass damping and stiffness matrices. The matrices

Mm, Dm, and Km are the masking matrices associated with the original mass, damping and

stiffness matrix. By definition, the masking matrix, Am, associated with matrix A is given

by

Am(i,j) = 1 if A(i,j) 0
(2.6.1)
Am(i,j) = 0 if A(i,j) = 0

In Eqs. (2.6.2), the operator 0 is the element-by-element (scalar) matrix multiplication. Let

B and C be two nxn matrices, then the element-by-element multiplication of B and C is given

by

S = B 0 C => S(i,j) = B(i,j) C(i,j) i,j = l,..., n (2.6.3)

At every iteration, the norms of the matrix differences between corresponding adjusted

and adjusted masked property matrices are computed. At a given iteration, if the three

computed norms are equal to zero or within user set limits, then the load paths of the original

three property matrices have been exactly achieved or achieved within user state guidelines.

Thus, the procedure is halted, and the refined model consists of the adjusted property

matrices computed at that particular iteration. It should be noted that there is no formal

guarantee of convergence in using this iterative procedure. Experience gained in using the

present algorithm indicates that if the experimental modal data are consistent with the

sparsity pattern, the procedure will converge. Consistent data means that there exist mass,

damping and stiffness matrices that have the same sparsity pattern as the original matrices

and also exhibit the measured test data. Otherwise, if the data are inconsistent, the original

sparsity patterns will not be exactly preserved. In this case, the added load path terms of Ma,

Da and Ka which should be zero will be closer to zero after application of the algorithm.















CHAPTER 3
INVERSE / HYBRID APPROACH FOR FINITE
ELEMENT MODEL REFINEMENT




3.1 Introduction

The inverse eigenvalue problem is concerned with the construction of the property

matrices (mass, damping or stiffness) of a dynamic model using experimentally measured

modal data. These techniques require complete modal properties. Thus, for an n-DOF

model, n natural frequencies, damping ratios and mode shapes (eigenvectors) must be

measured, and the identified mode shapes must be of dimension n. Due to practical testing

limitations, this is rarely accomplished for typical large structures. In this chapter, the

application of the inverse problem is extended to model updating by combining

experimental measurements and original analytical FEM modal information. Again,

refinement implies correlating the measured and analytical modal properties.


3.2 Theoretical Formulation

The dynamic structure under consideration is assumed to be successfully modelled by

an n-DOF proportionally damped nongyroscopic and noncirculatory (symmetric property

matrices) FEM. The free vibration motion of such a dynamic structure can be analytically

represented by a differential equation of the form

Mw(t) + Dw(t) + Kw(t) = 0 (3.2.1)

where the variables M, D, and K are nxn real symmetric matrix models of the mass, damping

and stiffness properties of the structure. The nxl time varying vector w(t) represents the n








displacements of the n-DOF model of the system. The overdots represent differentiation

with respect to time. The eigenvalue problem associated with the differential equation

shown in Eq. (3.2.1) is of the form

Mvy. X + Dv.i X + Kvy = 0 (3.2.2)

where Xi and yi are, respectively, the eigenvalue and eigenvector of the ith mode of vibration.

In this problem, it is assumed that the accuracy of the original FEM is lacking and, hence,

needs improvement. Furthermore, it is assumed that the inaccuracy of the original FEM is

solely due to modeling errors in the stiffness and damping properties.

The model refinement, proposed herein, exploits the cross-orthogonality relations that

arise from the symmetric nature of the property matrices and the proportional damping

assumption. As discussed in Chapter 2, by proper normalization of the eigenvectors these

cross-orthogonality relations have the form

VT M V = Inxn (3.2.3a)

VT D V = diag(2tlow, ,2tn n) = 5: (3.2.3b)

VTKV = diag(01 2, ,on2) = Q (3.2.3c)


V = [ ... V]

where wi and i are the natural frequency and damping ratio, respectively, of the ith mode

of the structure. The matrix Inxn is the nxn identity matrix. It is important to recognize that

Eq. (3.2.3a) represents necessary and sufficient conditions for conserving symmetry and

damping proportionality when updating the stiffness and damping properties of

proportionally damped systems.

Suppose that p (p << n) modes of an existing structure have been experimentally

identified (mode shapes or eigenvectors, frequencies and damping ratios). Assume that the

dimension of measured eigenvectors is equal to the dimension of the FEM, i.e., all n

components of the measured eigenvectors are available. It is widely accepted that in the








absence of specific experimental measurements a good approximation to the unmeasured

modes is their corresponding analytical modal information. With that in mind, a complete

hybrid set of modal data is generated by combining experimental and analytical information

Vea = [Ve Va]

[Qe 0
0 Qa (3.2.4)

le 0
lea = 0 a



where V is the eigenvector matrix; 2 and I are diagonal matrices of frequencies squared and

damping ratios, respectively. The subscripts e and a denote, respectively, experimental and

analytical sets. At this point, the complete "hybrid" set of eigenvectors, Vea, does not satisfy

the cross-orthogonality conditions defined in Eq. (3.2.3a); thus the conditions in Eqs.

(3.2.3b) and (3.2.3c) are not met. One possible solution to this problem is to modify all of

Vea in an optimal way to comply with the orthogonality requirement. This approach treats

all parts of Vea equally, and thus overlooks the fact that the experimental modes are known

with a higher confidence. Naturally, a technique that assigns a higher credibility to the

experimental eigenvectors, Ve, is preferable. This can be achieved by incorporating Vea into

the orthogonalization process group by group in the order of their descending credibility

(experimental then analytical). If the problem is set such that the experimental modes, Ve,

are corrected first, it is clear that the analytical modes, Va, will be subject to larger correction

when incorporated into the orthogonalization process since they will be subject to more

constraints.

The experimental modes, Ve, are orthogonalized by using the orthogonalization

technique formulated by Baruch and Bar Itzhack (1978). The general formulation of this

orthogonalization technique is discussed in Chapter 2. For this particular application, the

problem consists of finding the matrix ,Veo, that satisfies the cross-orthogonality condition,








Vo M Veo = Ipxp (3.2.5)

and that minimizes the weighted Euclidean norm,

(D = |N (Veo Ve) I
(3.2.6)
where N = M1/2

The solution to this problem, as reported in Chapter 2, is


Veo= Ve(V M Ve)-1/2 (3.2.7)

The next step is to invoke the orthogonality requirement on the analytical eigenvector

matrix, Va, by searching for a matrix Vao that satisfies the following two conditions,

Vao M Vao = I(n-p)x(n-p) (3.2.8)


and Vao M Veo = [0] (3.2.9)

while minimizing the objective function,

F =1 N (Vao Va) IIF
(3.2.10)
where N = M1/2


A similar problem was also treated by Baruch (1980) in a different context. A brief

discussion of the solution approach can be found in Chapter 2. The set of eigenvectors, Vao

that satisfies this problem is given by


Vao = Q(QT M Q)-/2
(3.2.11)
Q = Va VeoV oMVa


Clearly, the resultant matrix, Veao = [ Veo Vao ], satisfies Eq. (3.2.3a). The corrected
stiffness and damping matrices are then computed using Eqs. (3.2.2b) and (3.2.2c).








Ka = MVeao Q ea VTaoM (3.2.12a)

Da = MVeao E eaVeLoM (3.2.12a)

where Oea = 2Ileagea


Note that the matrices Ka and Da computed from Eqs. (3.2.12a,b) will be symmetric.

The above formulation suggests that the system modelled by the original mass matrix (M)

and the stiffness (Ka) and damping (Da) matrices computed from Eqs. (3.2.12a,b) will have

eigenvectors Veao, frequencies 4ea, and damping ratios Hea. Some structures exhibit rigid

body modes of vibration. Commonly, it is desirable to preserve these rigid body modes in the

updated model. The above formulation also suggest that the updated model will preserve the

original rigid body modes.

The procedure developed above can be easily contracted to address the case when the

system model does not account for the effects of damping (undamped model). The

contraction can be obtained by setting to zeros matrices D, Oea, and FIea in Eqs. (3.2.2-4).

The computational burden limits the size of the FEM which can be updated (order of

200). Essentially, the limiting factor is that all mode shapes of the structure that are not

available from experimental measurements must be calculated analytically by solving the

eigenvalue / eigenvector problem.


3.3 Numerical Illustration

The system addressed in this investigation is the commonly used eight degrees of

freedom model shown in Figure 3.1. This model was developed by Kabe (1985) to give a

common testbed for the evaluation of the performance of model refinement algorithms. An

original undamped analytical model of the system was generated by using the mass and

stiffness properties shown in Figure 3.1. The elements of the original analytical stiffness

matrix are displayed in the second column of Table 3.1. In this problem, it is assumed that the

original stiffness matrix of the model is incorrect. These inaccuracies were simulated by

























mi =0.001 m8= 0.002 mj=1.0 j=2,..., 7

k1 = 1000 k2 = 10 k3 = 900 k4 = 100 k5 = 1.5 k6 = 2.0


Figure 3.1 Kabe's Problem: Analytical Test Structure.


using incorrect stiffness constants for most of the springs. The elements of the actual correct

stiffness matrix are given in the last column of Table 3.1. Note that the present model

refinement problem is very challenging and because of the large difference between the

stiffness matrices of the original and exact model. In this study, two cases will be considered.

In the first case, it is assumed that only the modal parameters (eigenvalue and eigenvector) of

the first mode were measured. In the other, the modal parameters of first three modes are

assumed to be available. In both cases, the measured eigenvectors are supposed to be full (all

degrees of freedom of the eigenvector(s) are measured).

The main objective of this investigation is to compare the performance of the

inverse/hybrid method to the algorithm proposed by Baruch and Bar Itzhack (1978). The

Baruch and Bar Itzhack model update technique is one of the most commonly used model

refinement algorithm. One of the main reasons for its common use is because the algorithm

is computationally efficient. The updated stiffness matrices generated by using the Baruch

and Bar Itzhack approach for the one mode and three mode cases are shown in the third and








fourth column of Table 3.1, respectively. The fifth and sixth columns of Table 3.1 display the

results of using the inverse/hybrid algorithm for the one mode and three modes cases. For

both cases, it is clear that the performances of both algorithms is lacking in predicting the

exact stiffness matrix. As expected, It can also be seen that both algorithms generate better

results when using three measured modes. A comparison of the results generated using the

Baruch and Bar Itzhack approach and the inverse/hybrid method shows that both algorithms

give the same type of results. This similarity of results was also encountered in other

example problems.


3.4 Summary

A model refinement approach that uses a hybrid set of experimental and analytical

modal properties was formulated. The developed approach, termed the inverse/hybrid

algorithm, was illustrated by using a computer simulated example. Part of the evaluation of

the proposed algorithm was the comparison of its performance with the performance of the

Baruch and Bar Itzhack algorithm. It was found that both algorithms give the same type of

results. However, the computations involved in the inverse/hybrid approach exceed those

involved in the Baruch and Bar Itzhack algorithm. Essentially, the inverse/hybrid approach

requires the computation of all eigenvalues and eigenvectors of the structure that are not

available from experimental measurements. As will be illustrated in the forthcoming

chapter, the Baruch type approaches are not suited for damage assessment applications. For

these reasons, further development of the Inverse/Hybrid algorithm was not investigated and

new formulations (Chapters 4 & 6) were considered.








Table 3.1 Kabe's Problem: Elemental Stiffness Components.


Baruch Inverse/Hybrid
Element # Original Mode 1 Modes 1-3 Mode 1 Modes 1-3 Exact
(1,1) 2.0 2.0 2.0 2.0 2.0 1.5
(1,2) -2.0 -2.0 -3.0 -2.0 -3.0 -1.5
(1,3) 0.0 0.0 -0.1 0.0 -0.1 0.0
(1,4) 0.0 0.0 -0.1 0.0 -0.1 0.0
(1,5) 0.0 0.0 -0.1 0.0 0.1 0.0
(1,6) 0.0 0.0 0.1 0.0 0.1 0.0
(1,7) 0.0 0.0 0.1 0.0 0.0 0.0
(1,8) 0.0 0.0 0.0 0.0 0.0 0.0
(2,2) 1512.0 1508.6 1024.2 1510.6 1024.3 1011.5
(2,3) -10.0 -31.5 -68.5 -21.5 -71.4 -10.0
(2,4) 0.0 -8.9 -9.0 -15.0 -8.4 0.0
(2,5) 0.0 -8.9 -20.9 -15.0 23.8 0.0
(2,6) 0.0 -21.6 38.5 -11.5 35.2 0.0
(2,7) 0.0 -3.9 -9.1 -1.6 -8.7 0.0
(2,8) 0.0 -0.1 0.1 -0.1 0.1 0.0
(3,3) 1710.0 1574.6 1560.8 1624.1 1612.2 1110.0
(3,4) 0.0 -44.9 -49.8 -75.6 -76.2 0.0
(3,5) -200.0 -244.8 -244.1 -275.6 -276.9 -100.0
(3,6) 0.0 -136.5 -123.4 -86.2 -74.5 0.0
(3,7) 0.0 -24.6 50.8 -13.1 46.8 0.0
(3,8) 0.0 -0.4 -0.1 -0.4 0.0 0.0
(4,4) 850.0 1083.1 1087.9 1102.1 1099.9 1100.0
(4,5) -200.0 3.27 25.1 51.8 48.0 -100.0
(4,6) -200.0 -254.3 -242.2 -276.5 -274.9 -100.0
(4,7) 0.0 -10.1 27.0 -17.2 30.8 0.0
(4,8) 0.0 -0.2 -0.1 -0.2 -0.2 0.0
(5,5) 850.0 1082.2 1089.5 1101.4 1101.1 1100.0
(5,6) 0.0 -45.2 -49.5 -76.5 -75.1 0.0
(5,7) 0.0 -10.1 -12.4 -17.2 -11.7 0.0
(5,8) 0.0 -0.2 -0.1 -0.2 -0.2 0.0
(6,6) 1714.0 1576.4 1565.0 1617.5 1610.8 1112.0
(6,7) -10.0 -34.8 -82.6 -23.1 -76.2 -10.0
(6,8) -4.0 -4.4 -4.7 -4.4 -4.5 -2.0
(7,7) 1512.0 1507.5 1027.9 1510.1 1028.0 1011.5
(7,8) -2.0 2.1 -4.1 -2.1 -4.2 -1.5
(8,8) 6.0 6.0 6.0 6.0 6.0 3.5














CHAPTER 4
SYMMETRIC EIGENSTRUCTURE ASSIGNMENT MODEL
REFINEMENT ALGORITHM




4.1 Introduction

Eigenstructure assignment is a control concept used to alter the transient response of

linear systems. This is done by forcing the system to have some predetermined eigenvalues

and eigenvectors. A detailed overview of eigenstructure assignment theories can be found in

the paper by Andry et al. (1983). Inman and Minas (1990), Zimmerman and Widengren

(1989, 1990), and Widengren (1989) have developed model refinement algorithms based on

the mathematical framework of eigenstructure assignment. The basic idea of these model

refinement techniques is to design the pseudo-controller which is required to produce the

measured modal properties (natural frequencies, damping ratios and mode shapes) with the

original finite element model (FEM) of the structure. The pseudo-controller is then

translated into matrix adjustments applied to the initial FEM.

In this work, the eigenstructure assignment based model refinement algorithm proposed

by Zimmerman and Widengren (1989, 1990) is extended to better approach the damage

assessment problem. A subspace rotation algorithm is developed to enhance eigenvector

assignability. Finally, the enhanced algorithm is tested and compared to other techniques on

both "simulated" and actual experimental data.


4.2 Problem Formulation

In this section, a review of the Zimmerman and Widengren (1989, 1990) refinement

technique, which is termed the symmetric eigenstructure assignment model refinement








algorithm (SEAMRA), is presented. This review is essential in order to properly introduce

and discuss the extension and improvement proposed in this work.

4.2.1 Standard Eigenstructure Assignment Formulation

Consider the standard differential equation of motion of an n degrees of freedom

damped, nongyroscopic and noncirculatory structure with control feedback,

Mw(t) + Dw(t) + Kw(t) = Bou(t) (4.2.1)

Again, M, D, and K are n-by-n real symmetric matrix models of the mass, damping and

stiffness properties of the structure. Assume that these matrices were generated using the

finite element method. The nxl time varying vector, w(t), represents the n displacements

of the n-DOF FEM of the system. The overdots represent differentiation with respect to

time. In control terminology, BO is the nxm (m << n) control influence matrix describing

the actuator force distributions and u(t) is the mx 1 vector of output feedback control forces

defined by

u(t) = Fy(t) (4.2.2)

In Eq. (4.2.2), F is the mxr feedback gain matrix and y(t) is the rxl output of sensor

measurements defined by

y(t) = Cow + C1w (4.2.3)

in which Co and C1 are the rxn output influence matrices corresponding to position and

velocity, respectively. A substitution of Eqs. (4.2.2) and (4.2.3) into Eq. (4.2.1) yields

Mw(t) + (D BoFC1)w(t) + (K BoFCo)w(t) = 0 (4.2.4)

It is clear, from Eq. (4.2.4), that the feedback controller results in residual changes, BoFCo

and BoFC1, to the stiffness and damping matrices, respectively. These changes can be

viewed as perturbations to the initial finite element model (FEM) such that the adjusted FEM

matches closely the experimentally measured modal properties. The adjusted FEM consists

of the original mass matrix and the adjusted stiffness and damping matrices given by







Ka K
Da = D


- BoFCo
- BoFC1


(4.2.5)


Assume that modal analysis of the structure under consideration has been performed and that
p modes (p eigenvalues Xei, and p eigenvectors Ve) have been identified. As discussed earlier
in Chapter 2, in practice p is typically much less than n. The feedback gain matrix F, such
that the adjusted FEM eigendata matches the experimental modal parameters, is computed
using standard eigenstructure assignment theories (Andry et al. 1983):


F = [Z A f[C C -1


where


(4.2.6)


A0 Inxn
A M-1K -M-iD


B = M-1K


T = [B P]


A = T-'AT = A]



V = T-1 W W
WA WA


STB = Imxm
0


A = diag(Xe, ,e2, ., ep)


w = [eaivea,, .-,Veap]


Z = S-1[ W A [ A
WA WA


The overbar in the above equations indicates the complex conjugate operator. The vectors

Yea, in matrix W are the expanded "best achievable" eigenvectors associated with the
experimentally measured eigenvectors ve,. An explanation of the concept of "best


T- S
T- =1
S2








achievable" eigenvectors is discussed in Section 4.2.3. The submatrix P of matrix T is

arbitrary as long as T is invertible.

At this point, the variables BO, Co and C1 are still arbitrary. A random selection of these

variables will usually result in nonsymmetric perturbation matrices and, consequently

nonsymmetric adjusted stiffness and damping matrices. This clearly conflicts with the

fundamental symmetry requirement of most structures' FEM. In the formulation of Inman

and Minas (1990), the resulting perturbation matrices from the pseudo-controller are forced

to be symmetric through a nonlinear unconstrained optimization problem. Zimmerman and

Widengren (1989, 1990) proposed a non-iterative and computationally more efficient

approach to satisfy the symmetry requirement. This approach in discussed in the following

section.

4.2.2 Symmetric Eigenstructure Assignment Formulation

The perturbation matrices are symmetric if the following conditions are met,

BoFCo = CF TBo (a)
(4.2.7)
BoFC1 = CF TB (b)

At this point, two additional assumptions are made. As a prerequisite to the existence of the

inverse of some matrices used in the computations, it is assumed that the number of pseudo

sensors and actuators is equal to twice the number of measured modes (m=r=2p). The other

assumption consists of restricting the matrices Co and C1 by the conditions

CO = GoBo (a)
0 (4.2.8)
C1 = GIBT (b)

where Go and G1 are mxm invertible matrices. A substitution of Eq. (4.2.8) into Eq. (4.2.7)

simplifies the symmetry conditions to the following relationships,

FGo = GTFT (a)
(4.2.9)
FGi = G TFT (b)








By using the conditions in Eq. (4.2.9), along with the expression for the feedback gain matrix

(Eq. (4.2.6)), a necessary but not sufficient condition on Go and GI, for symmetric

perturbation matrices, is expressed in the form of a generalized algebraic Riccati equation,

A1X + XA2 + XA3X + A4 = [0] (4.2.10)

where

X-GG'G0

A, = ) T]

A2 =*- 1- 1 *


A3 = *-1-1(a ** a )~-*

A4 = o*a-10-la* Imxm


[w*BO' "AW*BO]
T= a = a= Z- AIV
W Bo AW Bo



The matrices Al, W, Z, and V are defined in Eqs. (4.2.6). The superscript ()-* indicates the

inverse of the complex conjugate transpose matrix. Equation (4.2.10) can be solved for X

by using the techniques described in the papers of Potter (1966) or Martensson (1971). In

general, there exist multiple solutions (X's) to this generalized algebraic Riccati equation.

With all solutions computed, the next step is to decompose these solutions into Go's and GI's.

It is shown in the paper by Zimmerman and Widengren (1989) that for a given solution X,

any selection of Go and GI satisfying X = G- 1Go results in the same adjusted damping

(Da) and stiffness (Ka) matrices. Hence, either GI (or Go) can be chosen arbitrarily, as long

as its inverse exists. Then, Go (or GI) is calculated from the relationship X = GI 'G0.

For each set (Go, G1), a feedback gain matrix F is calculated from Eq. (4.2.6), and the

corresponding adjusted damping (Da) and stiffness (Ka) matrices are computed using Eqs.








(4.2.5). At this point, a rationale is proposed to choose the most meaningful adjusted

damping and stiffness matrices. Among all computed sets (Da, Ka), it is apparent that only

the ones that are real and symmetric are acceptable. When dealing with a model refinement

problem, among all acceptable solutions, the final selection could be made by choosing the

set (Da, Ka) that minimizes the cost function,

J = q |I K Ka IIF + I| D Da IIF (4.2.11)

where q K IF



Clearly, this process selects the set (Da, Ka) that results in a minimum change from the

original set (D, K). The scale factor q in Eq. (4.2.11) is used to give equal weight to the

changes in D and K. For the damage detection problem, there is no unique rationale to

choose the "best" set (Da, Ka). A physically intuitive approach is to use engineering

judgement in selecting the "best" updated model. Thus, all acceptable "adjusted" sets of

solutions should be inspected to determine which best provides information concerning the

state of damage.


4.2.3 Best Achievable Eigenvectors


From standard eigenstructure assignment theory (Andry, et al. 1983), it is shown that the

measured eigenvectors are not always exactly assignable to the adjusted finite element

model. In fact, it can be shown that the measured eigenvectors are assigned exactly if and

only if they lie in their respective achievable subspace. The achievable subspace associated

to the ith mode is defined by

L, = (MX? + DXi + K)-B0 (4.2.12)


where 1i is the measured eigenvalue of the ith mode. When all n components of the

experimental eigenvectors are available, the ith best achievable eigenvectors is defined as the








least square projection of the ith experimental eigenvector Ve, on the ith achievable subspace

Li. This projection is schematically illustrated in Figure 4.1.




Achievable Subspace








Figure 4.1 Best Achievable Eigenvector Projection.



This best achievable eigenvector is given by


yea = LL L] L ve (4.2.13)


When only a subset s of the eigenvector components are measured, s < n, the least square

projection discussed above can be used to simultaneously expand and project the measured

eigenvectors. In this case, the ith expanded best achievable experimental eigenvector is

given by


Yea, = Li i Li LiVe (4.2.14)


where Li are the rows of Li which correspond to the measured eigenvector components.

Notice that the calculation of the p achievable subspace using Eq. (4.2.13) requires p

inversions of an nxn matrix. Although the matrix to be inverted is typically banded, this

may present a practical computational burden when dealing with large FEMs. The next

section discusses an approach that does not require the actual computation of the achievable

subspaces and hence avoids this computational burden.








4.2.4 Selection of Bu : The Subspace Rotation Method

So far, the control influence matrix has not yet been completely defined. The preceding

formulation suggest that different Bo may possibly result in different adjusted FEM. Hence,

it is essential to develop a physically meaningful rationale to select B0.

Zimmerman and Widengren (1989, 1990) proposed an approach, termed the mode

selection method, that consists in selecting Bo such that the unmeasured modes of the

structure are nearly unchanged. In other word, BO is selected such that only the measured

modes of the structure are corrected. This selection technique fixes the achievable subspaces

in which the eigenvectors must lie, and hence places a limitation on the assignment process.

In most studied cases, the experimental eigenvectors were not assigned exactly since their

assignment "success" depends on the locations vis-a-vis the achievable subspaces set by the

selection of B0.

In this work, a new method of selecting BO, termed the subspace rotation method, is

proposed. The subspace rotation method is based on selecting BO such that the measured

eigenvectors lie exactly in the achievable eigenvectors subspaces. This procedure is

illustrated in Figure 4.2 and is accomplished by setting BO as

B0 = [br, br2, ... ,brp I bi, bi, ... bi] (4.2.15)

where br = real [(MX2 + DXk + K)ve]


bi = imaginary [(MX2 + DX, + K)v

where ve. is the eigenvector associated with the jth experimentally measured mode, and it

is assumed that all n components of the experimental eigenvectors are available. This could

be accomplished by any of the procedures discussed in Chapter 2. Clearly, when Bo is

selected as shown in Eq. (4.2.15), the measured expanded eigenvectors lie exactly in the

achievable subspaces defined in Eq. (4.2.12). Hence, there is no need for the projection

operations defined in Eq. (4.2.13). This eliminates the required p inverses of nxn matrices








involved in computing the achievable subspaces. As will be seen in Chapters 5 and 6, the

elements of BO, as defined by Eq. (4.2.15), give an indication to the pseudo-controller about

the extent of modification of each DOF in order for the structure to exhibit the jth measured

eigenvalue and eigenvector.






Rotated Subspace --


Nominal Subspace







Figure 4.2 Rotation of the Achievable Subspace.





4.3 Numerical Illustrations



In this section, the characteristics of the proposed enhancement to the symmetric

eigenstructure assignment model refinement algorithm (SEAMRA) are evaluated and

compared to other refinement techniques for two example problems. The first problem is a

widely-used spring-mass computer simulated example (Kabe, 1985). It is used here for the

purposes of illustrating model refinement for a large local discrepancy, analogous to a

damage detection situation. The phenomena of global/local mode switching and load path

preservation are examined in this problem. The second problem is used to illustrate the

characteristics of the enhanced SEAMRA in updating the finite element model (FEM) of a

laboratory cantilever beam using actual measured modal parameters.








4.3.1 Damage Detection: Kabe's Problem


Kabe's eight degree of freedom spring-mass system is shown in Figure 4.3. The mass

and stiffness properties of the system are included in the figure. This problem presents a

challenging situation for damage detection in that stiffness values of various magnitudes are

included. The model exhibits closely-spaced frequencies and both local and global modes of

vibration.


m =0.001 m = 0.002 m = 1.0

ki = 1000 k2 = 10 k3 = 900 k4 = 100


j =12,.. 7


Figure 4.3 Kabe's Problem.




A variation of Kabe's original problem is used here. Rather than the standard initial

model commonly used, which has incorrect values for all of the connecting springs, only a

single spring constant is changed. This is reflective of the fact that damage may occur as a

large local change in the stiffness of a structural member.








4.3.1.1 Local to Global Mode Change


In the first problem, Kabe's initial model is only incorrect for the spring between

masses 3 and 5. A value of 500, five times that of the exact spring, is assumed in this problem.

Changing the spring value from 500 to 100 also causes a local mode of vibration to be

replaced by a global mode, thus presenting a difficult challenge for damage detection.

Figure 4.4 presents element-by-element stiffness matrix results for applying the Baruch

and Bar Itzhack update (1978) and the symmetric eigenstructure assignment model

refinement algorithm. Baruch Damage indicates that the update was made using Baruch and

Bar Itzhack's algorithm. SEA-M Damage indicates that the update was made using the

SEAMRA with BO selected by using the modal (M) selection method. SEA-SR Damage

indicates that the update was made by using the SEAMRA with BO selected using Subspace

Rotation (SR) method. The x-coordinate on all plots are the indices of a column vector

constructed by storing the upper triangular portion of the stiffness matrix in a column vector.

The y-coordinate on each plots consists of the difference between the updated stiffness

matrix elements and the original stiffness matrix.

In the first case, as shown in Figure 4.4, it is assumed that only the fundamental mode of

vibration is measured, but all eigenvector components have been measured. Thus, no

expansion of eigenvectors is required. It is evident from Figure 4.4 that the Baruch update is

unable to discern the damage, but that both the SEA-M and SEA-SR are able to clearly locate

the damage. In fact, the SEA-SR was able to exactly reproduce the correct stiffness matrix.

This was true independent of which mode was used in the update. Also, it should be noted

that the Baruch update tends to focus elemental changes in the third and fifth row of the

stiffness matrix, indicating the possibility of damage between these degrees of freedom, but

certainly giving no clear indication to the extent of damage. As is evident from the plot, the

Baruch update has spread errors over several elements. Using the algorithm of Lin (1990),

the damage vector is given as a = [1.0 0.93 0.72 0.83 0.70 0.90 0.97 1.0]T, where the

element number corresponds to the structural DOFs and a number less than 1 indicates the



















-500


Indices


Baruch Damage


JUULI I I i


0 10 20
Indices


0 10 20 30
Indices


0 10 20 30
Indices


Figure 4.4 Results for Kabe's Problem using the 1st Mode, Full Eigenvector.



possibility of damage affecting that DOF It is obvious that DOFs 3 and 5 are affected by

damage, but the results also indicate strong damage of DOF 4.

In the second case, as shown in Figure 4.5, it is assumed that the first three modes of

vibration have been measured, but only the first three components of the eigenvectors have

been measured. The eigenvectors components were expanded for the Baruch update using

dynamic expansion (Berman and Nagy, 1983) with subsequent orthogonalization (Baruch

and Bar Itzhack, 1978). The least squares expansion was used for the SEA-M update. The

SEA-SR update utilized the orthogonal Procrustes expansion (Chapter 2). In comparing

Figure 4.5 to Figure 4.4, it is clear that the damage detection capabilities of all three

algorithms have been degraded when using expanded mode shapes, even though more

modes have been measured. However, both the SEA-M and SEA-SR updates give a clear


enr\




















-500W '
0 10 20
Indices

500 __ SEA-M Dan
5001- 1


-500
0


Baruch Damage


3UUr j I "


0-


I I I


10 20 30
Indices


0


10 20
Indices

SEA-SR Dar


0 10 20 30
Indices


Figure 4.5 Results for Kabe's Problem using Modes 1, 2, 3,
and Eigenvectors Components 1, 2, 3.




indication to both the location and extent of damage. Using Lin's algorithm, the damage

vector is given as a = [1.0 0.81 0.75 0.83 0.82 0.79 0.85 1.0]1. It is difficult from

inspection of a to determine the location of damage.

The effect of applying the iterative load path algorithm described in Chapter 2 in the

update procedure is shown in Figure 4.6. For the Baruch update, 100 iterations were

performed. For the SEA-M and SEA-SR updates 2 and 3 iterations respectively, were

performed. The iterations were halted early for both SEA updates because the discrepancy

between the eigenstructure before and after masking was within the numerical precision of

the symmetric eigenstructure assignment software. It is seen that the load path enforcement

further enhances the damage detection capability of both SEA updates.


r f


0n


-









_500 Actual Damage




0----- ---,----


-500
0


10 20 30
Indices


0d


0-


-500 L
0


0 10 20 30
Indices


Baruch Damage


10 20
Indices


20
Indices


Figure 4.6 Results for Kabe's Problem using Load Path Preservation,
Modes 1, 2, 3, and Eigenvectors Components 1, 2, 3.



4.3.1.2 Consistent Modes


In the second problem, the initial model is only incorrect for the spring between masses 4

and 6. A value of 200, two times that of the exact spring, is assumed in this problem. In this

problem, all global and local modes remain global and local modes respectively after

damage. It should be noted that finding a problem with this feature was difficult.

In the first case, as shown in Figure 4.7, it is assumed that only the fundamental mode of

vibration is measured, but all eigenvector components have been measured. It is evident

from Figure 4.7 that the Baruch and SEA-M update are unable to discern the damage, but

that the SEA-SR is able to clearly locate damage. In fact, the SEA-SR was able to exactly

reproduce the correct stiffness matrix. Again, this was true independent of which mode was









Actual Damage Baruch Damage

t100- n 100-


0 0------^-- 0 0 -----jU


-100- -100-

0 10 20 30 0 10 20 30
Indices Indices

SEA-M Damage SEA-SR Damage

100 100





-100 -100-
I II
0 10 20 30 0 10 20 30
Indices Indices

Figure 4.7 Results for Kabe's Problem using the 1st Mode, Full Eigenvectors.



used in the update. It should be noted that the Baruch and SEA-M update tends to focus

elemental changes in the fourth and sixth rows of the stiffness matrix, indicating the

possibility of damage between these degrees of freedom, but certainly giving no clear

indication to the extent of damage. Using Lin's algorithm, the damage vector is given as a =

[1.0 0.98 0.92 0.85 0.87 0.84 0.95 1.0]T. This algorithm does not clearly identify the

damage location.

In the second case, as shown in Figure 4.8, it is assumed that the first three modes of

vibration have been measured, but only DOFs 1,3, and 6 of the eigenvectors have been

measured. In comparing Figure 4.8 to Figure 4.7, it is clear that the damage detection

capability of all three algorithms has again been degraded when using expanded mode

shapes. Only the SEA-SR update gives a clear indication to the location of damage, but is








unable to predict the exact extent. Using Lin's algorithm, the damage vector is given as a=

[0.99 0.18 0.55 0.34 0.52 0.33 0.41 1.0]T. Again, it is difficult from inspection of a to

determine the location of damage. In fact, inspection of a indicates that DOF 2 is the most

likely damaged DOF.


100


0


-100


0 10 20 30
Indices

SEA-M Damage


10 20 30
Indices


01


Baruch Damage


-100l


D 10 20 30
Indices

SEA-SR Damage


100K


0


-100


0 10 20 30
Indices


Figure 4.8 Results for Kabe's Problem using Modes 1, 2, 3,
and Eigenvectors Components 1, 2, 6.





It should be noted that in this problem, it was critical to have the proper DOFs measured.

When the second test case was run with the first three DOFs measured, no algorithm was able

to locate damage. In this case, the eigenvectors components were relatively unaffected by

damage, thus causing substantial error in the eigenvector expansion process.








4.3.2 Model Refinement of a Cantilever Beam: Experimental Study


4.3.2.1 Modal Test Description

The structure used in this investigation is the aluminum cantilevered beam shown in

Figure 4.9. The dimensions and material properties of the beam are given in Table 4.1.

Experimental modal analysis of the beam was performed to measure its modal properties.

Six equally spaced translational degrees of freedom shown in Figure 4.9 were selected as

measurement locations. The modal properties of the first three modes of vibration were

determined using frequency domain techniques and single degree of freedom curve fitting

algorithms. The excitation source used in the testing was an impact hammer and the driving

point measurement was an accelerometer mounted at the free end of the beam. Impact and

exponential windows were utilized to improve frequency response calculations. At each

measured degree of freedom, five frequency response measurements were made and

averaged to reduce the effects of measurement noise. Natural frequencies, damping ratios

and mode shapes of the beam's first three modes of vibration were identified and are reported

in Tables 4.2 and 4.3.


6 5


4 3 2 1


Figure 4.9 Experimental Cantilever Beam.








Table 4.1 Structural Properties of the Cantilever Beam.


Table 4.2 Measured Natural Frequencies and Damping Ratios of the
Cantilever Beam.

Mode # Natural Frequency (Hz) Damping Ratio
(%)
1 7.25 4.41
2 45.55 0.68
3 127.01 0.33


Table 4.3 Measured Mode Shapes of the Cantilever Beam.

Mode 1 2 3
DOF
1 1.00 1.00 1.00
2 0.95 0.16 -0.30
3 0.65 -0.53 -0.61
4 0.36 -0.72 0.20
5 0.15 -0.52 0.75
6 0.03 -0.12 0.28


4.3.2.2 Finite Element Model Description


A twelve DOF undamped finite element model (FEM) of the beam was generated using

six equal length beam elements as shown in Figure 4.9. The beam element has two degrees of

freedom (DOFs) at each node; bending and rotation. This model was then reduced using


Length 0.84 m
Mass/Length 2.364 kg/m
Moment of Inertia 3.02e-9 m4
Youngs Modulus 70 GPa








Guyan reduction (1965) eliminating the rotational degrees of freedom. There are several

possible errors affecting the accuracy of this FEM. The most obvious is the fact that a perfect

cantilever condition is assumed. In addition, an artificial error was purposely introduced by

selecting the Young's Modulus higher than that often assumed for aluminum.

4.3.2.3 Application of the SEAMRA.

Because the "true" finite element model is unknown, a comparison between the "true"

and updated structural matrices is not possible. Besides such comparison, a fair and useful

judgement on the quality of SEAMRA updating capability can be obtained by comparing

actual experimental frequency response functions with those predicted by the initial and

updated FEMs. Figure 4.10 shows a comparison of frequency response functions measured

between degrees of freedom 1 and 3 (i.e. sensor measurement at DOF 1 and impact excitation

at DOF 3). The solid curve corresponds to the experimental data, while the dotted line

corresponds to that predicted by the original analytical FEM. It is apparent that the

discrepancies between the frequency response function increases as the frequency of

excitation increases. This is in part due to the fact that the assumption of a perfect cantilever

condition affects the higher modes of vibration to a greater extent. The dashed lines in this

figure corresponds to the SEA-SR updated finite element model. It is clear from this

comparison that the SEA-SR provided a great deal of improvement to the original analytical

FEM. Inspection of the updated stiffness matrix indicates that changes occur throughout the

matrix, indicating that the discrepancy between the original and refined FEMs was due to

degradation of some global structural property (Youngs Modulus), as opposed to some form

of local damage, as seen in the previous problem.


4.4 Discussion of the SEAMRA's Formulation

In some problems SEAMRA in conjunction with either the subspace rotation or the

modal selection method, failed to find a symmetric updated FEM (symmetric Da and Ka).

This shortcoming was especially encountered in practical situations when the experimental









SEAMRA
10-2 ,
Solid Experimental Measurement
Dash Modified Analytical Model
10-3 Dotted Original Analytical Model


10-4


10.-5 ,


10-6 .,


10-7


10-8 -----
0 20 40 60 80 100 120 140 160 180 200

FREQUENCY (Hz)


Figure 4.10 Experimental and Analytical Frequency Response
Functions of the Cantilever Beam.



modal properties were corrupted by measurement errors. This can be attributed to the fact

that there is no symmetric updated FEM that is consistent with the present SEAMRA's

formulation. Recall that SEAMRA's modifications to the initial stiffness and damping

matrices in its updating process are given by

AK = BoHoBT with Ho = FGo
(4.5.1)
AD = BoHIBo with H, = FG1

Clearly, from Eq. (4.5.1), the perturbations (i.e. modifications) to the initial stiffness and

damping matrices are constrained by the relationship


range(AK) = range(AKT) = range(AD) = range(ADT)


(4.5.2)








This relationship can always be satisfied for the cases when the experimental modal

properties are consistent with an update in which (i) either AK or AD is zero or (ii) AK is

proportional to AD (AK = P AD, P3 is a scalar). For all other cases, SEAMRA might fail to

produce symmetric AK and AD. A more flexible and general formulation that accounts for

such shortcoming is proposed in Section 6.5 of Chapter 6. The formulation as presented in

Chapter 6 is more elegant, efficient and guarantees, for all situations, a symmetric updated

FEM.


4.5 Summary

A previously developed model refinement algorithm based on the general mathematical

framework of eigenstructure assignment theory (Zimmerman and Widengren, 1989, 1990)

has been extended and improved. A technique to enhance eigenvector assignability of the

algorithm has been developed. The method consists of rotating the achievable eigenvector

subspaces into the experimentally measured eigenvectors. The subspace rotation method,

used in conjunction with one of eigenvector expansion techniques discussed in Chapter 2,

results in both a decrease in the computational burden as well as an increase in the accuracy of

the assigned eigenvectors. Finally, the improved algorithm (SEA-SR) was tested for its

suitability for model refinement and structural damage assessment. The performances of

SEA-SR in a damage assessment problem on a challenging simulated structure was

presented and compared to other algorithms. The results acquired using the SEA-SR were

superior.















CHAPTER 5
DAMAGE LOCATION: THE SUBSPACE ROTATION ALGORITHM



5.1 Introduction

In this chapter, a computationally attractive algorithm is proposed to provide an insight

to the location of structural damage. The proposed algorithm is similar to the Modal Force

Error Criteria proposed by several researchers ( Ojalvo and Pilon 1988; Gysin, 1990).

However, a greater insight of the Modal Force Error criteria is provided. Further, a new

viewpoint which allows for the reduction of the effects of measurement errors in the

experimental modal parameters for a certain class of structures is also discussed. As will be

shown in the next sections, the proposed damage location algorithm requires only

matrix-scalar and matrix-vector multiplication.


5.2 The Subspace Rotation Algorithm: The Direct Method

Assume that an n-DOF finite element model of the "healthy" (undamaged) structure

exists. As seen in the earlier chapters, the standard differential equation governing the

dynamic motion of such structures is given by

Mw + Dw + Kw = 0 (5.2.1)

where M, D, and K are the n x n analytical mass, damping, and stiffness matrices, w is a n

x 1 vector of positions and the overdots represent differentiation with respect to time. The

eigenvalue problem associated with Eq. (5.2.1) in second order (lambda) form is given as

(2M + XhiD + K)vhi = 0 (5.2.2)








where kh and vh denote the ith eigenvalue and eigenvector, respectively, of the pre-damaged

"healthy" structure. It is assumed that Eq. (5.2.2) is satisfied for all measured "healthy"

eigenvalues/eigenvectors. This can be enforced by correlating the original FEM (M,D, and

K), possibly through the use of a model refinement procedure.

Next, consider that the p eigenvalues and eigenvectors, 4i and vvi, of a post-damage

modal survey of the structure are available, in which Xd e X i, Vdi VYhi. In the present

formulation, it is assumed that the dimension of the measured eigenvector is the same as the

analytical eigenvector. As discussed in Chapter 2, this is true (i) when all FEM DOFs are

measured (ii) after the application of an eigenvector expansion algorithm, or (iii) after the

application of a finite element model reduction algorithm. The ideal situation would be to

measure all FEM DOFs since the eigenvector expansion process would introduce additional

errors in the "expanded" eigenvectors and the model reduction process would introduce

errors in the FEM. It should be noted that in both cases the additional errors may become

significant as the ratio of measured to unmeasured DOFs become smaller.

Let AMd, ADd, and AKd be the exact perturbation matrices that reflect the nature of the

structural damage. For localized damage, the exact perturbation matrices are sparse matrices

with the nonzero elements reflecting the state of damage. The lambda equation for the

damaged structure is, defined by

(12(M AMd) + d(D ADd) + (K AKd))Vd = 0 (5.2.3)

Although only p of the n eigenvalues/eigenvectors are assumed measured, p << n, Eq. (5.2.3)

holds for any particular eigenvalue and eigenvector of the damaged structure because the

perturbation matrices are assumed to be exact. Grouping all perturbation matrices on the

right-hand side defines a damage vector di,

d. m Zdd (5.2.4a)

= (X2AMd + XdADd + AKd)Vd (5.2.4b)

where Zd = M + d dD + K








Although Eqs. (5.2.4a, b) yield the same damage vector, it should be noted that Zdj and

the coefficient matrix of Eq. (5.2.4b) are not equal. When the measured eigendata are not

corrupted by noise, an inspection ofdi in terms of the Eq. (5.2.4b) reveals that the jth element

ofdi will be zero when thejth rows of the perturbation matrices are zero, i.e. the finite element

model for thejj* degree of freedom is not directly affected by damage. Conversely, a degree

of freedom whose finite element model has been affected by damage will result in a nonzero

entry in di. Thus, the degrees of freedom which have been affected by damage can be

determined by inspecting the elements of di. Vector di as defined in Eq. (5.2.4b) also reveals

that only a single mode of vibration needs to be measured exactly to determine the damage

locations exactly. This is true in even multiple member damage situations. More

importantly, the vector di can be determined from the original finite element model (M,D,K)

and the measured eigenvalues and eigenvectors, Xd and vdi, using Eq. (5.2.4a). Thus, there is

no need to use a model refinement algorithm to attempt to estimate the exact perturbation

matrices in order to locate the damage. If the damping term is ignored, Eq. (5.2.4a) is

essentially the Modal Force Error criteria as proposed by Ojalvo (1988) for use as a

diagnostic "tool" to locate modelling errors in FEMs. A physical interpretation of Eq.

(5.2.4a) provided by Ojalvo was "... di is the applied harmonic force error distribution,

applied at frequency Xdi, which is necessary to cause the analytical model to vibrate with

mode shape vd ...". However, Eq. (5.2.4b) provides a much clearer interpretation of the

damage vector di for the damage location problem in which the perturbation matrices are

sparse.

In practice, the perfect zero/nonzero pattern of the damage vector di rarely occurs due to

errors present in the experimentally measured eigenvalues and eigenvectors. This scenario

was studied and discussed for the undamped case by Gysin (1990) in the context of

eigenvector expansion techniques. Gysin (1990) observed that in certain specific cases of

eigenvector expansion errors, the damage vector defined by Eq. (5.2.1a) may lead to

incorrect conclusions concerning the location of damage. The next section proposes a new








viewpoint which allows for the reduction of the effects of measurement errors for certain

classes of structures.


5.3 The Subspace Rotation Algorithm: The Angle Perturbation Method

In order to provide an alternative view of the state of damage, Eq. (5.2.4a) is rewritten as

d 4 = 1 IId II dcos(W) (5.3.1)


where d! is the jth component (or j* DOF) of the ith damage vector, zJ is the j* row of the

matrix Zd and 01 is the angle between the vectors zJ and vd*

In the case when the measurements are free of error, a zero dJ corresponds to a 06 of
I I

ninety degrees, whereas a nonzero d. corresponds to a OQ different from ninety degrees.
1 I

Errors in the experimental measurements of modal parameters will cause slight

perturbations in the angles 06 that destroy the zero/nonzero pattern of the damage vector.
I

One would initially expect that the components of di corresponding to the damaged DOFs

would be substantially larger than the other elements. However, by inspecting Eq. (5.3.1), a

large di component could be due to a zJ row norm substantially larger than other rows of Zd.,
di ,

coupled with a slight deviation of O9 from ninety degrees due to measurement noise. Hence,

when dealing with a structure whose FEM results in zJ row norms of different order of
-di

magnitude, it is more reasonable to use the deviation of the angles, W1, from ninety degrees

for damage location,

ai = 18 900 (5.3.2)


d!
where Oi = cos 1 -I d



The angle 01 is determined from Eq. (5.3.2) and a1 is the jth component of ai.
I I








5.4 Practical Issues


5.4.1 Cumulative Damage Location Vectors

The discussion in the previous section suggests that for a given mode, the damage is

locatable if the perturbation to the angle, OQ due to the measurement error is less than the
I
angle perturbation due to the damage. Hence, modes that are highly affected by the damage

are expected to provide better assessment to the location of the damage when errors are

present in the measured eigendata. As will be seen in example problems (Chapter 7), certain

modes are more susceptible to a given state of damage than others. This is mainly due to the

fact that different elements of the structure have different levels of contribution to the total

strain energy of a given mode (Kashangaki 1992). Furthermore, a particular part of the

structure usually has different orders of strain energy contribution for different modes.

Usually, if the damage occurs in a region of high strain energy for a given mode, that mode

would be highly susceptible to the damage and, hence would reflect the state of damage. To

accommodate this type of problem, when the number of measured modes p is greater than

one, two different composite damage vectors may be defined as


d i= 1 di


a 1 i (5.4.2)
-p
i=l

In Eq. (5.4.1), the damage vectors, di, are normalized with respect to their corresponding

eigenvectors, Vdi. The reason for this normalization is such that the composite vector ignores

the inherent "weighting" of I| yd 1i, which is usually of different orders of magnitudes for

different measured modes. It should be noted that in the multi-mode measurement case, Eq.

(5.4.2) is preferable when the values of 1| zJ || are of different orders of magnitude for
different measured modes. Again, in practice,
different measured modes. Again, in practice, the DOFs affected by the damage are expected








to have substantially larger d or _q. Finally, the damaged areas of the structure can then be

located using the knowledge of the "damaged" DOFs and the connectivity of the FEM.

It is interesting to note that Eqs. (5.2.4a, b) reveal an interesting relationship between

various model refinement algorithms. Model refinement techniques attempt to approximate

the exact perturbation matrices by using limited modal data, but do so in different manners.

However, Eq. (5.2.4a, b) indicates that if the model refinement technique has satisfaction of

the eigenproblem as an equality constraint, the calculated perturbation matrices AM, AD,

and AK are constrained to be related to the original finite element model M, D and K and the

measured eigendata by Eqs. (5.2.4a, b).

5.4.2 Eigenvector Filtering Algorithm

In a modal survey, the errors associated with the measured eigenvectors are typically

greater than the error associated with the measured eigenvalues. In addition, in the case of

incomplete eigenvector component measurements, these measurement errors are often

compounded with eigenvector expansion induced errors. A simple eigenvector noise

filtering algorithm is proposed assuming the measured eigenvalues to be correct. From the

cumulative damage vector defined in Eqs. (5.2.4a) or (5.2.4b), and the original FEM

connectivity, the engineer can deduce which DOFs have been damaged. It is reasonable to

assume that nonzero elements in each damage vector di associated with "undamaged DOFs"

are due to eigenvector errors. These elements can then be set to zero. In addition, the

magnitude of the elements ofdi at the "damaged" DOFs can be adjusted by using knowledge

of DOF connectivity and the properties of the element property matrices connecting the

"damaged DOFs." The element property matrices provide constraints relating the effect of

damage on each element DOF. The noise filtering algorithm consists simply of replacing the

di vectors by df, where df is obtained from di as described above. The ith filtered

eigenvector, vdf, can then be obtained from solving

(dXM + ,dD + K)vdf = df (5.4.3)








using Gaussian elimination. In this calculation, the bandedness of typical FEM matrices

should be exploited. Essentially, the filtered eigenvector is just the eigenvector that if

measured would have produced the damage vector df. Experience gained in using the

eigenvector filtering algorithm indicates that it is best to use structural matrix properties (M,

D, K) that (i) are finite element consistent, and (ii) have not been "corrupted" by

measurement noise. By finite element consistent, it is meant that the property matrices can

be achieved by a finite element program. Note that measurement noise can be introduced

in the property matrix through a FEM refinement algorithm. Hence, the property matrices

that should be used are the original property matrices (unrefined). The present eigenvector

filtering algorithm can be useful in improving the damage extent assessment. This concept

is discussed in Chapters 6 and 7.


5.5 Summary

A computationally attractive algorithm to determine the location of damage in

structures was developed. The algorithm completely bypasses the general framework of the

model refinement problem and involves only matrix-scalar and matrix-vector

multiplications. The effect of measurement error in the eigendata was discussed and

techniques to reduce these effect were presented. Furthermore, a simple eigenvector

filtering algorithm was developed. Practical example problems to illustrate and evaluate the

performance of the developed algorithm will be presented in Chapter 7.














CHAPTER 6
THE MINIMUM RANK PERTURBATION THEORY



6.1 Background

The theory developed in Chapter 5 is limited to determining the location of structural

damage. In a practical situation, it is essential to determine the extent of the damage to get a

good estimate about the overall integrity of the structure. In general, the extent problem, as

discussed in Chapter 1, has been approached by several researchers using existing model

refinement algorithms. The formulations of these algorithms were obviously based on the

"model refinement philosophy": minimum change made to the original FEM. The minimum

change constraint has a clear tendency to smear the changes throughout the entire FEM.

However, in most cases, this philosophy is not consistent with the effect of structural damage

on FEMs. In fact, the effects of structural damage on FEMs are usually "non-minimal"

localized perturbations. Structural damage often occurs at discrete locations. The effect of

damage on the analytical model is often restricted to just a few elements of the finite element

model. The rank of each element mass, damping or stiffness matrix is dependent on the

number of degrees of freedom defined by the element and the shape functions utilized.

However, it should be noted that in general the element matrices are not of full rank. For

example, the rank of the 6x6 element stiffness matrix of a three dimensional truss element is

just one. Thus, instead of using the matrix Frobenius norm minimization formulation to

arrive at unique perturbation matrices, minimum rank perturbation constraints are enforced.

In this Chapter, a computationally attractive damage extent algorithm is proposed. The

proposed damage extent algorithm is a minimum rank perturbation, which is consistent with

the effects of many classes of structural damage on a FEM.








Assume that "p damaged" eigenvalues and eigenvectors have been measured and that

the original FEM has been corrected such that its modal properties match the measured

modal properties of the healthy model. The eigenvalue problem of a damaged structure

shown in Eq. (5.2.3), for all p measured modes, can be written in matrix form, as

MVdAd + DVdAd + KVd = AMdVdAd + ADdVdAd + AKdVd B (6.1.1)

where Ad = diag(Xd,. d2' ... dp)

Vd= [d,,Y.d2, ..., Y-d

B = [dl,d2 ... dp]


where all variables have the same definitions as in the previous chapter. Note that matrix

B can be determined from the FEM (M, D, K) and the "p" measured eigenvalues and

eigenvectors. As discussed earlier, the damage extent problem consist of finding the

perturbation matrices, AMd, ADd, and AKd, such that Eq. (6.1.1) is satisfied. As already

discussed in Chapter 2, structures can be modelled using either undamped, proportionally

damped or nonproportionally damped finite element models. The proposed extent algorithm

is formulated to accommodate all three types of structural models. For each type of model,

several scenarios of damage effects are considered. Practical issues that can be used to

improve the damage extent estimate are also presented.


6.2 The Minimum Rank Perturbation Theory: Theoretical Background

In this section, the theoretical foundation of the Minimum Rank Perturbation Theory

(MRPT) is derived. This theory will be extensively used throughout the remainder of this

chapter.

PROPOSITION 6.1 Suppose that X, Y E RnxP are given where p < n and

rank(X)=rank(Y)=p.

Define % to be the set of matrices A in R"' that satisfy,








AX = Y with AT = A


Then,

(1.a) If the set % is nonempty, the minimum rank of any matrix, A, in % is p.

Next, define %P to be a subset of % comprised of all A such that rank(A) =p.
Then

(1.b) If the matrix YTX is symmetric, then one member of %P is given by


AP = YHYT


(6.2.2)


with H = (yTX)


and


(1.c) The matrix defined by Eq. (6.2.2) is the unique member of %9P.



Proof:
To prove Proposition (6.1.a), note that Eq. (6.2.1) is exactly satisfied if and only if
range(Y) is included in range(A), which is also the range(AT) by symmetry. This implies that
rank(Y) = p : rank(A).

///

To investigate Proposition (6.1.b), assume that the expanded singular value

decomposition of one member, AN, of 9%P to be of the form


(6.2.3)


where


Uj = [u{, u ... ..]

XJ = diag(oJ{, oJ, ... i)


(6.2.1)


APJ = Uj:jUjw








where the superscript j indicates the jth family member of 9%P, the uj are the left and right

singular vectors and the ao are the non-zero singular values of Ap. In the expanded singular

value decomposition, the (p+l) to n singular vectors are not shown in the factorization

because of their corresponding zero singular values. Note that the left and right singular

vectors are the same because AN is restricted to be symmetric. For Eq. (6.2.1) to be satisfied,

the range of Y, AP'N, ApJT and Ui must be equal. Thus, any column of Y can be written as

a linear combination of the u.'s. The matrices Y and Ui are then related by a unique pxp

invertible matrix QJ,

Y = UJ Q (6.2.4)

Substituting Eq. (6.2.4) into Eq. (6.2.3) gives

AP = Y(Qj-12jQJ-T)YT = yHJyT (6.2.5)


Thus, each family member is uniquely defined by the factorization of Eq. (6.2.3). From Eq.

(6.2.5), it is evident that HJ is of full rank because its inverse exists (HI-' = QjTFji~Q).
Inspection of Eq. (6.2.5) reveals that the only unknown term in the factorization is Hi.

By using the factorization of AN as defined by Eq. (6.2.5), Eq. (6.2.5) can be rewritten as

Y = A-PX = (YHjiYX = Y(HjYTX) (6.2.6)


Equation (6.2.6) is satisfied if and only if HJYTX = Ipxp, where Ipxp is the pxp identity

matrix. This is true because Y and X are of full column rank. Thus, Hi is uniquely calculated
to be

Hi = (yTX)1 (6.2.7)

///
Proposition (6.1.c) follows immediately by inspecting the right hand side of Eq. (6.2.7).
Inspection reveals that HJ is the same for all members of %P. This fact, in conjunction with








Eq. (6.2.5) leads to the conclusion that AN'j is the unique member of the set %P. This member

is given by Eq. (6.2.2).





At this point, the MRPT as defined in Proposition 6.1 assumes that the matrices X and Y

are of full rank. In practical uses, as will be seen later, matrix X is usually of full rank. The

rank requirement on matrix Y can be of some concern since it is directly related to the rank of

matrix A. The next proposition addresses the case in which matrix Y is rank deficient.


PROPOSITION 6.2 Suppose that X, Y E Rnxm are given and rank(X)=m and rank(Y)=p,

where p < m < n. Further, suppose that the matrix YTX is symmetric.

Define cl, to be the set of matrices A in Rnxn that satisfies the problem,

APXP = YP with (AP)T = AP (6.2.8)

where the superscript p indicates a rank p matrix. In Eq. (6.2.8), XP, YP E Rnxp are

corresponding full rank submatrices of X and Y.

Then

The set U contains a single member, AP, that can be calculated from Eq. (6.2.2) using

any corresponding XP and YP.

Proof:


The jth member of the set cM is given by

AP = YPHjYPj (6.2.9)

with HJ = (YPJTXPj)

where the additional superscript ( )-J indicates the jth member of c. Note that Hi is

symmetric since YTX is symmetric. The range of any YP is equal to the range of Y, thus

the YP'1 and YP'J are related by








YP- = YPj Qi' (6.2.10)

where Qid C Rpxp and rank(QiJ) = p. By utilizing Eq. (6.2.10), Eq. (6.2.9) can be written

for the ith member of cU as

Ap, = YPJ(Qi'HiQiT)YPJT (6.2.11)

with H' = (QijTYPTp,i) 1

Again, Hi is symmetric because YTX is symmetric. In comparing Eqs. (6.2.9) and (6.2.11),

it is seen that AN = Ap' if

HJ = Q'J H' QjiT (6.2.12)

or equivalently,

Hij- = Qi'-THi-TQ'ji' (6.2.13)

where Eq. (6.2.13) makes use of the symmetry of Hi. By using the definitions of Hi and H9,

Eq. (6.2.13) can be written as

YPjTXP' = Qi,-TXPiYPJ (6.2.14)

Pre-multiplying Eq. (6.2.14) by Qi'j and utilizing the relation in Eq. (6.2.10) gives the

condition such that AN = AP", namely

YPirXPj = XpirYPj (6.2.15)

This condition is clearly satisfied since YTX is symmetric.



The actual uses of Proposition 6.1 and Proposition 6.2 will be clearly seen in the next

sections. The practical implication of Proposition 6.2 is discussed in detail in Section 6.6.


6.3 Damage Extent: Undamped Structures

In some cases, the damping of the system under consideration is assumed to be

negligible. For this type of system, MRPT based algorithms will be developed assuming that








the structural damage affects (i) only the mass properties, or (ii) only the stiffness properties,

or (iii) simultaneously the mass and stiffness properties.

6.3.1 Damage Extent: Mass Properties

In this case, it is assumed that the effect of damage on the stiffness properties of the

structure is negligible. With this assumption, Eq. (6.1.1) can be rewritten as
MVdA2 + KVd = AMVdA2 = B (6.3.1)


Note that the eigenvectors are real and the eigenvalues are purely imaginary. Further, the

eigenvectors are linearly independent, which implies that the matrix product VdA2 is of full

column rank if rigid body modes are not included. Assume, for the moment, that B is of full

rank (rank(B) = p). Then, Proposition 6.1 can be applied to determine the perturbation

matrix, AMd, as

AMd = B(BTVdA) BT (6.3.2)

by letting Y=B and X=VdA Note that the required inversion is that of a pxp matrix, where

"p" is the number of measured modes. As discussed in Proposition 6.1, this inversion is

feasible if matrix B is of full rank and the rigid body modes of the system are omitted in the

computations. When matrix B is rank deficient, Proposition 6.2 should be used to render

the computation possible.

The properties associated with AMd as computed in Eq. (6.3.2) are as follows:

PROPERTY 6.3.1 The perturbation matrix, AMd, defined in Eq. (6.3.2) will be symmetric if

the eigenvectors, Vd, are stiffness orthogonal, i.e., the eigenvectors are orthogonal with

respect to the original stiffness matrix, K.

Proof:

Proposition (6.1.c) in conjunction with Proposition (6.1.b) implies that the existence of

the unique symmetric rank p AMd requires the symmetry of the matrix product BTVdA 2

The symmetric equivalence associated with this matrix product is








BTVdA -= A VB (6.3.3)

Substituting the expression for B from Eq. (6.3.1) into Eq. (6.3.3) gives

AdV MVdAd + V KVdA2 AdVdMVdAd + A2VTKVd (6.3.4)

where the symmetry of M, K and A2 has been used in writing Eq. (6.3.4). From Eq. (6.3.4),

it is clear that the equivalence is true if

(VKVd)A == Ad(V KVd) (6.3.5)

Equation (6.3.5) will obviously be satisfied if the measured eigenvectors are stiffness

orthogonal. Baruch (1978) treated one approach to mass orthogonalize the measured

eigenvectors. A similar approach can be used to orthogonalize the measured eigenvectors

with respect to the stiffness matrix.




PROPERTY 6.3.2 The updated finite element model (FEM) defined by the original mass

and stiffness matrix along with the perturbation mass matrix computed using Eq. (6.3.2)

preserve the rigid body characteristics of the original FEM.

Proof:

This is apparent in that the original stiffness matrix is unchanged and that the rigid body

modes are defined as modes whose eigenvectors lie in the null space of the stiffness matrix.

6.3.2 Damage Extent: Stiffness Properties

Here, it is assumed that the effect of damage on the mass properties of the structure is
negligible. With this assumption, Eq. (6.1.1) can be rewritten as

MVdA2 + KVd = AKdVd -= B (6.3.6)

For this problem, the eigenvectors are real and the eigenvalues are purely imaginary. The
eigenvectors are also linearly independent, which implies that matrix Vd is of full rank. If








matrix B is assumed to be of full rank (rank(B)=p), Proposition 6.1 can be used to determine

the perturbation to the original stiffness matrix,

AKd = B(BTVd) IB (6.3.7)

This expression for AKd is determined by setting Y=B and X=Vd in Eq. (6.1.3).

The properties associated with AKd as computed by Eq. (6.3.7) are as follow.

PROPERTY 6.3.3 The matrix AKd will be symmetric if the eigenvectors are mass

orthogonal, i.e., the eigenvectors are orthogonal with respect to the original mass matrix.

The proof of Property 6.3.3 follows very much the same pattern as the one presented for AMd

(Property 6.3.1).

PROPERTY 6.3.4 The updated FEM defined by the original mass and stiffness matrices and

the perturbation stiffness matrix, AKd, preserves the rigid body mode characteristics if the

measured eigenvectors and the rigid body modes are mass orthogonal.



The original rigid body modes of an undamped system are defined by the eigenvalue

problem,

KVr = XrMVr = 0 (6.3.8)

where the subscript r denotes the rigid body mode(s) and Xr is equal to zero. Thus, the rigid

body modes lie in the null space of the original stiffness matrix. The rigid body modes of

the system will be preserved in the updated model if the original rigid body modes lie in the

null space of the updated stiffness matrix,

e = (K AKd)Vr (6.3.9a)
= AKdYr (6.3.9b)

where vector e is zero if the the rigid body modes are preserved. Equation (6.3.8) has been

used to arrive at the expression shown in Eq. (6.3.9b). Substituting Eq. (6.3.7) into (6.3.9b)

gives








e = B(BTVd) BTr (6.3.10)

By utilizing the symmetry of the original mass and stiffness matrices, along with Eq. (6.3.6),
Eq. (6.3.10) can be expanded as

e = B(BTVd)- [VTKKv + AdVdMv,] (6.3.11)

The first term in the parenthesis is zero because the matrix-vector product Kvr is zero by
definition. The second term will be zero if the rigid body modes and the measured mode
shapes are mass orthogonal.



6.3.3 Damage Extent: Mass and Stiffness Properties

In this case, it is assumed that the structural damage affects simultaneously the mass and
stiffness properties of the structure. With this assumption, Eq. (6.1.1) can be rewritten as

MVdAd + KVd = AMdVdA2 + AKdVd =- B (6.3.12)

6.3.3.1 Application of The MRPT

Assume that Eq. (6.1.12) can be decoupled as follows,

AMdVd = Bm (6.3.13a)
AKdVd = Bk (6.3.13b)

Then, the MRPT, as formulated in Proposition 6.1, can be applied to determine the
perturbation matrices AMd and AKd, as

AMd=Bm (BLVd) -1 B (6.3.14a)

AKd = Bk (B Vk) Bk (6.3.14b)


Note that the matrices BLVd and B Vd are invertible if Bm and Bk are of full rank. When

these rank requirements are not met, Proposition 6.2 can be used to make the computations
possible.








6.3.3.2 Decomposition of Matrix B

The decomposition problem as illustrated in the previous section is equivalent to the

problem of solving for the matrices Bm and Bk. So far, the only constraint that these

unknown matrices must satisfy is given by the expression,

B = Bm A2 + Bk (6.3.15)

which results from Eqs. (6.3.12) and (6.3.13). Naturally, there is an infinite set of solutions

(Bm, Bk) that satisfy Eq. (6.3.15). To arrive at a unique solution, additional physically

meaningful constraints can be enforced. The decomposition proposed herein exploits the

cross-orthogonality relations that arise from the symmetric nature of the property matrices

and the undamped assumption. By measuring mass normalized "damaged" eigenvectors

(which is possible if a driving point measurement is made), the cross-orthogonality relations

associated with the damaged structure can be written as

V(M AMd)Vd = Ipxp (6.3.16a)

VTK AKd)Vd = diag( d,2, 2) = 2 d (6.3.16b)

in which Odi is the natural frequency of the ith mode of the "damaged" structure. Matrix

Ipxp is the pxp identity matrix. A rearrangement of Eq. (6.3.16) yields

VT AMd Vd = VT M Vd Ipxp VT Bm (6.3.17a)

VdT AKd Vd = V K Vd d = V Bk (6.3.17b)

Clearly, the matrices Bm and Bk can be computed from Eqs. (6.3.17). In the rare situation

that the number of measured modes is equal to the number of DOFs in the FEM (p = n), these
can be computed by simply inverting matrix Vd. Unfortunately, as discussed earlier, the

number of measured modes is usually much less than the number of FEM DOFs (p << n).
In this case, the solution that naturally comes to mind is to use the pseudo-inverse of matrix

VT. The inconvenience of this approach is that the sparsity pattern of matrix B will not be








reflected in the computed matrices Bm and Bk. Remember that the sparsity pattern of B, as

discussed in Chapter 5, indicates the location of the damage affecting the structure. A more

physically intuitive approach is to constrain Bm and Bk to exhibit the same sparsity pattern

as matrix B. This is done by casting B in an equation similar to the expressions of Eqs.

(6.3.17). The problem in question is then to find an nxp matrix P that satisfies

P(VTB) = B (6.3.18)

Matrix P can be computed as


P = B (VTB) (6.3.19)

The inverse involved in this computation is that of a pxp matrix which is invertible if matrix

B is of full rank. Now that P is computed, Bm and Bk can be computed using Eq. (6.3.19)

as


B = P(V M Vd Ipxp) (6.3.20a)

Bk = P K Vd Qd) (6.3.20b)

It is clear from Eq. (6.3.19) that P will have the same sparsity pattern as matrix B. Hence

Bm and Bk will also reflect the important sparsity pattern of B. The computed matrices Bm

and Bk can also be used to determine the effect of the damage, respectively, on the mass, and

stiffness properties. As in Chapter 5, cumulative vectors associated to Bm and Bk can also be

defined when more than one measured mode is available.

d 1 IdrI
Am = I --i (6.3.21a)


dk Idl (6.3.21b)
i=111 di


where dm. and dki are, respectively, the ith column of matrix Bm and Bk.








PROPERTY 6.3.5 The perturbation matrices (AMd, AKd) computed from the MRPT using
the Bm and Bk resulting from the decomposition discussed above will be symmetric.

Proof:
The perturbation matrices AMd and AKd will be symmetric because they satisfy the
relationships in Eqs. (6.3.17) and the right hand sides of these equations are symmetric.


6.4 Damage Extent: Proportionally Damped Structures

Since many structures have non-negligible damping, it is of practical interest to extend
the MRPT to address damped structures. In this analysis, the structure under consideration is
assumed to exhibit proportional damping.

6.4.1 Damage Extent: Stiffness and Damping Properties

It is assumed that the effect of the structural damage on the mass properties is negligible.
In this context, Eq. (6.1.1) is rewritten as

MVdA2 + DVdAd + KVd = ADdVdAd + AKdVd B (6.4.1)

The complex conjugate of Eq. (6.4.1) is

ADdVdAd + AKdVd = B (6.4.2)

where the overbar indicates the complex conjugate operator, and the fact that ADd, AKd and
Vd are real has been used in writing Eq. (6.4.2). Subtracting Eq. (6.4.2) from Eq. (6.4.1)
gives

ADdVd(Ad Ad) = (B B) (6.4.3)

If (B B) is assumed to be of full rank, Proposition 6.1 can be applied to determine the

perturbation matrix, ADd, as

ADd = (B B)Hd(B BT

with Hd = [(B TVd(Ad d)-1 (6.4.4)







Note that ADd as defined by Eq. (6.4.4) is real. Post-multiplying Eq. (6.4.1) by Ad and
Eq. (6.4.2) by Ad and subtracting the two equations leads to

AKdVd(Ad Ad) = (BAd BAd) (6.4.5)

where the fact that Ad and Ad are diagonal matrices has been used in writing Eq. (6.4.5).

If (BAd BAd) is assumed to be of full rank, Proposition 6.1 can also be applied to
determine the perturbation matrix, AKd, as

AKd = (BAd BAd)Hk(BAd BAd)
-1 (6.4.6)
with Hk = [(BAd BAd)TVd(Ad Ad)] (6.4.6)

Note that AKd as defined by Eq. (33) is also real.

PROPERTY 6.4.1 The perturbation matrices ADd and AKd, as computed above, will be
symmetric if the measured eigenvectors, Vd, are mass orthogonal; i.e., the eigenvectors are
orthogonal with respect to the original unperturbed mass matrix.


Proof:
Matrix ADd is symmetric if Hd is symmetric or, equivalently, if Hd 1 is symmetric.
Hence, to get a symmetric ADd, the following equivalence must be satisfied.

(B B)TVd(Ad Ad) = (Ad AVT(B B) (6.4.7)


Substituting the expressions for B and B, from Eqs. (6.4.1) and (6.4.2) respectively, into
Eq. (6.4.7) yields

(A2VTM + AdVTD AdV2TM AdVD Vd Ad)
d d d (6.4.8)
(Ad Ad)V (MVdA2 + DVdAd MVdd + DVdAd







Note that in Eq. (6.4.8) the terms involving matrix AKd canceled out. A further expansion
and simplification of Eq. (6.4.8) yields

(Ad Ad) VMMV (Ad d) (Ad d) VMVd (A Ad) (6.4.9)


which is clearly satisfied if the measured "damaged" eigenvectors, Vd, are mass orthogonal.


Likewise, the perturbation matrix AKd as computed in Eq. (6.4.6) is symmetric if Hk is

symmetric or, equivalently, if Hk 1 is symmetric. This symmetry requirement yields the
following equivalence.

(BAd BAd)TVd(Ad Ad) (Ad Ad)VT(BAd BAd) (6.4.10)

Substitution of the expressions for B and B into Eq. (6.4.10) yields

(AdA2VM + AdVK AdAdV M AdVK) Vd(Ad Ad)
T 2 2 (6.4.11)
S(Ad Ad)Vd (MVdA dAd + KVdAd dAd KVdAd)

in which the terms involving matrix ADd cancel. Manipulating and simplifying Eq. (6.4.11)
yields

(?dAdA A dAd)VMVd(Ad Ad) = (Ad Ad)V dMVd(XdAd AdAd) (6.4.12)

This equivalence is obviously satisfied if the eigenvectors are mass orthogonal.
///




PROPERTY 6.4.2 The updated FEM, defined by the original FEM and the perturbation
matrices ADd and AKd computed from Eqs. (6.4.4) and (6.4.6), preserves the original rigid
body modes if the measured eigenvectors and the rigid body modes are mass orthogonal.








Proof:
As discussed earlier, a rigid body mode is defined as a mode whose eigenvalue is equal to
zero and whose eigenvector lies in the null space of the FEM stiffness matrix. Hence, the
rigid body modes of the original system are preserved in the updated FEM if they lie in the
null space of the perturbed stiffness matrix. Consider the relationship

e = (K AKd)Vr (6.4.13)

where Yr is a rigid body mode eigenvector. Clearly, the rigid body mode associated to
eigenvector Yr is preserved if e = Q. By definition, _r is a rigid body eigenvector of the
original system, hence Eq. (6.4.13) can be simplified as

e = AKdYr (6.4.14)

Substituting the expression for AKd as defined in Eq. (6.4.6), into Eq. (6.4.14) gives

e = (BAd BAd)Hk(BAd BAd) r (6.4.15)

Substitution of the expressions for B and B into this equation yields

e = (BAd Ad)Hk[ AdAVM + A K AdAdVdM AdVK ]vr (6.4.16)

By using the fact that vr is a rigid body eigenvector of the original system (i.e. Kvr = 0), Eq.
(6.3.16) can be simplified as

e = (BAd BAd)Hk( AdA2 Ad )VMvr (6.4.17)

It is clear from Eq. (6.4.17) that e = 0 if the rigid body mode Yr and the measured

eigenvectors Vd are mass orthogonal (i.e. VTMvr = 0).



6.4.2 Damage Extent: Mass and Damping Properties

In this case it is assumed that the effect of the structural damage on the stiffness
properties is negligible. In this context, Eq. (6.1.1) is rewritten as







MVdAd + DVdAd + KVd = AMdVdA2 + ADdVdAd = B (6.4.18)

By using an approach similar to one used in the preceding section, Eq. (6.4.18) and its
complex conjugate can be manipulated to yield the following decomposition

AMdVd(Ad Ad) (BAd BA) (6.4.19)

ADdVd(AdAd Ad) = (BAd dBA) (6.4.20)

Again by applying the MRPT to the preceding equations, AMd and ADd are determined to
be

AMd = (BAd Ad)Hm(BAd BAd)T
S-I (6.4.21)
with Hm = (BAd- 2BA) V d Ad


ADd = (BAd RAd)Hd(BAd TA2

with Hd ( BA) Vd(Add dA(6.4.22)

Clearly, the perturbation matrices AMd and ADd as defined by Eqs. (6.4.21) and (6.4.22) are
real.

PROPERTY 6.4.3 The perturbation matrices AMd and ADd, as computed above, will be
symmetric if the measured eigenvectors, Vd, are stiffness orthogonal; i.e., the eigenvectors
are orthogonal with respect to the original unperturbed stiffness matrix.

PROPERTY 6.4.4 The updated FEM, defined by the original FEM and the perturbation
matrices, AMd and ADd, preserves the original rigid body modes.
The proof of Property 6.4.4 is straightforward since the original stiffness matrix is
unchanged (see Property 6.3.2). The proof of Property 6.4.4 follows very much the same
pattern as the proof of Property 6.4.2.







6.4.3 Damage Extent: Mass and Stiffness Properties

In this problem, it is assumed that the effect of the structural damage on the damping
properties is negligible. For this situation, the general eigenvalue problem defined in Eq.
(6.1.1) associated to this case can be simplified as

MVdA2 + DVdAd + KVd = AMdVdA2 + AKdVd B (6.4.23)

Algebraic manipulations of Eq. (6.4.23) and its complex conjugate yield the following
decomposition

AMdVd(A Ad) = (B B) (6.4.24)

AKdVd(A A) = (BAd dBA) (6.4.25)

The perturbation matrices AMd and AKd can then be computed using the MRPT.

AMd = (B )Hm(B B)T

with Hm = (B BTVdA2 1 (6.4.26)

T
AKd = (BAd BA)HA BAS d
T (6.4.27)
with Hk (BA A) Vd( (6.4.27)

Note that AMd and AKd as defined by Eqs. (6.4.26) and (6.4.27) are real.

PROPERTY 6.4.5 The perturbation matrices AMd and AKd, as computed above, will be
symmetric if the measured eigenvectors, Vd, are damping orthogonal; i.e., the eigenvectors
are orthogonal with respect to the original unperturbed damping matrix.


PROPERTY 6.4.6 The updated FEM, defined by the original FEM and the perturbation
matrices, AMd and AKd, preserves the original rigid body modes if the measured
eigenvectors and the rigid body modes are damping orthogonal.







These proofs of the above two properties are not reported here. They follow very much the
same pattern as the proofs in Section 6.4.2.

6.4.4 Damage Extent: Mass. Damping and Stiffness Properties

The eigenvalue problem of a proportionally damped system with all property matrices
simultaneously affected by damage can be rearranged into the form

MVdAd + DVdAd + KVd = AMdVdA + ADdVdAd + AKdVd = B (6.4.28)

The theory developed in Section 6.3.3 can be expanded to address this particular problem.
The cross-orthogonality relationships associated with this type of structures are

Vj(M AMd)Vd = Ipxp (6.4.29a)

VJ(D ADd)Vd = diag(2Odd1, ,2E ddp) = (6.4.29b)

VT(K AKd)Vd = diag(wd2, .. ap2) = d (6.4.29c)

Notice that the cross-orthogonality relationships in Eqs. (6.4.29a) and (6.4.29c) are exactly
the same as the ones associated with undamped systems reported in Eqs. (6.3.16a) and
(6.3.16b). As before, these cross-orthogonality conditions can also be rearranged as

V AMd Vd = V M Vd pxp = V Bm (6.4.30a)

VT ADd Vd = V D Vd d VT Bd (6.4.30b)

VT AKd Vd = V K Vd d = V Bk (6.4.30c)

Following the exact same argument discussed for undamped systems in Section 6.3.3, an nxp
matrix P that satisfies the relation,

P(V dB) = B (6.4.31)

is sought, where B is computed using Eq. (6.4.18) and (vdTB) is a pxp matrix. Although B

is a complex matrix, the nxp matrix P is real, since Vd is real. Hence, for computational
efficiency, matrix P can computed from








P = Br(V Br) (6.4.32)


where Br is the real part of B. In Eq. (6.4.32), it is assumed that matrix (VTB) is invertible.

With P computed, the next step is to determine the decomposed damage vectors that indicate
the effects of the damage on the mass, damping and stiffness matrices,

Bm = AMVd = P (V M Vd pxp) (6.4.33a)

Bd = ADVd = P (V D Vd Yd) (6.4.33b)

Bk = AKVd = P (V' K Vd d) (6.4.33c)

The minimum rank perturbation theory (MRPT), as formulated in Proposition 6.1, can again
be applied to determine the perturbation matrices, AMd, ADd and AKd, as

AMd = Bm (BTVd)-1 Bm (6.4.34a)

ADd =Bd (B Vd)-1 BT (6.4.34b)

AKd = k (BVk) -1 Bk (6.4.34c)


Note that the matrices BTVd, BJVd and B Vd are pxp matrices that are invertible if Bm,

Bd and Bk are of full rank. As in all other cases already studied, Proposition 6.2 can be used
to deal with the situation when any one of these matrices are rank deficient. The cumulative
damage location vector associated to Bm and Bk, defined in Eqs. 6.3.21, are also applicable
to this problem. An additional cumulative damage vector associated to the perturbations in
the damping properties can be similarly defined as


dd = I (6.4.35)
pi=l 1 d


where dd is the ith column of matrix Bd.