ZhouShiyu air force PI meeting 2016 final
AFOSR Project FA9550-14-1-0384
Progressive Fault Identification and Prognosis in Aircraft Structure Based on Dynamic Data Driven Adaptive Sensing and Simulation
Shiyu Zhou University of Wisconsin-Madison
Jiong Tang University of Connecticut
Yong Chen University of Iowa Air Force DDDAS Program PI Meeting, 01/27/2016
1
Outline • Introduction and project objectives • Research results accomplished the first project year Identifiability analysis for finite element based dynamic damage detection Establishment of piezoelectric impedance-based damage identification algorithms and test-beds
• Summary and on-going efforts
Air Force DDDAS Program PI Meeting, 01/27/2016
2
Introduction • Structural health monitoring and management (SHM2) is of critical importance for air force applications • Dynamic responses of an aerospace structure contain rich information • Limitations of the current dynamic response based SHM2 systems – Sensory system is fixed and re-active – Applications modeling is either based on mechanistic principles or simple data driven models using historical data – Sensor-structure interaction dynamics are not exploited; Uncertainties not fully addressed
Air Force DDDAS Program PI Meeting, 01/27/2016
3
Research objectives Overarching Goal: To create a new methodology of progressive structural fault identification and prognosis for Air Force applications based on the framework of dynamic data-driven applications systems (DDDAS). Applications measurements systems and methods Piezoelectric transducer
Vi
i
Applications modeling
Mathematical and statistical algorithms Failure time prognosis
Data collection
Data-driven Sensor tuning
Structural weakness estimation
Aircraft structure Tunable piezoelectric impedance sensor systems
Air Force DDDAS Program PI Meeting, 01/27/2016
Data driven sensor tuning and structural weakness estimation
Structural weaknesses
Applications
Failure threshold
Time
Structural weakness growth modeling and prognosis
4
First-year Progress: Identifiability
Systematic identifiability study on damage detection using system dynamic responses
Air Force DDDAS Program PI Meeting, 01/27/2016
5
Existing Work on Identifiabiliy in Structural Systems • Liu, P. L. (1995) – show the solution of element properties (e.g., mass and stiffness) to minimize error function of a truss structure is global minimum for sufficient data-identifiability in the sense of optimization
• Franco, G., Betti, R., & Longman, R. W. (2006) – discuss the relationship between identifiability of shear-type linear structural system and the number of actuators and sensors. Provide the minimum requirement in order to identify all system parametersidentifiability in the sense of experimental setups
• Chang, J. D., & Guo, B. Z. (2007) – prove that the density and the flexural rigidity of the Euler-Bernoulli beam in class 𝐶𝐶 4 can be uniquely determined from input and output functionsidentifiability in the sense of continuous beam structure
Limited work has been done in investigating the identifiability of the method Finite Element Model (FEM).
Air Force DDDAS Program PI Meeting, 01/27/2016
6
Finite Element Model and Linear Time Invariant System • Finite Element Model (FEM) is a broadly used numerical tool in structural response analysis. Beam Structure
𝑴𝑴𝒙𝒙̈ + 𝑪𝑪𝒙𝒙̇ + 𝑲𝑲𝑲𝑲 = 𝒇𝒇(𝑡𝑡)
Approximated FEM-based Beam
123
…
n
• A FEM can be formulated into LTI by setting 𝒛𝒛 =
𝑻𝑻
𝒙𝒙𝑻𝑻 , 𝒙𝒙𝑻𝑻̇ , 𝐀𝐀 =
𝟎𝟎 𝟎𝟎 𝐈𝐈 and 𝐂𝐂 = [𝐄𝐄𝐂𝐂 , 𝟎𝟎𝐓𝐓 ], where 𝐄𝐄𝐁𝐁 and 𝐄𝐄𝐂𝐂 are matrices that , 𝐁𝐁 = −𝟏𝟏 −𝟏𝟏 𝐄𝐄 𝐌𝐌 −𝐌𝐌 𝐊𝐊 𝟎𝟎 𝐁𝐁 specify the input and output locations
𝑴𝑴𝒙𝒙̈ + 𝑪𝑪𝒙𝒙̇ + 𝑲𝑲𝑲𝑲 = 𝒇𝒇(𝑡𝑡)
Air Force DDDAS Program PI Meeting, 01/27/2016
𝑑𝑑𝒛𝒛 𝑡𝑡 𝑑𝑑𝑑𝑑
= 𝑨𝑨 𝜽𝜽 𝒛𝒛 𝑡𝑡 + 𝑩𝑩 𝜽𝜽 𝒖𝒖 𝑡𝑡
𝒚𝒚 𝑡𝑡 = 𝑪𝑪 𝜽𝜽 𝒛𝒛 𝑡𝑡
7
Parameter Identifiability in Linear Time Invariant System 𝑼𝑼(𝑠𝑠)
𝒀𝒀(𝑠𝑠)
𝑯𝑯(𝑠𝑠)
𝒀𝒀 𝑠𝑠 = 𝑯𝑯 𝑠𝑠 𝑼𝑼 𝑠𝑠 ,
𝑯𝑯 𝑠𝑠 = 𝑪𝑪(𝜽𝜽) 𝑠𝑠𝑰𝑰 − 𝑨𝑨(𝜽𝜽)
Basic Definition:
−1
𝑩𝑩(𝜽𝜽) + 𝑫𝑫(𝜽𝜽)
Let 𝑨𝑨, 𝑩𝑩, 𝑪𝑪, 𝑫𝑫 𝜽𝜽 be a parametrization of the system matrices 𝑨𝑨, 𝑩𝑩, 𝑪𝑪, 𝑫𝑫 . This parametrization is said to be parameter-identifiable if: 𝑪𝑪 𝜽𝜽𝟏𝟏 𝑠𝑠𝑰𝑰 − 𝑨𝑨 𝜽𝜽𝟏𝟏
= 𝑪𝑪 𝜽𝜽𝟐𝟐 𝑠𝑠𝑰𝑰 − 𝑨𝑨 𝜽𝜽𝟐𝟐
for all 𝑠𝑠 ∈ ℂ, implying 𝜽𝜽𝟏𝟏 = 𝜽𝜽𝟐𝟐 . Air Force DDDAS Program PI Meeting, 01/27/2016
−𝟏𝟏
𝑩𝑩 𝜽𝜽𝟏𝟏 + 𝑫𝑫 𝜽𝜽𝟏𝟏
−𝟏𝟏
𝑩𝑩 𝜽𝜽𝟐𝟐 + 𝑫𝑫 𝜽𝜽𝟐𝟐 8
Problem Setup • Damage Severity Identifiability: For a fixed damage location 𝒑𝒑, whether the change of damage severity 𝜸𝜸 will lead the change of the transfer function matrix 𝑯𝑯 𝑠𝑠 , i.e., 𝑯𝑯| 𝑝𝑝, 𝛾𝛾1 ? = 𝑯𝑯|(𝑝𝑝, 𝛾𝛾2 ) for 𝛾𝛾1 ≠ 𝛾𝛾2 . Log Amplitude
Transfer Function
Frequency
• Damage Location Identifiability:
For different pairs of (𝒑𝒑, 𝜸𝜸), whether or not the transfer function matrix 𝑯𝑯 𝑠𝑠 will be the same, i.e., 𝑯𝑯| 𝑝𝑝𝟏𝟏 , 𝛾𝛾1 = ? 𝑯𝑯|(𝑝𝑝𝟐𝟐 , 𝛾𝛾2 ) for 𝑝𝑝𝟏𝟏 , 𝛾𝛾1 ≠ 𝑝𝑝𝟐𝟐 , 𝛾𝛾𝟐𝟐 . Log Amplitude
Transfer Function
Frequency
Air Force DDDAS Program PI Meeting, 01/27/2016
9
LTI System Formulation Based on FEM: Beam Illustration The transfer function of a beam structure:
𝑯𝑯 𝑠𝑠 = 𝑬𝑬𝒄𝒄
−𝑠𝑠 𝟐𝟐 𝑴𝑴 +
Element
𝑣𝑣1
𝑴𝑴𝒊𝒊 = 𝜌𝜌
𝑠𝑠𝑠𝑠𝑠𝑠
11𝑙𝑙 2 210 𝑙𝑙 3 105
9𝑙𝑙 70 13𝑙𝑙 2 420 13𝑙𝑙 35
13𝑙𝑙 2 420 𝑙𝑙 2 − 140 11𝑙𝑙 2 − 210 𝑙𝑙 3 105
𝑬𝑬𝑩𝑩 ,
𝜽𝜽 = (𝒑𝒑, 𝜸𝜸)
Damage location
Node 2
−
𝑲𝑲𝒊𝒊 = 𝐸𝐸𝐸𝐸
12 𝑙𝑙 3
𝑠𝑠𝑠𝑠𝑠𝑠
Air Force DDDAS Program PI Meeting, 01/27/2016
6 𝑙𝑙 2 4 𝑙𝑙
Damage severity
Joint nodes
𝑑𝑑2
Node 1
13𝑙𝑙 35
𝑲𝑲(𝜽𝜽)
−𝟏𝟏
12 𝑙𝑙 3 6 − 2 𝑙𝑙 12 𝑙𝑙 2
−
6 𝑙𝑙 2 2 𝑙𝑙 6 − 2 𝑙𝑙 4 𝑙𝑙
1 𝑴𝑴 𝑜𝑜𝑜𝑜 𝑲𝑲 =
4×4
2
𝟎𝟎
⋱
𝟎𝟎
10
Properties of Transfer Function 𝑯𝑯(𝒔𝒔)
Transfer function 𝑯𝑯 𝑠𝑠 is the inverse of a block tri-diagonal matrix −𝑠𝑠 𝟐𝟐 𝑴𝑴 + 𝑲𝑲 𝜽𝜽 𝑹𝑹1 −𝑸𝑸𝑻𝑻2 −𝑸𝑸2 𝑹𝑹2 −𝑸𝑸𝑻𝑻3 𝟐𝟐 ⋱ Define 𝑽𝑽 = −𝑠𝑠 𝑴𝑴 + 𝑲𝑲 𝜽𝜽 = ⋱ ⋱ 𝑻𝑻 𝑹𝑹 −𝑸𝑸𝒏𝒏 𝒏𝒏 −𝑸𝑸𝒏𝒏+1 −𝑸𝑸𝒏𝒏+1 𝑹𝑹𝒏𝒏+1
Define recursive functions 𝜟𝜟 and 𝜮𝜮 For single damage at the same location in collocated setup
𝑽𝑽 =
Input
𝑗𝑗 𝑡𝑡𝑡
𝜟𝜟1 = 𝑹𝑹1 , � 𝑻𝑻 𝜟𝜟𝒊𝒊 = 𝑹𝑹𝒊𝒊 − 𝑸𝑸𝒊𝒊 𝜟𝜟−1 𝒊𝒊−1 𝑸𝑸𝒊𝒊 ,
𝑹𝑹𝒋𝒋 −𝑸𝑸𝑻𝑻𝒋𝒋+𝟏𝟏 −𝑸𝑸𝒋𝒋+𝟏𝟏 𝑹𝑹𝒋𝒋+𝟏𝟏 𝐻𝐻𝑘𝑘𝑘𝑘 ≠ 𝐻𝐻𝑘𝑘𝑘𝑘
= 𝑽𝑽−1
𝑖𝑖 < 𝑗𝑗 𝜴𝜴𝒊𝒊 = 𝜟𝜟𝒊𝒊 , 𝑖𝑖 > 𝑗𝑗 , �𝜴𝜴𝒊𝒊 = 𝜮𝜮𝒊𝒊 , 𝜮𝜮𝒏𝒏+1 = 𝑹𝑹𝒏𝒏+1 , 𝜴𝜴𝒋𝒋 = 𝜟𝜟𝒋𝒋 + 𝜮𝜮𝒋𝒋 − 𝑹𝑹𝒋𝒋 � 𝑻𝑻 −1 𝜮𝜮𝒊𝒊 = 𝑹𝑹𝒊𝒊 − 𝑸𝑸𝒊𝒊+1 𝜮𝜮𝒊𝒊+1 𝑸𝑸𝒊𝒊+1 . For damages at different locations in collocated setup
Output damage
⋱
, then 𝑯𝑯 𝑠𝑠
𝑯𝑯𝟏𝟏𝒏𝒏 ⋯ 𝑯𝑯𝟏𝟏𝟏𝟏 𝑯𝑯𝟏𝟏𝟏𝟏 𝑯𝑯𝟐𝟐𝟐𝟐 𝑯𝑯𝟐𝟐𝟐𝟐 = ⋮ ⋱ ⋱ ⋱ 𝑯𝑯𝒏𝒏𝒏𝒏 𝑯𝑯 ⋮ 𝒏𝒏𝒏𝒏+𝟏𝟏 𝐇𝐇𝐧𝐧+𝟏𝟏𝟏𝟏 𝑯𝑯 𝑯𝑯𝒏𝒏+𝟏𝟏𝟏𝟏 𝒏𝒏+𝟏𝟏𝐧𝐧+𝟏𝟏
Air Force DDDAS Program PI Meeting, 01/27/2016
⋱
𝑽𝑽 =
𝑯𝑯𝒋𝒋𝒋𝒋 = 𝜴𝜴−𝟏𝟏 𝒋𝒋
𝑻𝑻 −𝟏𝟏 𝑻𝑻 −𝟏𝟏 𝑯𝑯𝒋𝒋𝒋𝒋−𝒍𝒍 = 𝚫𝚫−𝟏𝟏 𝐣𝐣−𝐥𝐥 𝑸𝑸𝒋𝒋−𝒍𝒍+𝟏𝟏 ⋯ 𝚫𝚫 𝐣𝐣−𝟏𝟏 𝑸𝑸𝒋𝒋 𝜴𝜴𝒋𝒋 −𝟏𝟏 −𝟏𝟏 𝑯𝑯𝒋𝒋𝒋𝒋+𝒍𝒍 = 𝚺𝚺𝒋𝒋+𝒍𝒍 𝑸𝑸𝒋𝒋+𝒍𝒍 ⋯ 𝜮𝜮−𝟏𝟏 𝒋𝒋+𝟏𝟏 𝑸𝑸𝒋𝒋+𝟏𝟏 𝜴𝜴𝒋𝒋
𝑗𝑗 𝑡𝑡𝑡
⋱
𝑙𝑙 𝑡𝑡𝑡
−𝑸𝑸𝑻𝑻𝒋𝒋+𝟏𝟏
𝑹𝑹𝒋𝒋 −𝑸𝑸𝒋𝒋+𝟏𝟏 𝑹𝑹𝒋𝒋+𝟏𝟏
⋱
𝑽𝑽 =
𝑇𝑇 𝑅𝑅𝑙𝑙 −𝑄𝑄𝑙𝑙+1 −𝑄𝑄𝑙𝑙+1 𝑅𝑅𝑙𝑙+1
Untraceable!
? 𝐻𝐻𝑘𝑘𝑘𝑘 𝐻𝐻𝑘𝑘𝑘𝑘 =
11
⋱
⋱
Identifiability of Damage Severity
Input Output
damage
Lemma 𝑯𝑯𝑗𝑗 |(𝑝𝑝, 𝛾𝛾1 ) ≠ 𝑯𝑯𝑗𝑗 |(𝑝𝑝, 𝛾𝛾2 ) if 𝛾𝛾1 ≠ 𝛾𝛾2 for 𝑝𝑝 = 1,2, … 𝑛𝑛 and 𝑗𝑗 = 1,2, , … 𝑛𝑛 + 1, where 𝑯𝑯𝑗𝑗 |(𝑝𝑝, 𝛾𝛾1 ) represents the matrix 𝑯𝑯𝑗𝑗 given parameters (𝑝𝑝, 𝛾𝛾1 ). Input 𝑛𝑛 2 Output
Corollary 𝑯𝑯𝑛𝑛 |(𝑝𝑝, 𝛾𝛾) = 𝑯𝑯𝑛𝑛 |(𝑛𝑛 + 1 − 𝑝𝑝, 𝛾𝛾) for an even number 𝒏𝒏 2
2
Air Force DDDAS Program PI Meeting, 01/27/2016
12
Identifiability of Damage Location Input 𝒑𝒑∗ , 𝜸𝜸∗ , 𝓛𝓛 and 𝒕𝒕
• Transfer Function can also be expressed as 𝐻𝐻 𝑠𝑠 =
𝑐𝑐1 s−𝑧𝑧1 )(s−𝑧𝑧2 )(s−𝑧𝑧3 )⋯(s−𝑧𝑧𝑚𝑚−1 )(s−𝑧𝑧𝑚𝑚 2 )(s−𝑓𝑓 2 𝑐𝑐2 s−𝑓𝑓12 )(s−𝑓𝑓22 )(s−𝑓𝑓32 )⋯(s−𝑓𝑓𝑛𝑛−1 𝑛𝑛
Update 𝑲𝑲
,
𝑓𝑓𝑖𝑖2 = 𝜆𝜆𝑖𝑖 , eigenvalue of 𝑴𝑴−𝟏𝟏 𝑲𝑲 𝒑𝒑, 𝜸𝜸
• Sufficient condition: the difference in natural frequencies lead to different transfer functions
Compute 𝒇𝒇𝒊𝒊
𝒑𝒑∗ ,𝜸𝜸∗
Compute 𝑆𝑆𝑘𝑘 = min ∑𝑖𝑖∈ℒ for all 𝑘𝑘
𝜸𝜸
Sort 𝑆𝑆𝑘𝑘 as 𝑺𝑺
𝟏𝟏
𝒑𝒑 ,𝜸𝜸 𝒑𝒑 ,𝜸𝜸 𝑓𝑓𝑖𝑖 ∗ ∗ −𝑓𝑓𝑖𝑖 𝒌𝒌 𝒑𝒑 ,𝜸𝜸 𝑓𝑓𝑖𝑖 ∗ ∗
≤ 𝑺𝑺
If 𝑺𝑺 𝟐𝟐 − 𝑺𝑺 𝟏𝟏 > 𝒕𝒕
𝟐𝟐
≤⋯
Yes
𝒑𝒑∗ , 𝜸𝜸∗ is identifiable
No
Air Force DDDAS Program PI Meeting, 01/27/2016
𝒑𝒑∗ , 𝜸𝜸∗ is NOT identifiable
13
Identifiability of Damage Location: Example
𝑆𝑆𝑘𝑘 = min � 𝜸𝜸
0.5
𝑖𝑖∈ℒ
𝑓𝑓𝑖𝑖
𝒑𝒑∗ ,𝜸𝜸∗
𝑓𝑓𝑖𝑖
− 𝑓𝑓𝑖𝑖
𝒑𝒑∗ ,𝜸𝜸∗
𝒑𝒑𝒌𝒌 ,𝜸𝜸
𝓛𝓛 = {𝟏𝟏: 𝟓𝟓}
0.8
𝑆𝑆𝑘𝑘 0
e.g. 𝑝𝑝∗ , 𝛾𝛾∗ = 23,0.4
𝑆𝑆𝑘𝑘 1
𝑘𝑘
• The minimum value of 𝑆𝑆𝑘𝑘 = 0 when 𝑘𝑘 = 23 and 38 for a beam with 60 elements in a fixedfixed experimental setup. • The parameter 𝑝𝑝∗ , 𝛾𝛾∗ = 23,0.4 can only be identified in pairs due to symmetry. Air Force DDDAS Program PI Meeting, 01/27/2016
60
0
1
𝑘𝑘
60
• The minimum value of 𝑆𝑆𝑘𝑘 = 0 when 𝑘𝑘 = 23 only for a beam with 60 elements in a fixed-free experimental setup. • The parameter 𝑝𝑝∗ , 𝛾𝛾∗ = 23,0.4 can be uniquely identified.
14
First-year Progress: Adaptive Sensor & Damage Identification
Structural damage identification using piezoelectric impedance active sensing
Air Force DDDAS Program PI Meeting, 01/27/2016
15
Piezoelectric Impedance Sensory System: Concept Impedance sensor and measurement: two-way electromechanical coupling Structural impedance coupled with transducer impedance Resistor
Transducer
200 Real Part of Impedance (Ohm)
Impedance with noticeable changes
160
120
80
40 5
10
15 20 Excitation frequency (kHz)
Transducer
25
30
35
Circuitry elements
New idea: impedance approach + piezoelectric circuitry
Integrate circuitry elements to enhance dynamics & resonant effects Fundamental advantage: circuitry can amplify signals and enrich measurements
Air Force DDDAS Program PI Meeting, 01/27/2016
Mathematical Formulation Finite element-based formulation Transducer
+ Cq + Kq + K 12Q = Mq 0 + RQ + K Q + K T q = LQ Vin c 12 • M,C,K – mass matrix, damping matrix, stiffness matrix • 𝐊𝐊 𝟏𝟏𝟏𝟏 – coupling vector
Circuitry elements
• q – displacement vector • Q – electrical charge
• L, R – circuitry inductance and resistance
Under harmonic excitation, the impedance measurement is Z=
Vin 1 T ( K − Mω 2 + Cωi ) −1 K 12 ] = [( −ω 2 L + iω R + K c ) − K 12 Q ωi
Structural damage impedance change: How to identify damage under impedance changes (at various frequencies)?
Damage Quantification and Influence to Piezoelectric Circuitry Impedance Assumption: damage is modeled as elemental stiffness reduction • 𝐋𝐋 – indicates how the elemental stiffness
= Kd
m
∑ L [K j =1
T j
ej
(1 − D j )]L j
𝐣𝐣
matrices are assembled • 𝐷𝐷𝑗𝑗 – damage index (percentage of stiffness reduction)
Impedance of damaged structure
m Vin 1 T Z d = = [ Z c − K 12 ( ∑ LTj [K ej (1 − D j )]L j − Mω 2 + Cωi ) −1 K 12 ] Q ωi j =1
• Taylor series expansion • Small damage - Disregard the higher order terms
∂Z Z d ( D1 , D2 ,...Dm ) ≈ Z ( D1 , D2 ,...Dm = 0) + ∑ d |D1 , D2 ,... Dm = 0 D j j =1 ∂D j m
Damage Sensitivity Matrix Impedance change expression: m
1 T K 12 ( K − Mω 2 + Cωi ) − ( LTj K ej L j )( K − Mω 2 + Cωi ) − K 12 }D j ωi
∆Z = Zd − Z (D = 0) = ∑{ j =1
Compare impedances at n frequency points
∆Z (ω1 ) D1 ∆Z (ωi ) = [ S nm ]n×m Di ∆Z (ω ) D n n×1 m m×1 Sensitivity matrix entries: 1 T 2 − T 2 − S= K K − M ω + C ω i L K L K − M ω + C ω i K12 ( ) ( )( ) nm 12 n n m em m n n ωn i Grand Challenge: Problem is severely under-determined, as the number of finite elements (m) in practical structure is very large, much greater than the number of impedance measurements (n)!
Damage Identification Approach – Exclusion Solution to challenge: pre-screening of possible regions Consider the most likely case of single damage ∆Z (ω1 ) D1 ∆Z (ωi ) = S n×m Di ∆Z (ω ) D n n×1 m m×1
∆Z (ω1 ) 0 ∆Z (ωi ) = S n×m Di = ∆Z (ω ) 0 n n×1 m×1
S1,i S D j ,i i S n ,i
Impedance change vector in theory similar to one column of the sensitivity matrix(𝐒𝐒:,𝑖𝑖 ) (multiplied with a scale factor(𝐷𝐷𝑖𝑖 )) The ratio of the two vectors is the damage severity
Eliminate impossible columns/elements Consider each column of S matrix D pi = mean(ΔZ. / S:,i ) All elements with D pi < 0 or D pi > 1 can be directly excluded This approach is applicable regardless the mesh size; Can be extended to ‘super-elements’ or ‘regions’.
Damage Identification Approach – Similarity Comparison Damage prediction: First estimate impedance change correlation ∆Zei = D pi . × S:,i Then compute and check similarity between estimated impedance change and the actual measurement SI i =∆ ( Z ei ⋅ ΔZ) /( ΔZ ⋅ ΔZ)
The closer 𝑆𝑆𝑆𝑆𝑖𝑖 to 1, the more likely damage occurs at the i-th element and 𝐷𝐷𝑝𝑝𝑖𝑖 is the severity. Define a similarity function for the i-the element Pi ( SI i ,1) =
1 / SI i − 1
sum(1 / SI i − 1)
Higher Pi indicates more probable damage location
New algorithm features: Avoid inversion of under-determined matrix equation Naturally fit the probabilistic and data-driven analysis framework
Case Illustration - Setup Cantilever plate
Piezoelectric transducer • Size: 15 mm×19.05mm×1.4 mm • Piezoelectric constant : 𝑑𝑑31 : −140 × 10−12 𝑚𝑚/𝑉𝑉 𝑔𝑔33 :25 × 10−3 𝑉𝑉𝑉𝑉/𝑁𝑁
• Size: 561 mm×19.05mm×4.7625mm • Young’s modulus: 6.89 GPa • Density: 2700 𝐾𝐾𝐾𝐾/𝑚𝑚3
4.7625mm
561 mm
Piezoelectric transducer 19.05mm
180 mm
15 mm Signal analyzer (Agilent 35670A)
Vin Vout
PZT Cantilever plate
Small resistor Rs
STRUCTURE AND SYSTEM DYNAMICS LAB
Case Illustration – Model Calibration Finite element model with convergence study and model updating To facilitate high detection sensitivity, high-frequency impedance is needed, requiring large number of elements. The plate is discretized into 12500 elements, grouped into 250 segments. FEA Model updating of healthy baseline is EXP (Hz) Original Updated (Hz) (Hz) conducted. Boundary conditions and 12.20 12.46 12.08 material properties are updated to 74.00 76.50 74.29 match FE predictions with experiments. 209.50 213.95 208.22 S151
S152
S153
S154
S225
S76
S77
S78
S79
S150
S1
S2
S3
S4
S75
414.70 682.80 1017.10 1893.18 2428.44 3033.20 3703.09
424.47 699.14 1040.69 1936.73 2482.46 3104.30 3781.89
413.34 681.33 1015.93 1893.58 2429.08 3034.05 3704.05
Impedance Measurements A small mass is attached onto the plate to emulate damage (Segment 103; 0.5% equivalent stiffness change) Impedances are measured around 14th, 16th , 20th resonance frequencies
(1) 14th resonance frequency Air Force DDDAS Program PI Meeting, 01/27/2016
(2) 16th resonance frequency
(3) 20th resonance frequency
Identification Result Damage Segment # Similarity severity Use measurements around the and 103 0.0786 0.005 20th resonances: 214 0.0622 0.003 64 0.0622 0.003 55 0.0522 0.002 205 0.0522 0.002 121 0.0521 0.005 161 0.0482 0.003 11 0.0482 0.003 139 0.0482 0.003 169 0.0464 0.002 19 0.0464 0.002 86 0.0429 0.003 28 0.0392 0.005 178 0.0392 0.005 94 0.0379 0.002 196 0.0364 0.005 Location 46 0.0364 0.005 Multiple possible damage locations are predicted; 130 0.0353 0.002 Actual damage location has the highest likelihood; 131 0.0090 0.003 56 0.0070 0.003 Damage severity prediction is accurate.
Similarity
16th
Air Force DDDAS Program PI Meeting, 01/27/2016
Identification Result
Similarity
Use measurements around one more resonance, the 14th:
Segment # Similarity 103 0.088 121 0.084 28 0.084 178 0.084 196 0.077 46 0.077 140 0.012 85 0.012 102 0.012 122 0.011
Damage severity 0.005 0.005 0.005 0.005 0.005 0.005 0.010 0.009 0.004 0.005
Location
Damage prediction result is improved: number of possible damage locations reduced drastically, decreasing the false-alarm possibility Perfectly suited for probabilistic-based inference analysis Air Force DDDAS Program PI Meeting, 01/27/2016
Probabilistic Inference – Another Study Another approach to avoiding direct matrix inversion Update the probability estimate of interested parameters for a hypothesis as additional information is acquired – take the effect of uncertainty into account Can be fully represented by Bayes’ rule p (θ | D) =
Incorporate data
p (D | θ) p (θ) ∫ p (D | θ) p (θ)dθ
Prior PDF
Posterior PDF
:prior distribution / hypothesis p (D | θ) :likelihood function / evidence p (θ | D) :posterior distribution p (θ)
θ: model parameters – Damage location & severity;
D: data/observations – piezoelectric impedance measurements; Likelihood function 𝑝𝑝(𝐷𝐷|𝜃𝜃): p ( D | θ ) =
1
2πε
e
Air Force DDDAS Program PI Meeting, 01/27/2016
− ( D − D ( θ ))2 2ε
𝐷𝐷:measured impedance � 𝐷𝐷:associated impedance parameterized by 𝜃𝜃�
Bayesian Inference – Challenge and Possible Solution Bayesian inference procedure Generate a full set of parameter samples subject to assumed distribution Call finite element analysis of each updated model (w.r.t each parameter sample) to derive the associated impedance prediction Use Bayesian inference to predict the parameter distribution
Challenge:
Size of parameter sample is huge — a lot of possible damage locations and severity levels Repeated FEA under large parameter space— huge Computational cost
Possible Solution:
Only perform finite element analysis under a reduced sample size Gaussian Process—enrich the FEA output under full set of parameters
Air Force DDDAS Program PI Meeting, 01/27/2016
Case Illustration – Problem Setup Bayesian inference with Gaussian Process
117 mm
PZT 134 mm
plate
PZT
circuitry
Clamped-clamped aluminum plate Size :610 mm×508 mm×4.826mm Young’s modulus: 73 GPa Density: 2780 𝐾𝐾𝐾𝐾/𝑚𝑚3 PZT transducer (PZT-4): Size :44.88 mm×44.88 mm×2.8mm Piezoelectric constant 𝑑𝑑31 : −140 × 10−12 𝑚𝑚/𝑉𝑉 𝑔𝑔33 :25 × 10−3 𝑉𝑉𝑉𝑉/𝑁𝑁
Finite model and parametrization: 2688 elements 9 severity levels
Computational costs: Full-size: 2688 X 9 rounds of finite element dynamic analysis; 145 hours Reduced size: 165 X 3 rounds; 3 hours Air Force DDDAS Program PI Meeting, 01/27/2016
Case Illustration – Result When impedance measurements around two resonances are used Location
Severity
Probability
1
39
6
3.1551
2
143
5
3.1550
3
24
3
3.1550
4
49
2
3.1547
5
129
5
3.1546
6
58
2
3.1546
7
72
5
3.1541
Posterior Distribution 4 3 2 1 0 6 4 2 0
5
0
15
10
20
25
30
35
Severity 5 4 3 2 1 0 5 4 3 2 1
0
10
20
30
Severity 10 Air Force DDDAS Program PI Meeting, 01/27/2016
40
The results indicate to a possible area The actual damage parameter has the 7th highest possibility
Case Illustration – Result When impedance measurements around a third resonance is added
Location Severity
Posterior Distribution 10 8 6 4 2 0 5 4 3 2 1
0
10
20
30
40
Severity 5
Probability
1
72
5
10.18
2
36
2
10.12
3
159
2
9.99
4
15
3
9.93
5
105
5
9.88
As more information is obtained and applied, we can get better result – the actual damage parameter now has the highest possibility The computational cost is still relatively high On-going effort: To combine pre-screening with Bayesian inference! Air Force DDDAS Program PI Meeting, 01/27/2016
Summary and Ongoing Work Brief Summary on Identifiability Issue Theoretically proved the identifiability of damage at a fixed location in the FEM-based beam structure. A numerical algorithm is proposed to numerically check and validate the location identifiability.
Brief Summary on impedance-based damage identification
A promising damage pre-screening scheme is developed.
Probabilistic inference framework has been established.
On-going work on the estimation of damage changing trend Existing studies focus on static detection, that is, one time estimation of the current health status with available data. Lack of the ability to inference the future trend of the health status. A Dynamic Data Driven Based Hierarchical Bayesian Degradation Model will be developed. As data are accumulated, the new method can recover the true trend quickly.
Air Force DDDAS Program PI Meeting, 01/27/2016
32
Publications •
•
• • •
Q. Shuai, K. Zhou, and J. Tang , 2015, Structural Damage Identification Using Piezoelectric Impedance and Bayesian Inference, Proceedings of SPIE, Smart Structures / NDE, V9435, 2015. Yuhang Liu , Shiyu Zhou, Jiong Tang, 2015, Identifiability Analysis of Finite Element Models for Vibration-Response Based Structure Damage Detection, ASME Transactions, Journal of Dynamic Systems, Measurement, and Control, in revision Hou, Y., Wu, J., and Chen, Y. 2015, “Online Steady State Detection Based on Rao-Blackwellized Sequential Monte Carlo”. IEEE Transactions on Automation Science and Engineering, in revision K. Zhou, Q. Shuai, and J. Tang, Damage identification using piezoelectric impedance and intelligent inference, to be submitted to Journal of Intelligent Material Systems and Structures. Liu, Y., Zhou, S., and Tang, J, Dynamic data-driven Bayesian structural failure prognosis, to be submitted to IIE Transactions.
All publications acknowledge the support by AFOSR DDDAS Program under grant FA9550-14-1-0384. Air Force DDDAS Program PI Meeting, 01/27/2016
33
Progressive Fault Identification and Prognosis in Aircraft Structure Based on Dynamic Data Driven Adaptive Sensing and Simulation
Shiyu Zhou University of Wisconsin-Madison
Jiong Tang University of Connecticut
Yong Chen University of Iowa Air Force DDDAS Program PI Meeting, 01/27/2016
1
Outline • Introduction and project objectives • Research results accomplished the first project year Identifiability analysis for finite element based dynamic damage detection Establishment of piezoelectric impedance-based damage identification algorithms and test-beds
• Summary and on-going efforts
Air Force DDDAS Program PI Meeting, 01/27/2016
2
Introduction • Structural health monitoring and management (SHM2) is of critical importance for air force applications • Dynamic responses of an aerospace structure contain rich information • Limitations of the current dynamic response based SHM2 systems – Sensory system is fixed and re-active – Applications modeling is either based on mechanistic principles or simple data driven models using historical data – Sensor-structure interaction dynamics are not exploited; Uncertainties not fully addressed
Air Force DDDAS Program PI Meeting, 01/27/2016
3
Research objectives Overarching Goal: To create a new methodology of progressive structural fault identification and prognosis for Air Force applications based on the framework of dynamic data-driven applications systems (DDDAS). Applications measurements systems and methods Piezoelectric transducer
Vi
i
Applications modeling
Mathematical and statistical algorithms Failure time prognosis
Data collection
Data-driven Sensor tuning
Structural weakness estimation
Aircraft structure Tunable piezoelectric impedance sensor systems
Air Force DDDAS Program PI Meeting, 01/27/2016
Data driven sensor tuning and structural weakness estimation
Structural weaknesses
Applications
Failure threshold
Time
Structural weakness growth modeling and prognosis
4
First-year Progress: Identifiability
Systematic identifiability study on damage detection using system dynamic responses
Air Force DDDAS Program PI Meeting, 01/27/2016
5
Existing Work on Identifiabiliy in Structural Systems • Liu, P. L. (1995) – show the solution of element properties (e.g., mass and stiffness) to minimize error function of a truss structure is global minimum for sufficient data-identifiability in the sense of optimization
• Franco, G., Betti, R., & Longman, R. W. (2006) – discuss the relationship between identifiability of shear-type linear structural system and the number of actuators and sensors. Provide the minimum requirement in order to identify all system parametersidentifiability in the sense of experimental setups
• Chang, J. D., & Guo, B. Z. (2007) – prove that the density and the flexural rigidity of the Euler-Bernoulli beam in class 𝐶𝐶 4 can be uniquely determined from input and output functionsidentifiability in the sense of continuous beam structure
Limited work has been done in investigating the identifiability of the method Finite Element Model (FEM).
Air Force DDDAS Program PI Meeting, 01/27/2016
6
Finite Element Model and Linear Time Invariant System • Finite Element Model (FEM) is a broadly used numerical tool in structural response analysis. Beam Structure
𝑴𝑴𝒙𝒙̈ + 𝑪𝑪𝒙𝒙̇ + 𝑲𝑲𝑲𝑲 = 𝒇𝒇(𝑡𝑡)
Approximated FEM-based Beam
123
…
n
• A FEM can be formulated into LTI by setting 𝒛𝒛 =
𝑻𝑻
𝒙𝒙𝑻𝑻 , 𝒙𝒙𝑻𝑻̇ , 𝐀𝐀 =
𝟎𝟎 𝟎𝟎 𝐈𝐈 and 𝐂𝐂 = [𝐄𝐄𝐂𝐂 , 𝟎𝟎𝐓𝐓 ], where 𝐄𝐄𝐁𝐁 and 𝐄𝐄𝐂𝐂 are matrices that , 𝐁𝐁 = −𝟏𝟏 −𝟏𝟏 𝐄𝐄 𝐌𝐌 −𝐌𝐌 𝐊𝐊 𝟎𝟎 𝐁𝐁 specify the input and output locations
𝑴𝑴𝒙𝒙̈ + 𝑪𝑪𝒙𝒙̇ + 𝑲𝑲𝑲𝑲 = 𝒇𝒇(𝑡𝑡)
Air Force DDDAS Program PI Meeting, 01/27/2016
𝑑𝑑𝒛𝒛 𝑡𝑡 𝑑𝑑𝑑𝑑
= 𝑨𝑨 𝜽𝜽 𝒛𝒛 𝑡𝑡 + 𝑩𝑩 𝜽𝜽 𝒖𝒖 𝑡𝑡
𝒚𝒚 𝑡𝑡 = 𝑪𝑪 𝜽𝜽 𝒛𝒛 𝑡𝑡
7
Parameter Identifiability in Linear Time Invariant System 𝑼𝑼(𝑠𝑠)
𝒀𝒀(𝑠𝑠)
𝑯𝑯(𝑠𝑠)
𝒀𝒀 𝑠𝑠 = 𝑯𝑯 𝑠𝑠 𝑼𝑼 𝑠𝑠 ,
𝑯𝑯 𝑠𝑠 = 𝑪𝑪(𝜽𝜽) 𝑠𝑠𝑰𝑰 − 𝑨𝑨(𝜽𝜽)
Basic Definition:
−1
𝑩𝑩(𝜽𝜽) + 𝑫𝑫(𝜽𝜽)
Let 𝑨𝑨, 𝑩𝑩, 𝑪𝑪, 𝑫𝑫 𝜽𝜽 be a parametrization of the system matrices 𝑨𝑨, 𝑩𝑩, 𝑪𝑪, 𝑫𝑫 . This parametrization is said to be parameter-identifiable if: 𝑪𝑪 𝜽𝜽𝟏𝟏 𝑠𝑠𝑰𝑰 − 𝑨𝑨 𝜽𝜽𝟏𝟏
= 𝑪𝑪 𝜽𝜽𝟐𝟐 𝑠𝑠𝑰𝑰 − 𝑨𝑨 𝜽𝜽𝟐𝟐
for all 𝑠𝑠 ∈ ℂ, implying 𝜽𝜽𝟏𝟏 = 𝜽𝜽𝟐𝟐 . Air Force DDDAS Program PI Meeting, 01/27/2016
−𝟏𝟏
𝑩𝑩 𝜽𝜽𝟏𝟏 + 𝑫𝑫 𝜽𝜽𝟏𝟏
−𝟏𝟏
𝑩𝑩 𝜽𝜽𝟐𝟐 + 𝑫𝑫 𝜽𝜽𝟐𝟐 8
Problem Setup • Damage Severity Identifiability: For a fixed damage location 𝒑𝒑, whether the change of damage severity 𝜸𝜸 will lead the change of the transfer function matrix 𝑯𝑯 𝑠𝑠 , i.e., 𝑯𝑯| 𝑝𝑝, 𝛾𝛾1 ? = 𝑯𝑯|(𝑝𝑝, 𝛾𝛾2 ) for 𝛾𝛾1 ≠ 𝛾𝛾2 . Log Amplitude
Transfer Function
Frequency
• Damage Location Identifiability:
For different pairs of (𝒑𝒑, 𝜸𝜸), whether or not the transfer function matrix 𝑯𝑯 𝑠𝑠 will be the same, i.e., 𝑯𝑯| 𝑝𝑝𝟏𝟏 , 𝛾𝛾1 = ? 𝑯𝑯|(𝑝𝑝𝟐𝟐 , 𝛾𝛾2 ) for 𝑝𝑝𝟏𝟏 , 𝛾𝛾1 ≠ 𝑝𝑝𝟐𝟐 , 𝛾𝛾𝟐𝟐 . Log Amplitude
Transfer Function
Frequency
Air Force DDDAS Program PI Meeting, 01/27/2016
9
LTI System Formulation Based on FEM: Beam Illustration The transfer function of a beam structure:
𝑯𝑯 𝑠𝑠 = 𝑬𝑬𝒄𝒄
−𝑠𝑠 𝟐𝟐 𝑴𝑴 +
Element
𝑣𝑣1
𝑴𝑴𝒊𝒊 = 𝜌𝜌
𝑠𝑠𝑠𝑠𝑠𝑠
11𝑙𝑙 2 210 𝑙𝑙 3 105
9𝑙𝑙 70 13𝑙𝑙 2 420 13𝑙𝑙 35
13𝑙𝑙 2 420 𝑙𝑙 2 − 140 11𝑙𝑙 2 − 210 𝑙𝑙 3 105
𝑬𝑬𝑩𝑩 ,
𝜽𝜽 = (𝒑𝒑, 𝜸𝜸)
Damage location
Node 2
−
𝑲𝑲𝒊𝒊 = 𝐸𝐸𝐸𝐸
12 𝑙𝑙 3
𝑠𝑠𝑠𝑠𝑠𝑠
Air Force DDDAS Program PI Meeting, 01/27/2016
6 𝑙𝑙 2 4 𝑙𝑙
Damage severity
Joint nodes
𝑑𝑑2
Node 1
13𝑙𝑙 35
𝑲𝑲(𝜽𝜽)
−𝟏𝟏
12 𝑙𝑙 3 6 − 2 𝑙𝑙 12 𝑙𝑙 2
−
6 𝑙𝑙 2 2 𝑙𝑙 6 − 2 𝑙𝑙 4 𝑙𝑙
1 𝑴𝑴 𝑜𝑜𝑜𝑜 𝑲𝑲 =
4×4
2
𝟎𝟎
⋱
𝟎𝟎
10
Properties of Transfer Function 𝑯𝑯(𝒔𝒔)
Transfer function 𝑯𝑯 𝑠𝑠 is the inverse of a block tri-diagonal matrix −𝑠𝑠 𝟐𝟐 𝑴𝑴 + 𝑲𝑲 𝜽𝜽 𝑹𝑹1 −𝑸𝑸𝑻𝑻2 −𝑸𝑸2 𝑹𝑹2 −𝑸𝑸𝑻𝑻3 𝟐𝟐 ⋱ Define 𝑽𝑽 = −𝑠𝑠 𝑴𝑴 + 𝑲𝑲 𝜽𝜽 = ⋱ ⋱ 𝑻𝑻 𝑹𝑹 −𝑸𝑸𝒏𝒏 𝒏𝒏 −𝑸𝑸𝒏𝒏+1 −𝑸𝑸𝒏𝒏+1 𝑹𝑹𝒏𝒏+1
Define recursive functions 𝜟𝜟 and 𝜮𝜮 For single damage at the same location in collocated setup
𝑽𝑽 =
Input
𝑗𝑗 𝑡𝑡𝑡
𝜟𝜟1 = 𝑹𝑹1 , � 𝑻𝑻 𝜟𝜟𝒊𝒊 = 𝑹𝑹𝒊𝒊 − 𝑸𝑸𝒊𝒊 𝜟𝜟−1 𝒊𝒊−1 𝑸𝑸𝒊𝒊 ,
𝑹𝑹𝒋𝒋 −𝑸𝑸𝑻𝑻𝒋𝒋+𝟏𝟏 −𝑸𝑸𝒋𝒋+𝟏𝟏 𝑹𝑹𝒋𝒋+𝟏𝟏 𝐻𝐻𝑘𝑘𝑘𝑘 ≠ 𝐻𝐻𝑘𝑘𝑘𝑘
= 𝑽𝑽−1
𝑖𝑖 < 𝑗𝑗 𝜴𝜴𝒊𝒊 = 𝜟𝜟𝒊𝒊 , 𝑖𝑖 > 𝑗𝑗 , �𝜴𝜴𝒊𝒊 = 𝜮𝜮𝒊𝒊 , 𝜮𝜮𝒏𝒏+1 = 𝑹𝑹𝒏𝒏+1 , 𝜴𝜴𝒋𝒋 = 𝜟𝜟𝒋𝒋 + 𝜮𝜮𝒋𝒋 − 𝑹𝑹𝒋𝒋 � 𝑻𝑻 −1 𝜮𝜮𝒊𝒊 = 𝑹𝑹𝒊𝒊 − 𝑸𝑸𝒊𝒊+1 𝜮𝜮𝒊𝒊+1 𝑸𝑸𝒊𝒊+1 . For damages at different locations in collocated setup
Output damage
⋱
, then 𝑯𝑯 𝑠𝑠
𝑯𝑯𝟏𝟏𝒏𝒏 ⋯ 𝑯𝑯𝟏𝟏𝟏𝟏 𝑯𝑯𝟏𝟏𝟏𝟏 𝑯𝑯𝟐𝟐𝟐𝟐 𝑯𝑯𝟐𝟐𝟐𝟐 = ⋮ ⋱ ⋱ ⋱ 𝑯𝑯𝒏𝒏𝒏𝒏 𝑯𝑯 ⋮ 𝒏𝒏𝒏𝒏+𝟏𝟏 𝐇𝐇𝐧𝐧+𝟏𝟏𝟏𝟏 𝑯𝑯 𝑯𝑯𝒏𝒏+𝟏𝟏𝟏𝟏 𝒏𝒏+𝟏𝟏𝐧𝐧+𝟏𝟏
Air Force DDDAS Program PI Meeting, 01/27/2016
⋱
𝑽𝑽 =
𝑯𝑯𝒋𝒋𝒋𝒋 = 𝜴𝜴−𝟏𝟏 𝒋𝒋
𝑻𝑻 −𝟏𝟏 𝑻𝑻 −𝟏𝟏 𝑯𝑯𝒋𝒋𝒋𝒋−𝒍𝒍 = 𝚫𝚫−𝟏𝟏 𝐣𝐣−𝐥𝐥 𝑸𝑸𝒋𝒋−𝒍𝒍+𝟏𝟏 ⋯ 𝚫𝚫 𝐣𝐣−𝟏𝟏 𝑸𝑸𝒋𝒋 𝜴𝜴𝒋𝒋 −𝟏𝟏 −𝟏𝟏 𝑯𝑯𝒋𝒋𝒋𝒋+𝒍𝒍 = 𝚺𝚺𝒋𝒋+𝒍𝒍 𝑸𝑸𝒋𝒋+𝒍𝒍 ⋯ 𝜮𝜮−𝟏𝟏 𝒋𝒋+𝟏𝟏 𝑸𝑸𝒋𝒋+𝟏𝟏 𝜴𝜴𝒋𝒋
𝑗𝑗 𝑡𝑡𝑡
⋱
𝑙𝑙 𝑡𝑡𝑡
−𝑸𝑸𝑻𝑻𝒋𝒋+𝟏𝟏
𝑹𝑹𝒋𝒋 −𝑸𝑸𝒋𝒋+𝟏𝟏 𝑹𝑹𝒋𝒋+𝟏𝟏
⋱
𝑽𝑽 =
𝑇𝑇 𝑅𝑅𝑙𝑙 −𝑄𝑄𝑙𝑙+1 −𝑄𝑄𝑙𝑙+1 𝑅𝑅𝑙𝑙+1
Untraceable!
? 𝐻𝐻𝑘𝑘𝑘𝑘 𝐻𝐻𝑘𝑘𝑘𝑘 =
11
⋱
⋱
Identifiability of Damage Severity
Input Output
damage
Lemma 𝑯𝑯𝑗𝑗 |(𝑝𝑝, 𝛾𝛾1 ) ≠ 𝑯𝑯𝑗𝑗 |(𝑝𝑝, 𝛾𝛾2 ) if 𝛾𝛾1 ≠ 𝛾𝛾2 for 𝑝𝑝 = 1,2, … 𝑛𝑛 and 𝑗𝑗 = 1,2, , … 𝑛𝑛 + 1, where 𝑯𝑯𝑗𝑗 |(𝑝𝑝, 𝛾𝛾1 ) represents the matrix 𝑯𝑯𝑗𝑗 given parameters (𝑝𝑝, 𝛾𝛾1 ). Input 𝑛𝑛 2 Output
Corollary 𝑯𝑯𝑛𝑛 |(𝑝𝑝, 𝛾𝛾) = 𝑯𝑯𝑛𝑛 |(𝑛𝑛 + 1 − 𝑝𝑝, 𝛾𝛾) for an even number 𝒏𝒏 2
2
Air Force DDDAS Program PI Meeting, 01/27/2016
12
Identifiability of Damage Location Input 𝒑𝒑∗ , 𝜸𝜸∗ , 𝓛𝓛 and 𝒕𝒕
• Transfer Function can also be expressed as 𝐻𝐻 𝑠𝑠 =
𝑐𝑐1 s−𝑧𝑧1 )(s−𝑧𝑧2 )(s−𝑧𝑧3 )⋯(s−𝑧𝑧𝑚𝑚−1 )(s−𝑧𝑧𝑚𝑚 2 )(s−𝑓𝑓 2 𝑐𝑐2 s−𝑓𝑓12 )(s−𝑓𝑓22 )(s−𝑓𝑓32 )⋯(s−𝑓𝑓𝑛𝑛−1 𝑛𝑛
Update 𝑲𝑲
,
𝑓𝑓𝑖𝑖2 = 𝜆𝜆𝑖𝑖 , eigenvalue of 𝑴𝑴−𝟏𝟏 𝑲𝑲 𝒑𝒑, 𝜸𝜸
• Sufficient condition: the difference in natural frequencies lead to different transfer functions
Compute 𝒇𝒇𝒊𝒊
𝒑𝒑∗ ,𝜸𝜸∗
Compute 𝑆𝑆𝑘𝑘 = min ∑𝑖𝑖∈ℒ for all 𝑘𝑘
𝜸𝜸
Sort 𝑆𝑆𝑘𝑘 as 𝑺𝑺
𝟏𝟏
𝒑𝒑 ,𝜸𝜸 𝒑𝒑 ,𝜸𝜸 𝑓𝑓𝑖𝑖 ∗ ∗ −𝑓𝑓𝑖𝑖 𝒌𝒌 𝒑𝒑 ,𝜸𝜸 𝑓𝑓𝑖𝑖 ∗ ∗
≤ 𝑺𝑺
If 𝑺𝑺 𝟐𝟐 − 𝑺𝑺 𝟏𝟏 > 𝒕𝒕
𝟐𝟐
≤⋯
Yes
𝒑𝒑∗ , 𝜸𝜸∗ is identifiable
No
Air Force DDDAS Program PI Meeting, 01/27/2016
𝒑𝒑∗ , 𝜸𝜸∗ is NOT identifiable
13
Identifiability of Damage Location: Example
𝑆𝑆𝑘𝑘 = min � 𝜸𝜸
0.5
𝑖𝑖∈ℒ
𝑓𝑓𝑖𝑖
𝒑𝒑∗ ,𝜸𝜸∗
𝑓𝑓𝑖𝑖
− 𝑓𝑓𝑖𝑖
𝒑𝒑∗ ,𝜸𝜸∗
𝒑𝒑𝒌𝒌 ,𝜸𝜸
𝓛𝓛 = {𝟏𝟏: 𝟓𝟓}
0.8
𝑆𝑆𝑘𝑘 0
e.g. 𝑝𝑝∗ , 𝛾𝛾∗ = 23,0.4
𝑆𝑆𝑘𝑘 1
𝑘𝑘
• The minimum value of 𝑆𝑆𝑘𝑘 = 0 when 𝑘𝑘 = 23 and 38 for a beam with 60 elements in a fixedfixed experimental setup. • The parameter 𝑝𝑝∗ , 𝛾𝛾∗ = 23,0.4 can only be identified in pairs due to symmetry. Air Force DDDAS Program PI Meeting, 01/27/2016
60
0
1
𝑘𝑘
60
• The minimum value of 𝑆𝑆𝑘𝑘 = 0 when 𝑘𝑘 = 23 only for a beam with 60 elements in a fixed-free experimental setup. • The parameter 𝑝𝑝∗ , 𝛾𝛾∗ = 23,0.4 can be uniquely identified.
14
First-year Progress: Adaptive Sensor & Damage Identification
Structural damage identification using piezoelectric impedance active sensing
Air Force DDDAS Program PI Meeting, 01/27/2016
15
Piezoelectric Impedance Sensory System: Concept Impedance sensor and measurement: two-way electromechanical coupling Structural impedance coupled with transducer impedance Resistor
Transducer
200 Real Part of Impedance (Ohm)
Impedance with noticeable changes
160
120
80
40 5
10
15 20 Excitation frequency (kHz)
Transducer
25
30
35
Circuitry elements
New idea: impedance approach + piezoelectric circuitry
Integrate circuitry elements to enhance dynamics & resonant effects Fundamental advantage: circuitry can amplify signals and enrich measurements
Air Force DDDAS Program PI Meeting, 01/27/2016
Mathematical Formulation Finite element-based formulation Transducer
+ Cq + Kq + K 12Q = Mq 0 + RQ + K Q + K T q = LQ Vin c 12 • M,C,K – mass matrix, damping matrix, stiffness matrix • 𝐊𝐊 𝟏𝟏𝟏𝟏 – coupling vector
Circuitry elements
• q – displacement vector • Q – electrical charge
• L, R – circuitry inductance and resistance
Under harmonic excitation, the impedance measurement is Z=
Vin 1 T ( K − Mω 2 + Cωi ) −1 K 12 ] = [( −ω 2 L + iω R + K c ) − K 12 Q ωi
Structural damage impedance change: How to identify damage under impedance changes (at various frequencies)?
Damage Quantification and Influence to Piezoelectric Circuitry Impedance Assumption: damage is modeled as elemental stiffness reduction • 𝐋𝐋 – indicates how the elemental stiffness
= Kd
m
∑ L [K j =1
T j
ej
(1 − D j )]L j
𝐣𝐣
matrices are assembled • 𝐷𝐷𝑗𝑗 – damage index (percentage of stiffness reduction)
Impedance of damaged structure
m Vin 1 T Z d = = [ Z c − K 12 ( ∑ LTj [K ej (1 − D j )]L j − Mω 2 + Cωi ) −1 K 12 ] Q ωi j =1
• Taylor series expansion • Small damage - Disregard the higher order terms
∂Z Z d ( D1 , D2 ,...Dm ) ≈ Z ( D1 , D2 ,...Dm = 0) + ∑ d |D1 , D2 ,... Dm = 0 D j j =1 ∂D j m
Damage Sensitivity Matrix Impedance change expression: m
1 T K 12 ( K − Mω 2 + Cωi ) − ( LTj K ej L j )( K − Mω 2 + Cωi ) − K 12 }D j ωi
∆Z = Zd − Z (D = 0) = ∑{ j =1
Compare impedances at n frequency points
∆Z (ω1 ) D1 ∆Z (ωi ) = [ S nm ]n×m Di ∆Z (ω ) D n n×1 m m×1 Sensitivity matrix entries: 1 T 2 − T 2 − S= K K − M ω + C ω i L K L K − M ω + C ω i K12 ( ) ( )( ) nm 12 n n m em m n n ωn i Grand Challenge: Problem is severely under-determined, as the number of finite elements (m) in practical structure is very large, much greater than the number of impedance measurements (n)!
Damage Identification Approach – Exclusion Solution to challenge: pre-screening of possible regions Consider the most likely case of single damage ∆Z (ω1 ) D1 ∆Z (ωi ) = S n×m Di ∆Z (ω ) D n n×1 m m×1
∆Z (ω1 ) 0 ∆Z (ωi ) = S n×m Di = ∆Z (ω ) 0 n n×1 m×1
S1,i S D j ,i i S n ,i
Impedance change vector in theory similar to one column of the sensitivity matrix(𝐒𝐒:,𝑖𝑖 ) (multiplied with a scale factor(𝐷𝐷𝑖𝑖 )) The ratio of the two vectors is the damage severity
Eliminate impossible columns/elements Consider each column of S matrix D pi = mean(ΔZ. / S:,i ) All elements with D pi < 0 or D pi > 1 can be directly excluded This approach is applicable regardless the mesh size; Can be extended to ‘super-elements’ or ‘regions’.
Damage Identification Approach – Similarity Comparison Damage prediction: First estimate impedance change correlation ∆Zei = D pi . × S:,i Then compute and check similarity between estimated impedance change and the actual measurement SI i =∆ ( Z ei ⋅ ΔZ) /( ΔZ ⋅ ΔZ)
The closer 𝑆𝑆𝑆𝑆𝑖𝑖 to 1, the more likely damage occurs at the i-th element and 𝐷𝐷𝑝𝑝𝑖𝑖 is the severity. Define a similarity function for the i-the element Pi ( SI i ,1) =
1 / SI i − 1
sum(1 / SI i − 1)
Higher Pi indicates more probable damage location
New algorithm features: Avoid inversion of under-determined matrix equation Naturally fit the probabilistic and data-driven analysis framework
Case Illustration - Setup Cantilever plate
Piezoelectric transducer • Size: 15 mm×19.05mm×1.4 mm • Piezoelectric constant : 𝑑𝑑31 : −140 × 10−12 𝑚𝑚/𝑉𝑉 𝑔𝑔33 :25 × 10−3 𝑉𝑉𝑉𝑉/𝑁𝑁
• Size: 561 mm×19.05mm×4.7625mm • Young’s modulus: 6.89 GPa • Density: 2700 𝐾𝐾𝐾𝐾/𝑚𝑚3
4.7625mm
561 mm
Piezoelectric transducer 19.05mm
180 mm
15 mm Signal analyzer (Agilent 35670A)
Vin Vout
PZT Cantilever plate
Small resistor Rs
STRUCTURE AND SYSTEM DYNAMICS LAB
Case Illustration – Model Calibration Finite element model with convergence study and model updating To facilitate high detection sensitivity, high-frequency impedance is needed, requiring large number of elements. The plate is discretized into 12500 elements, grouped into 250 segments. FEA Model updating of healthy baseline is EXP (Hz) Original Updated (Hz) (Hz) conducted. Boundary conditions and 12.20 12.46 12.08 material properties are updated to 74.00 76.50 74.29 match FE predictions with experiments. 209.50 213.95 208.22 S151
S152
S153
S154
S225
S76
S77
S78
S79
S150
S1
S2
S3
S4
S75
414.70 682.80 1017.10 1893.18 2428.44 3033.20 3703.09
424.47 699.14 1040.69 1936.73 2482.46 3104.30 3781.89
413.34 681.33 1015.93 1893.58 2429.08 3034.05 3704.05
Impedance Measurements A small mass is attached onto the plate to emulate damage (Segment 103; 0.5% equivalent stiffness change) Impedances are measured around 14th, 16th , 20th resonance frequencies
(1) 14th resonance frequency Air Force DDDAS Program PI Meeting, 01/27/2016
(2) 16th resonance frequency
(3) 20th resonance frequency
Identification Result Damage Segment # Similarity severity Use measurements around the and 103 0.0786 0.005 20th resonances: 214 0.0622 0.003 64 0.0622 0.003 55 0.0522 0.002 205 0.0522 0.002 121 0.0521 0.005 161 0.0482 0.003 11 0.0482 0.003 139 0.0482 0.003 169 0.0464 0.002 19 0.0464 0.002 86 0.0429 0.003 28 0.0392 0.005 178 0.0392 0.005 94 0.0379 0.002 196 0.0364 0.005 Location 46 0.0364 0.005 Multiple possible damage locations are predicted; 130 0.0353 0.002 Actual damage location has the highest likelihood; 131 0.0090 0.003 56 0.0070 0.003 Damage severity prediction is accurate.
Similarity
16th
Air Force DDDAS Program PI Meeting, 01/27/2016
Identification Result
Similarity
Use measurements around one more resonance, the 14th:
Segment # Similarity 103 0.088 121 0.084 28 0.084 178 0.084 196 0.077 46 0.077 140 0.012 85 0.012 102 0.012 122 0.011
Damage severity 0.005 0.005 0.005 0.005 0.005 0.005 0.010 0.009 0.004 0.005
Location
Damage prediction result is improved: number of possible damage locations reduced drastically, decreasing the false-alarm possibility Perfectly suited for probabilistic-based inference analysis Air Force DDDAS Program PI Meeting, 01/27/2016
Probabilistic Inference – Another Study Another approach to avoiding direct matrix inversion Update the probability estimate of interested parameters for a hypothesis as additional information is acquired – take the effect of uncertainty into account Can be fully represented by Bayes’ rule p (θ | D) =
Incorporate data
p (D | θ) p (θ) ∫ p (D | θ) p (θ)dθ
Prior PDF
Posterior PDF
:prior distribution / hypothesis p (D | θ) :likelihood function / evidence p (θ | D) :posterior distribution p (θ)
θ: model parameters – Damage location & severity;
D: data/observations – piezoelectric impedance measurements; Likelihood function 𝑝𝑝(𝐷𝐷|𝜃𝜃): p ( D | θ ) =
1
2πε
e
Air Force DDDAS Program PI Meeting, 01/27/2016
− ( D − D ( θ ))2 2ε
𝐷𝐷:measured impedance � 𝐷𝐷:associated impedance parameterized by 𝜃𝜃�
Bayesian Inference – Challenge and Possible Solution Bayesian inference procedure Generate a full set of parameter samples subject to assumed distribution Call finite element analysis of each updated model (w.r.t each parameter sample) to derive the associated impedance prediction Use Bayesian inference to predict the parameter distribution
Challenge:
Size of parameter sample is huge — a lot of possible damage locations and severity levels Repeated FEA under large parameter space— huge Computational cost
Possible Solution:
Only perform finite element analysis under a reduced sample size Gaussian Process—enrich the FEA output under full set of parameters
Air Force DDDAS Program PI Meeting, 01/27/2016
Case Illustration – Problem Setup Bayesian inference with Gaussian Process
117 mm
PZT 134 mm
plate
PZT
circuitry
Clamped-clamped aluminum plate Size :610 mm×508 mm×4.826mm Young’s modulus: 73 GPa Density: 2780 𝐾𝐾𝐾𝐾/𝑚𝑚3 PZT transducer (PZT-4): Size :44.88 mm×44.88 mm×2.8mm Piezoelectric constant 𝑑𝑑31 : −140 × 10−12 𝑚𝑚/𝑉𝑉 𝑔𝑔33 :25 × 10−3 𝑉𝑉𝑉𝑉/𝑁𝑁
Finite model and parametrization: 2688 elements 9 severity levels
Computational costs: Full-size: 2688 X 9 rounds of finite element dynamic analysis; 145 hours Reduced size: 165 X 3 rounds; 3 hours Air Force DDDAS Program PI Meeting, 01/27/2016
Case Illustration – Result When impedance measurements around two resonances are used Location
Severity
Probability
1
39
6
3.1551
2
143
5
3.1550
3
24
3
3.1550
4
49
2
3.1547
5
129
5
3.1546
6
58
2
3.1546
7
72
5
3.1541
Posterior Distribution 4 3 2 1 0 6 4 2 0
5
0
15
10
20
25
30
35
Severity 5 4 3 2 1 0 5 4 3 2 1
0
10
20
30
Severity 10 Air Force DDDAS Program PI Meeting, 01/27/2016
40
The results indicate to a possible area The actual damage parameter has the 7th highest possibility
Case Illustration – Result When impedance measurements around a third resonance is added
Location Severity
Posterior Distribution 10 8 6 4 2 0 5 4 3 2 1
0
10
20
30
40
Severity 5
Probability
1
72
5
10.18
2
36
2
10.12
3
159
2
9.99
4
15
3
9.93
5
105
5
9.88
As more information is obtained and applied, we can get better result – the actual damage parameter now has the highest possibility The computational cost is still relatively high On-going effort: To combine pre-screening with Bayesian inference! Air Force DDDAS Program PI Meeting, 01/27/2016
Summary and Ongoing Work Brief Summary on Identifiability Issue Theoretically proved the identifiability of damage at a fixed location in the FEM-based beam structure. A numerical algorithm is proposed to numerically check and validate the location identifiability.
Brief Summary on impedance-based damage identification
A promising damage pre-screening scheme is developed.
Probabilistic inference framework has been established.
On-going work on the estimation of damage changing trend Existing studies focus on static detection, that is, one time estimation of the current health status with available data. Lack of the ability to inference the future trend of the health status. A Dynamic Data Driven Based Hierarchical Bayesian Degradation Model will be developed. As data are accumulated, the new method can recover the true trend quickly.
Air Force DDDAS Program PI Meeting, 01/27/2016
32
Publications •
•
• • •
Q. Shuai, K. Zhou, and J. Tang , 2015, Structural Damage Identification Using Piezoelectric Impedance and Bayesian Inference, Proceedings of SPIE, Smart Structures / NDE, V9435, 2015. Yuhang Liu , Shiyu Zhou, Jiong Tang, 2015, Identifiability Analysis of Finite Element Models for Vibration-Response Based Structure Damage Detection, ASME Transactions, Journal of Dynamic Systems, Measurement, and Control, in revision Hou, Y., Wu, J., and Chen, Y. 2015, “Online Steady State Detection Based on Rao-Blackwellized Sequential Monte Carlo”. IEEE Transactions on Automation Science and Engineering, in revision K. Zhou, Q. Shuai, and J. Tang, Damage identification using piezoelectric impedance and intelligent inference, to be submitted to Journal of Intelligent Material Systems and Structures. Liu, Y., Zhou, S., and Tang, J, Dynamic data-driven Bayesian structural failure prognosis, to be submitted to IIE Transactions.
All publications acknowledge the support by AFOSR DDDAS Program under grant FA9550-14-1-0384. Air Force DDDAS Program PI Meeting, 01/27/2016
33
