建築物微振量測之最佳感測器配置

(1)

國立臺灣大學工學院土木工程學系博士論文

Department of Civil Engineering College of Engineering National Taiwan University

Doctoral Dissertation

建築物微振量測之最佳感測器配置

Optimal Sensor Placements for Ambient-Vibration Monitoring of Building

楊晏瑜 Yang, Yen-Yu

指導教授：呂良正教授 Advisor: Prof. Leu, Liang-Jenq

中華民國 112 年 1 月

January 2023

(2)

(3)

(4)

(5)

致 ^謝

這本論文能順利完成，首先要非常感謝恩師呂良正老師，他多年來悉心指導以及照顧，讓我受益匪淺，由於呂老師諄諄教誨以及不斷的督促指導，使得我能在原本不熟悉的領域上順利完成此本論文，有幸能跟隨呂老師多年，在休學當助教期間，呂老師也不曾放棄過我，從老師身上學到非常多知識以及做研究的態度，對我未來影響非常深遠，在此謹致上由衷感謝。

論文口試承蒙羅俊雄教授、張國鎮教授、宋裕祺教授、韓仁毓教授以及黃仲偉教授指導，並給予適當建議使其論文更臻完備，在此感謝。

感謝台灣大學醫學院附設癌症中心、國家地震研究中心提供場地及數據，讓我在實驗上能順利完成。

在台大多年認識非常多的師長們，在此特別感謝歐昱辰教授在研究上給予指導與協助；也感謝大師兄黃仲偉學長、唐瑢書學姊在生活上的協助以及解惑問題；

組辦的宛香姊謝謝您常常聽我抱怨跟我聊天以及幫我許多事務，讓我在博士生涯中不會感到無聊；感謝可葳韶珍夫婦為我點綴博士班生活的色彩；也感謝助教室的小江、憶玲、怡穎、屏先大哥，曉梅以及離職的高宇、紫芳、榆峻，在我當助教時，能協助我；感謝閔中從碩士時期與我一起奮鬥，感謝鈺翔夫婦對我的協助，

時常聽我在那邊抱怨；感謝嘉玟、靚芳、祐年、昌翰為研究室帶來歡樂；感謝一起奮鬥的至伸、翔佑、俊穎、郁珊、子欽、奕汎、少驊、盧易、享樑、銘姍、孟笙、琮堯、玟慧、佳逢、德鈞、彬祐、家萱、承憲、紹庭、穎君、俊廷、懷寬、

牧軒以及鈺庭、采霏、智傑、權恩等學弟妹們。

在此特別感謝女友辛玟錡陪我一起度過博士班最後三年，有妳在我才能走完博士學位的最後一哩路，好友煥勝、峻輔、志偉以及承偉能陪我聚餐哈拉，讓我不會感到無聊。

最後謹以本文獻給我最愛的父母親，感謝你們這些年來不懈地支持我、照顧以及關愛我，讓我在生活上無後顧之憂，我才得以完成博士學位；謝謝弟弟晏誠能在台中照顧爸媽，讓我不用擔心家裡情況，在台北專心地完成學業。

由衷感謝這一路上曾幫助過我的每個人，謹將本論文獻給所有關心我的人，

願所有人心想事成。

楊晏瑜謹誌於國立台灣大學土木研究大樓 511 室 2022 年 6 月

(6)

(7)

摘要

台灣位屬多地震區域，由於地震頻繁，為了結構物之安全性，需於定期進行結構物健康檢測，尤其是災害發生後能確保結構物之安全。設置感測器極為耗成本且耗時，故本文所提出之最佳化配置為減少再次檢測同棟結構物時所需之感測器顆數，並且提高判斷結構物模態頻率的精度。為此目的，本文結合了四種方法進行最佳化配置，分別為隨機子空間識別法、三次樣條內插法、K-means 演算法與基因演算法等方法判斷出加速度感測器配置之最佳樓層位置。提出各種不同模型，並與模態頻率解析解相互比較，驗證此法能以少數層樓之時間歷時，得到整體結構物的模態頻率，以確保結構物之安全性。為了驗證此方法之可行性，於數值模擬上，本研究建立幾個不同勁度折減之模型，再次驗證結構物災害發生前以及災害發生後，最佳配置差異；在現地實驗方面，本研究以六棟真實結構物進行實驗，分別為國立台灣大學土木研究大樓、國立台灣大學醫學院附設癌醫中心、

淡水施工中20 層樓建築物、板橋兩棟 20 與 21 層樓建築物及國家地震工程研究中心等建築，並以傅立葉快速轉換法來驗證此本文方法之可行性，由各種數值模擬以及實驗上證明。故本文最主要的貢獻在於提出之感測器最佳配置法，可以大幅減少第二次檢測時所使用之感測器數量與獲得較為精準之中高模態頻率，作為日後檢測之重要相關依據，當建築物受到損害時，即時檢測以瞭解建築物之受損情形。

關鍵字：系統識別、結構健康監測、隨機子空間識別法、頻率域分解法、三次樣條內插法、K-means 分群法，感測器最佳配置

(8)

(9)

ABSTRACT

Taiwan is located in the seismic zone with high frequency earthquake occurrences.

In order to increase structure safety, it needs to monitor the structural health before and after disaster occurs. This study proposes a method to obtain the optimal sensors placements (OSP), which could reduce the number of sensors for building monitoring.

In additions, the method could find out the higher modal frequencies for structures.

First, collect the real time-histories and use Cubic spline interpolation method to obtain simulated time-histories for each floor. Second, use Stochastic Subspace Identification to generate stabilization diagrams. Third, K-means clustering method is used to obtain modal frequencies. Finally, use Genetic Algorithm method to find OSP. In this study, created some model to proved this method that can obtained the high accuracy modal frequency, and then created some model which stiffness decrease to proved that used OSP method to obtain modal frequency before and after stiffness decrease which similar with real modal frequency. There are six in-situ experiments for the method verifying. The six in-situ building are NTU Civil Engineering Research Building, NTU Cancer Center, Tamsui building (under construction), two apartment complex in Banqiao and National Center for Research on Earthquake Engineering, respectively.

The OSP method can obtained great result from these in-situ experiments.

keyword：System identification, Structural health monitoring, Stochastic subspace identification, Frequency domain decomposition, Cubic-spline interpolation, K-means clustering, Optimal sensor placements

(10)

(11)

圖目錄

圖2-1. 所有樓層時間歷時分析之穩態圖 ... 18

圖2-2. 單一層樓分析之穩態圖 ... 19

圖2-3. Newmark-β 法之時間歷時與其內插之時間歷時比較 ... 23

圖2-4. Newmark-β 法之時間歷時之 FFT 頻率域圖 ... 23

圖2-5. 內插法之時間歷時之 FFT 頻率域圖 ... 23

圖2-6. 樓層數與參考感測器關係圖 ... 25

圖2-7. 以台大土木研究大樓內插比較圖 ... 26

圖2-8. 以台大土木研究大樓外插比較圖 ... 27

圖2-9. 真實量測速度時間歷時圖與 FFT 圖 ... 28

圖2-10. 內插速度時間歷時圖與 FFT 圖 ... 28

圖3-1(a). 穩態圖中所有資料點 ... 32

圖3-1(b). 隨機灑取 n 個點位作為群心 ... 32

圖3-1(c). 計算各點到群心之絕對距離 ... 33

圖3-1(d). 將各點歸群到最近的群心上 ... 33

圖3-1(e). 重新計算群心 ... 33

圖3-1(f). 再次分群直到收斂為止 ... 34

圖3-2. 基因演算法流程圖 ... 38

圖3-3. 最佳感測器配置流程圖 ... 41

(15)

圖4-1. 真實感測器與虛擬感測器示意圖 ... 44

圖4-2. 10 層模型示意圖 ... 45

圖4-3. 10 層樓之自由振動穩態圖 ... 47

圖4-4. 5 層樓最佳配置之穩態圖 ... 49

圖4-6. 方形窗 ... 52

圖4-7. 梯形窗 ... 53

圖4-8. 原始資料與重新取樣後資料之比較圖 ... 54

圖4-9. 10 層樓修正後之穩態圖 ... 55

圖4-10. 20 層未修正最佳配置之穩態圖 ... 57

圖4-11. 20 層穩態圖與 FDD 圖 ... 63

圖4-12. 20 層樓修正後最佳配置之穩態圖 ... 63

圖4-13. 15 層樓最佳配置穩態圖 ... 64

圖4-14. 30 層樓最佳配置穩態圖 ... 66

圖4-15. 10 層樓取三層之穩態圖 ... 68

圖4-16. 5 層樓加權配置法之穩態圖 ... 71

(16)

圖4-21 模態振型配置法流程圖 ... 79

圖5-1. CV-374V 速度感測器 ... 93

圖5-2. CV-374V 詳細規格表 ... 93

圖5-3. 土木研究大樓 ... 94

圖5-4. 感測配置位置以及平面圖 ... 95

圖5-5. 土木研究大樓最佳配置穩態圖 ... 97

圖5-6. 土木研究大樓單獨 3 樓分析穩態圖 ... 98

圖5-7. 土木研究大樓配置 3 樓與 7 樓分析穩態圖 ... 98

圖5-8. 台灣大學醫學院附設癌醫中心 ... 103

圖5-9. 台灣大學醫學院附設癌醫中心感測器設置位置及方向 ... 104

圖5- 10 台灣大學醫學院附設癌醫中心感測器設置樓層 ... 104

圖5-11. 台灣大學醫學院附設癌醫中心最佳配置穩態圖 ... 105

圖5-12. 淡水 21 層樓感測器配置樓層 ... 108

圖5-13. 淡水 21 層樓感測器設置平面示意圖 ... 109

圖5-14. 淡水 21 層樓穩態圖 ... 111

圖5-15. 板橋區 A 棟與 B 棟感測器設置平面圖 ... 112

圖5-16. 板橋區兩棟社區大樓感測器配置樓層 ... 112

圖5-17. A 棟穩態圖 ... 114

圖5-18. B 棟穩態圖 ... 114

(17)

圖5-19. 國家地震研究中心 ... 118

圖5-20. 國家地震研究中心感測器配置位置以及檢測方向 ... 119

圖5-21. 國家地震研究中心最佳配置穩態圖 ... 119

圖5-22. 國家地震研究中心 5 樓 FFT 頻率域圖 ... 121

圖5-23. 國家地震研究中心 13 樓 FFT 頻率域圖 ... 122

(18)

表目錄

表2-1. 各層數值模型參數 ... 19

表2-2. 5 層數值模型之頻率解析解(Hz) ... 19

表2-3. 1 層至 5 層分析(Hz) ... 19

表2-4. 單獨 5 層分析(Hz) ... 20

表2-5. 5 層數值模型之 5 個模態頻率平均誤差比較(Hz) ... 20

表2-6. 參考感測器與樓層數模態頻率平均誤差關係表 ... 24

表4-1. 5 個模型參數 ... 46

表4-2. 10 層樓模型理論模態頻率 ... 46

表4-3. 10 層樓自由振動最佳配置頻率與理論頻率誤差表 ... 48

表4-4. 5 層樓參考模態頻率與理論模態頻率誤差表 ... 50

表4-5. 5 層樓最佳配置模態頻率與理論模態頻率誤差表 ... 50

表4-6. 10 層樓未修正之最佳配置模態頻率誤差表 ... 51

表4-7. 10 層樓修正後之最佳配置模態頻率誤差表 ... 56

表4-8. 20 層樓未修正最佳配置模態頻率誤差表 ... 58

表4-9. 20 層樓已修正最佳配置模態頻率誤差表 ... 62

表4-10. 15 層樓最佳配置模態頻率誤差表 ... 65

表4-11. 30 層樓最佳配置模態頻率誤差表 ... 67

表4-12. 不使用 CSI 法之 1、5、10 層分析模態頻率 ... 69

(19)

表4-13. 有效模態質量表 ... 71

表4-14. 5 層樓加權配置模態頻率誤差表 ... 74

表4-19 理論模態振型 ... 80

表4-20 10 層樓參考模態振型 ... 81

表4-21 10 層樓參考模態振型與理論模態振型 MAC 值 ... 81

表4-22 10 層樓最佳配置位置之模態振型 ... 83

表4-23 10 層樓參考模態振型與理論模態振型 MAC 值 ... 83

表4-24. 10 層樓模型總比較表 ... 85

表4-25. 各模型勁度折減後數值模型參數 ... 86

表4-26. 10 層樓 case1 勁度折減最佳配置與原先配置比較表 ... 88

表5-1. 土木研究大樓參考模態頻率 ... 96

表5-2. 土木研究大樓最佳配置樓層與其他配置樓層模態頻率表 ... 96

(20)

表5-3. 土木研究大樓最佳配置與文獻及 FFT 頻率比較表 ... 100

表5-4. 不同時間檢測之土木研究大樓最佳配置模態頻率 ... 102

表5-5. 國立臺大醫學院附設癌醫中心最佳配置與 FFT 比較表 ... 106

表5- 6 淡水 20 層樓建築物之參考模態頻率與最佳配置頻率表 ... 110

表5-7. 板橋區 A 棟最佳配置頻率與 FFT 頻率比較表 ... 115

表5-8. 板橋區 B 棟最佳配置頻率與 FFT 頻率比較表 ... 116

表5-9. 國家地震中心最佳配置模態頻率表 ... 120

(21)

第一章緒論

1.1 ^研究背景

台灣位於環太平洋地震帶上，地處歐亞板塊與菲律賓板塊的交會處，地震發生的相當頻繁，基於此先天因素，其建築物耐震設計之考量就顯的相當重要。微振量測實驗(Ambient vibration test)與系統識別(System identification) 的技術進步，

使用微振量測求得結構物的結構特性越來越普遍。當地震發生時，會影響地震對建築物是否造成損傷的主要因素為建築物之基本振動周期與阻尼比，對於建築物基本振動周期的探討已有不少研究，如葉祥海等人(2000、2001)。陳易瑞(2002)利用 ETABS 數值結構模型找出有無非結構牆之基本振動周期比值，進而修正微振量測所得結果，定義設計周期上限值與周期建議公式，並提出結構物空構架之基本振動周期。劉醇宇(2003) 以建築物結構模型建立建築物勁度與強度的關係式，

並以微振量測與系統識別技術所得到的建築物基本周期，用於建築物耐震能力快速評估的依據。

模態參數可區分為自然頻率(Natural frequency)、阻尼比(Damping ratio)與模態振型(Mode shape)，其中自然頻率為本文所探討之分析重點。建築物的自然頻率需先檢測結構物之動力反應，目前檢測方式大致可分為直接量測法與間接量測法，前者是將感應器直接安裝在建築物上，透過感應器直接量測並紀錄建築物的反應，此過程之反應稱為時間歷時，直接量測法是目前最多人會採用的方式，但直接量測法會因根據所需的模態參數決定建築物量測位置而需要大量佈點，並在實驗施作時要花費大量時間與金錢，也會給當地住戶帶來不便等缺點；間接量測則為將一個被測量轉化為若干可直接測量的量加以測量，而後再依據由定義或規律導出的關係式(即測量式)進行計算或作圖，從而間接獲得測量結果的測量方法。

(22)

為了確保建築物經災害發生後是否安全，需對災害發生前與災害發生後之結構物進行結構健康監測或檢測，從而判斷出結構物是否損害及勁度是否下降。由於採用直接量測法進行監測或檢測時，架設實驗儀器及量測非常耗時的，若需多層樓同步量測更是個艱困的工作，為求簡化設置實驗儀器，本文針對實驗儀器設置樓層進行最佳化配置，利用最少感測器達到量測結構物之模態頻率的最大成效。

結構物配置感測器時，為了得到各層樓的時間歷時，Limongelli (2003)曾提出利用三次樣條內插法(Cubic spline interpolation Method 簡稱 CSI)得到未設置感測器的樓層之時間歷時，並最佳化感測器之位置，得到該棟結構物之模態形狀。

Kodera et al.(2018)使用少數樓層真實量測資料並搭配 CSI 製作了 1/3 比例的鋼結 構模型，放置於振動台上以不同的地震進行振動測試與模擬，，將取得之新的量測資料再與內插模擬之樓層資料進行比較，以驗證此內插法之準確性。Limongelli (2004)也利用 CSI 於橋梁檢測，對未設置感測器之量測點進行內插模擬，也有不錯之成效。

過去已有學者提出最佳化感測器配置方法，如Kammer (1991)提出一種有效之演算法，基於各感測器配置位置對識別之模態皆有線性獨立貢獻，可由初始擺設位置快速篩選出最佳配置位置。Guo, et al. (2004)先找出結構損傷參數，利用此參數配合簡易基因演算法去求得感測器最佳配置位置，因簡易的基因演算法較不精準，故提出改進基因演算法得到最佳感測器配置位置。Meo & Zumpano (2005) 運用有效獨立驅動點殘差法(Effective independence driving-point residue method，

簡稱 EFI-DPR)進行最佳感測器配置位置，此法可強化訊號強度和抗噪能力，藉此提高獲取橋梁模態頻率之精準度。

本研究採用之系統識別為隨機子空間識別(Subspace stochastic identification，

簡稱SSI)，於 1991 年 SSI 由 Van Overschee (1996)在第 30 屆 IEEE(The Institute of Electrical and Electronics Engineers)大會上首次提出，對於目前各種時間域輸出系統識別法，以 SI 為基礎推導理論背景、得系統參數識別結果呈現合理且較佳之結論，為近年來振動分析領域熱烈討論之識別方法。根據量測資料組成空間矩

(23)

陣並從中分離出觀測矩陣、控制矩陣及噪訊矩陣等子空間矩陣，以取得系統參數矩陣，進而求得對應之模態頻率、阻尼比、與量測位置之振態，代表結構特定狀態之動態參數。

基於隨機子空間之時間域唯輸出系統識別方法可概分為協方差型 (SSI- COV)與資料型(SSI-DATA)兩種方法。Peeters (1993)首次提出 SSI-DATA 法並用此法於土木工程模態分析，而 SSI-COV 法首次於 2000 年 Peeters 的博士論文中提及，這兩種方法如同其他唯輸出系統識別法具有基本假設，像是輸入外力為假定成零平均值之白噪音，不同時間點之相關性為零，同時假設之輸入系統狀態與輸出反應並無相關性，但隨機子空間法之狀態空間僅適用於線性非時變系統，

故有許多學者研究成果發表提出了改良方法；Fan, et al.(2004)利用移動窗格對 SSI-COV 法處理過程進行改善並應用於香港九龍橋颱風氣候下之監測，以減少識別所需時間。Chang et al.(2007)提高穩態圖中模態參數識別精準度。

隨機子空間識別會因輸入樓層檢測資料的不同所得到的分析也不盡相同，且提供的樓層時間歷時資料越齊全，)，所得到的模態頻率、阻尼及振態形狀也越精確。因此若要量測一棟大樓，希望得到該棟大樓各模態頻率，則需各層樓都架設感測器進行同時量測，經SSI 識別所得到之結構物模態頻率相較於一般擺在樓頂上量測識別出來之模態頻率來的精確，但安裝感測器成本較高，若各層樓皆裝設感測器，是非常耗時且高成本，故本研究希望找出最佳感測器配置位置，能同時兼顧經濟性與準確性。

1.2 ^{研究動機與目的}

感測器配置是極為麻煩且複雜的事情，像是配置感測器點位、擺設感測器數量多寡、了解感測器之精度及高精度的感測器花費成本高等都是需探討的。本研究主要目的為用於再次對同棟建築物進行量測時採用最佳化配置，要得知一棟建築物是否損壞而需進行補強皆經由感測器進行量測並判定建築物是否受損，此過

(24)

程須進行多次量測，無法就單次量測便能判斷建築物有無受損。為減少後續量測工作，而發展出最佳化感測器配置法。

由於高模態頻率之建築物不易觀測，且廣泛使用之快速傅立葉轉換(Fast fourier transform)中高頻率易被雜訊影響，因此本研究將隨機子空間法為主體，發展最佳化感測配置，以最少之感測器得到較為精準的結構物模態頻率，得以判斷結構物是否受損需進行補強等後續工作。

1.3 ^文章架構

本文使用隨機子空間識別法、三次樣條內插法、K-means 演算法、基因演算法等四大方法，構建出了整個最佳化感測器配置法，下述分別介紹本研究各章節內容。

第一章：本文的研究背景與動機，各方法的發展內容與過程。

第二章：介紹微振與隨機子空間識別法、穩態圖、三次樣條內插法等方法，用簡

單的例子證明結構系統不足與找到之最後結構參數的差異，以及三次樣條內插法運用於結構系統上之可行性。

第三章：介紹非監督式機器學習中的K-means 演算法，使穩態圖中的極點分群得到各模態頻率值。並介紹基因演算法跟其理論和適應函數，最後再講述整體最佳化配置流程。

第四章：用五個剪力構架進行模擬，並搭配兩種外力：衝擊外力及白噪音。衝擊

外力對單一10 層剪力構架進行分析，而白噪音則為五個模型進行分析，並探討一些遇到之問題且如何解決，最後再提出若不用 CSI 法情形與若對低模態頻率進行加權後感測器配置的情況。

第五章：提出了五個實驗分別為台灣大學土木研究大樓、台灣大學醫學院附設癌

醫中心、淡水1 棟 20 層施工中建築物、板橋 2 棟 21 層與 22 層社區大樓進行現地實驗以及取得改建後國家地震研究中心時間歷時等6 棟建築物進行實驗分析。

(25)

第六章：整體最佳化感測器配置法做總結並提出改善以及未來之發展性。

(26)

第二章系統識別方法

本章節探討連續與離散狀態空間之模型，而後延伸成隨機狀態空間之模型，

並說明隨機子空間識別法(Stochastic subspace identification，簡稱 SSI)之假設及定義，從而探討SSI 之理論與使用流程。

2.1 ^微振量測

結構物的動力特性如自然振動周期、阻尼比與模態振型等，經由檢測結構物的振動行為，可由頻譜分析或系統識別等方法加以求得，而檢測資料包含強迫振動、地震、微振量測等各種方式，其中以微振檢測最為方便。微振檢測是利用結構物周遭的干擾作為外力，其外力源包含常態風力、車輛經過時產生振動與人員於建築物中所引起的微振動等，感測器紀錄建築物在微振動下之反應訊號，並分析求得建築物之動力參數。

周遭外力對結構物來說非常微小，可視之為微振( Ambient vibration )，微振之振幅約為 10^-5~10^-3g。近年來因為科技進步，高靈敏度的振動計(如速度規、加速度規等)快速發展並廣泛被使用，可精確量測極微弱之振動訊號。以國內對微振之研究發展來說，建築物之微振量測實驗早期如劉文(1976)、曾紫銘(1977)、馬俊強(1979)等人對數棟建築物進行微振量測，並分析其結構動力參數。發展至今感應器本身與可攜帶型主機的尺寸也不斷縮小，甚至不需帶主機與尋找電源即可進行檢測，因此使微振實驗之技術被廣泛地使用而逐漸趨於普遍，再加上檢測技術與分析方法純熟，使研究人員利用分析得到之結構動力參數，進行相關研究與統計分析。

微振檢測利用系統識別等方法求取結構物之動力特性，由以上假設可知，利用微振試驗以識別結構物的動力特性應具有以下優點(楊永斌等人，1995):

1. 時間歷時擷取長度越長可增進識別結果的精確度。

(27)

2. 結構物進行檢測時，微振情況與地震狀況或強迫振動三種情況下，微振檢測較具線性化，即在單一量測記錄中結構物之振動為線性彈性行為。

另一方面，假設外力為平穩過程(Stationary process)，因此根據輸出反應所求得之結構動力參數並不會因時間不同而改變。

3. 分析過程中，假設輸入為白躁音(White noise)訊號，即涵蓋各頻率範圍以及各頻率點之能量相同的訊號，故只需對輸出訊號作分析而毋須探討輸入訊號，因此可節省相當多的計算時間。

2.2 ^{快速傅立業轉換} (Fast Fourier Transform)

利用微振量測資料識別結構系統動態特性，大概可分為頻率域與時間域分析兩大類，而頻率域分析發展較早也較簡便。一般來說以微振檢測求取結構系統動態特性最簡單的方法為頻譜分析(Spectral analysis)，本節介紹快速傅立業轉換分析。根據頻譜分析中之單一輸入/單一輸出模式(Single input / Single output)，其自相關能量譜(Auto-correlation spectrum)的關係式為(Bendat & Piersol，2011)

(2.1) 其中Svv為量測物理量的自相關能量譜，Spp為所有外力與擾動的自相關能量譜，H 為系統轉換函數。假設在量測過程中所有的外力與擾動為白噪音，即Spp在每一個頻率之下均為常數，Svv與H²為倍數關係。從上述可知，得到觀測物理量之自相關能量譜，以獲得系統的轉換函數。並可藉由轉換函數的峰值，找出水平向以及扭轉向振動模態所相對應的頻率。然而由自相關能量譜與傅立業轉換(Fourier transformation)的關係為(Bendat & Piersol，2011)

(2.2)

其中E[]代表期望值運算子，V代表v(t)的傅立業轉換，V*代表V的共軛複數，T 為時間長度。自相關能量譜可以由傅立業轉換的振幅近似之。因此資料分析的流程如下:首先獲得時間歷時資料 (位移、速度或加速度)，然後對所記錄之微振歷

2

vv pp

S = H S

* 2

1 1

lim lim

vv T T

S E VV E V

T T

®¥ ®¥

é ù é ù

= êë úû= êë úû

(28)

時資料作FFT，求得傅立業振幅。假定系統無雜訊的干擾，則FFT之頻率域圖峰值所對應之頻率即代表結構物之自然振動頻率。

離散化的傅立業轉換即 FFT 可寫成下列形式：

(2.3) 其中f(k)為離散化的時間域函數，F(n)為離散化的傅立業轉換。N值最常選取值為

N = 2γ，其中γ為正整數。本文取樣頻率之200z，Δt=0.005 sec，量測時間為5分鐘，

共60000筆資料點。在N值的選取方面，當 N 取1024時，也就是說每 1024筆資料 即為一個周期，則60000筆資料可以視為58個周期的組合；當N 取2048時，則可以視為29個周期的組合，以此類推。但依據頻譜分析的誤差理論可知，當使用上述方法來求的自相關能量譜時，其正規化平方誤差(Normalized mean square error) 可表示為

(2.4)

其中Be為自相關能量譜在頻率域上的解析度(resolution bandwidth)，與N值大小有 關，即Be=1/N∆t。∆t 為採樣間隔時間，本研究均固定為0.005秒; 而Tr為總量測

時間，當N值過小時，則Be會過大，從式(2.4)中可以看出，此時第一項的隨機誤

差(random error)雖然會減小，但第二項的偏見誤差項(bias error)會隨之擴大，此結果反應在傅立業頻率域圖上的曲線有較為平滑之現象，但是頻率域之解析度太

差，即刻度過於疏鬆；反之，N值過大則會使第一項的隨機誤差過大，反應在傅

立業頻率域圖上的曲線則較為鋸齒狀，使得尖峰的判讀較顯不易，因此N值的選取必須折衷，並無統一之標準。

2.3 ^狀態空間

建立結構系統之參數識別，需先建立其結構模式，並根據空間領域分析、結構之物理行為、外力形式與訊號來源等特徵運用數學模型描述之。本節以狀態空間模型為基礎，先介紹運動方程式之連續時間狀態空間方程式，由此方程式推出

( )

¹

( )

⁽² ⁾^/

0

0,1, , 1

N kn N

k

F n ^- f k e^- ^p n N

=

å

= ^! -

2 " 2

2 1

24

e pp

e r pp

B S

B T S

e = + ç^æçè ^ö÷÷ø

(29)

離散時間下之狀態空間方程式，並根據馬可夫參數說明系統之確定性。最後整合外力項，修正成隨機離散時間之狀態空間方程式，此狀態空間方程式即為隨機子空間識別法之基礎方程式。

2.3.1 連續時間狀態空間方程式

考慮一個多自由度結構系統，其質心座標在地震運動作用下之運動方程式可表示如下：

(2.5) 或

(2.6)

上兩式中M、C 與 K 分別為結構系統之質量、阻尼與徑度矩陣，L 為影響

係數矩陣，U 為樓層相應之結構固定端相對位移，!̈^!為樓層之絕對位移，!̈_"為外力引起之絕對加速度。

令， (2.7)

將(2.6)式表示成狀態空間方程式，其方程式如下：

(2.8)

上式中x(t)為 n 階系統之狀態向量，n=2n0為系統階數，n0為模態數，u(t) 為系統之輸入向加速度向量。一般量測之時間歷時訊號為絕對加速度，(2.6) 可推得相關輸出狀態方程式。

(2.9) 其中 y(t)為裝設感測器絕對加速度量測之向量，Ca 為輸出位置影響矩陣 (Output Location Influence Matrix)，#̈^!為相應質心絕對加速度。

( ) ( ) ( )

⁰

U tt + U t + U t = M!! C ! K

( ) ( ) ( )

g

( )

U t + U t + U t = - U t M!! C ! K ML!!

( ) ( )

( )

t U t U t

é ù

= ê ú

ë û

x ! u t

( )

= !!U tg

( )

^t ⁼^é_ê_- ⁰_-₁ _- ^I_-₁ ^ù_ú

( )

^t ⁺^{é ù}_{ê ú}_-⁰ ^{u t}

( )

ë û ë û

x x

M K M C L

!

( )

t = aU^t

( )

t y C !!

(30)

(2.10)

由(2.8)與(2.10)兩式重新整理，可得連續時間之狀態方程式：

(2.11) (2.12)

其中、、

Ac為狀態矩陣，Bc為輸入影響矩陣(Input influence matrix)，Cc為輸出影響矩陣(Output influence matrix)，Dc為直接傳導矩陣(Direct transmission matrix)。

2.3.2 離散時間狀態空間方程式

上章節提到連續狀態空間方程式(2.11)，將其離散化

(2.13)

(2.14)

其中Δ 為取樣時間，令 t=kΔ，將(2.14)改寫為：

(2.15)

(2.16)

由於連續時間系統成立，定義其狀態方程式之解：

( ) ( )

1 1

( )

a

t U t

U t

- - é ù

é ù

= ë- - ûê ú

ë û

y C M K M C !

( )

t = c

( )

t + cu

( )

t x! A x B

( )

t = c

( )

t + cu

( )

t y C x D

1 -1

0 I 0

c - c

é ù é ù

=êë- - úû =êë- úû

A B

M K M C、 L C_c = -éë M K^-¹ -M C^-¹ ùû

[ ]

⁰

c = D

( ) ( ) ( )

lim0 x t x t x t D®

+ D -

= D

!

( ) ( ) ( )

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )

c c

c

c c

x t cx t t

x t x t x t u t

= + u

Þ + D - = D + D Þ + D = + D + D

A B

!

( ) ( ) ( ) ( )

( ) ( ) ( )

+ = I+ _c _c

c c

x k x k u k

y k k k

D D D D + D D

ìïí D = D + D

ïî

A B

C D

( )

( ) ( ) ( ) ( ) ( ) ( ) ( )

1 = I+ _c _c

c c

x k x k u k

y k k k

ì + D D D + D D

Þ íï

D = D + D

ïî

A B

C D

(31)

(2.17)

令t0=kΔt、t=(k+1) Δt，k∈N，上述式子可以改寫成

(2.18) '()) = {'(-∆/)|-∆/ ≤ ) < - + 1)∆/} (2.19) 將上述兩式子整合一下改寫成

(2.20) 其中

(2.21)

(2.22)

x(k)為第 k 點之離散時間狀態向量，y(k)為第 k 點的時間輸出向量，Ad為離

散狀態矩陣，Bd為離散輸入影響矩陣，Cd為離散輸出影響矩陣，Dd為直接傳導矩陣。

2.3.3 ^{隨機狀態空間方程式}

在量測上，環境條件具有不確定性，無法準確的掌握輸入與輸出反應，因此考慮系統之不確定因素，將式(2.21)與式(2.22)改寫成：

(2.23)

(2.24)

其中w(k)為過程噪訊向量，主要適用來表示系統與擾動所產生之誤差；v(k)為量

( ) ( ) ( )

( ⁰)

( )

^{( )}

( )

0 0

:

t t t t

t

k k

e t e u d

t

t t t

- -

= +

= D =

ò

c c

A A

x Bc

x x x

( )

⁽ ⁰⁾

( )

⁽⁽ ⁾ ⁾

( )

0

1 ^{t t} Δ ^t ^k 1 Δ^t

k t e ^- k t te ⁺ ^-^t u t td

é + D =ù +

ë û ^A^c

ò

^A^c ^c

x x B

( )

⁽ ⁰⁾

( ) ( )

0

1 Δ ^{t t} Δ ^t Δ

k t e ^- k t ^é te d^t^¢ t ^ù u k t

é + ù= + ¢

ë û ^A^c êë

ò

^A^c úû ^c

x x B

(

^k ^{1 Δ}

)

^t

t¢= + -t

Δ 1 2 1 3

I Δ ( Δ ) ( Δ )

2! 3!

t

d =e^A^c = + _c t+ _c t + _c t +…

A A A A

[ ]

0 c c

t c

t

d =

ò

e^A^c^t^¢dt¢ = - ^-1

B B A I A B

(

k ¹

)

d

( )

k du

( )

k Þx + =A x +B

( )

k = d

( )

k + du

( )

k

y C x D

(

^k^{+ =}¹

)

d

( )

^k ⁺ d^{u k}

( ) ( )

⁺ ^k

x A x B w

( )

^k ⁼ d

( )

^k ⁺ d^{u k}

( ) ( )

⁺ ^k

y C x D v

(32)

測噪訊向量，其與量測儀器之不準確性相關。此兩項無法從量測結果得知。假設噪訊為零平均值之白噪，其協方差矩陣可以表示成：

(2.25)

E 為期望值因子，Q、S 與 R 為協方差矩陣； 為Kronecker delta。

(2.26)

即w(k)與 v(k)之相關性為 0。

因真實外力不易求得，且外擾噪訊無法移除，因此將輸入項設為一隨機過程，

合併輸入向量及噪訊向量將式(2.23)與式(2.24)改寫為：

(2.27)

(2.28) 此兩式稱為隨機狀態空間方程式，外力向已含在噪訊項中，且w(k)與 v(k)仍 假設為零平均值之白噪音。

2.4 ^隨機過程

假設隨機過程為穩態(static state)：

(2.29) 為狀態協方差矩陣，在任意時間點下，狀態協方差矩陣均相同。且w(k)與 v(k)皆為零平均值，並假設與狀態向量 x(k)無相關性，因此：

(2.30) 用Lyapunov 方程式表示協方差矩陣：

( ) ( )

^p

( ( )

^T

( )

^T

)

^pq

E q q

p d

éæ ö ù æ ö

êç ÷ ú=ç ÷

è ø

êè ø ú

ë û ^T

w Q S

w v

v S R

dpq

1 0

pq

p q p q

d d

= Þ =

ìí ¹ Þ = î

(

^k^{+ =}¹

)

d

( ) ( )

^k ⁺ ^k

x A x w

( )

^k ⁼ d

( ) ( )

^k ⁺ ^k

y C x v

( ) ( ^{( )}

[ )^T , 0

E x k x k ùû= S Eéëx k ùû= S

( ) (

⁾^T ^0,

^{( )} (

⁾^T ⁰

Eéëx k w k ùû= Eéëx k v k ùû= S

(33)

(2.31)

再來定義協方差輸出矩陣：

(2.32) Ri為m x m 階；i 為時間延遲(Arbitrary Time Lag)，由下式推導得：

(2.33)

定義矩陣G 為狀態向量與輸出量測向量之協方差可得：

(2.34)

將上述三式展開可得：

(2.35)

(2.36) 上兩式可知輸出量測協方差矩陣與系統矩陣Ad與Cd之函數。

2.5 ^{隨機子空間識別法}

隨機子空間識別法由Peeters (1993)提出，是由 y(t)是唯一可經由儀器所量測 到的資訊，Ad與 Cd矩陣欲求之系統特性，因 Ad矩陣之特徵值中含有模態參數

( ) ( ) ( ) ( )

( ) ( ^{( )} ^{( )} )

( ) ( ) ( ) ( )

1 1 ^T

T T

T

E k k

E k k k k

E k k E k k

é ù

S = ë + + û

é ù

= êë + + úû

é ù é ù

= ë û× + ë û

= S +

d d

x x

A x w A x w

A x x A w w

A A Q

( ) (

[ 1 )^T

Ri =E y k+ y k ùû

( ) (

( ) ( )

( ) ( ^{( ) ( )} )

( ) ( ) ( ) ( )

0 [ )^T

T

T T T

T

R E k k

E k k k k

E k k E k k

= ùû

é ù

= êë + + úû

é ù é ù

= ë û× + ë û

= S +

d d

y y

C x v C x v

C x x C v v

C C R

( ) ( ) ( ) ( )

( ) ( ^{( ) ( )} )

( ) ( ) ( ) ( )

1 ^T

T

T T T

T

G

E k k k k

E k k

k

E k

k E

k é

û S

+ ù

é ù

= êë + + úû

é ù é ù

= ë û× + ë

=

= û

+ ë

d d

A x w C x v

A x x C w v

A C

x y

S

1 i

Ri =C A G_d _d^- ( 1)

T i T T

R_-i =G A_d^- C

(34)

要進行隨機子空間識別法，先需進行量測輸出值的重組，假設有 m 個輸出 量測數目，每個點位共有N 個時間反應歷時，由起始點擷取 j 點資料，依序每格 一點取一次共取 2i 次，對此取之資料矩陣 Y 稱一列區塊，並組成 Hankel 矩陣 Y0|2i-1，其表示如下：

(2.37)

量測資料同乘^#

$%，以利之後的計算，y(0)、y(1)、…、y(2i+j-2)為不同時間點 下之m x 1 輸出反應向量。建立量測輸出矩陣 Y0|2i-1，需要有N=2i+j-2 個時間點之量測資料，即矩陣Y0|2i-1行數j=N-2i+2，以達到統計上之穩定性。i 為自行定義參數，需符合矩陣階數要求，且需滿足mi >= n (系統階數)，即 i >= n/m。

2.5.1 協方差隨機子空間識別法

協方差隨機子空間識別法(Covariance-driven stochastic subspace identification，

簡稱SSI-COV)首次出現則在 Peeters (2000)的博士論文上。本研究採用此法，利

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )

0|2 1

0 1 1

1 2

1

1 1

1 2

2 1 2 2 2

i

j j

i i i j

j

i i i j

-

é - ù

ê ú

- + -

ê ú

= ê + + - ú

ê ú

+ + +

ê ú

- + -

ê ú

ë û

y y y

Y y y y

y y y

!

" " # "

!

" " # "

!

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )

( ) ( ) ( ) ( ) ( ) ( )

( ) ( ) ( )

0 1 1

1 2

1 1

1

1 2

2 1 2 2 2

j j

i i i j

j

i i i j

é - ù

ê ú

- + -

ê ú

ê + + - ú

ê ú

= ê + + + ú

ê ú

- + -

ê ú

ë û

y y y

!

" " # "

!

" " # "

!

(35)

用量測資料關係計算協方差矩陣，從中找出觀測矩陣，進而求得系統矩陣 Ad之識別方法。由上節提之輸出協方差矩陣組成Teoplitz 矩陣。

(2.38)

輸出協方差矩陣可改寫成：

(2.39) N 為量測反應資料點數，由上兩式可得之矩陣 T 可以改寫為：

(2.40)

Oi為觀測矩陣、為控制矩陣，將T 矩陣進行奇異值解耦(Singular Value Decomposition，簡稱 SVD)，可得：

(2.41)

上述兩式整理可得：

(2.42)

1 1

1 2

2 1 2 2

i i

i i i

R R R

- +

- -

é ù

ê ú

=ê ú

ê ú

ë û

!

" " # "

! T

( ) ( )

¹

( )

0

l m 1 1 (

1 ^T i ^N )^T

i N

k

R k

E k k N ^- k

®¥ =

é ù=

= ë^y + ^y û

å

^y + ^y

1 2

1

2 2 2 3 1

i i

i i i

- -

-

- - -

é ù

ê ú

=ê ú

ê ú

ë û

d d d d d

d d d d d d

C A G C A G C G

C A G C A G C A G

T

C A G C A G C A G

!

" " # "

!

(

¹ ²

)

1

i i

i

- -

-

æ ö

ç ÷

=ç ÷ … =

ç ÷

è ø

d d d

d d

C

C A A G A G G O Γ

C A

!

Γ_i

[

1 2

]

¹ ¹ 1

2 2

1 1

0

0 0

T

T T

T

é ù

= = êë » úû ëê úû=

S V

T USV U U U S V

S V

1 1 12 1

i

i-

æ ö

ç ÷

=ç ÷=

ç ÷

è ø

d d d

d d

C

C A U S

C A O !

(

ⁱ ¹ ⁱ ²

)

^1/2¹ ¹^T

i

- -

= _d _d … =

Γ A G A G G S V

(36)

從觀測矩陣Oi求得系統矩陣Ad：

(2.43)

其中為虛擬反矩陣。

2.6 ^{系統模態參數}

前一節已由輸出向量y(k)求得系統矩陣 Ad，系統模態參數隱含於Ad和Cd

矩陣之中。將Ad矩陣進行解耦，其特徵方程是可寫成：

(2.44) 其中γ 為特徵值矩陣(Eigenvalue Matrix)，Ψ 為特徵向量(Eigenvector)，分別表示 如下：

, (2.45)

Ad與Ac其關係如下：

(2.46)

由可解特徵值：

(2.47) ζ 與 ω 會成對出現，因此狀態數其值應保持為偶數，故階數需設定為兩倍 系統自由度。

( ( ) )

^†

( )

d = _i 1:l i-1 ,: _i l+1: ,:li

A O O

( )

! ^†

g

= AΨ Ψ

[

1 n

]

=

Ψ Ψ ! Ψ

1

n

l g

l

é ù

ê ú

= ê ú

ê ú

ë û

!

( )

ln

c t

d

d c

e

t

= D

= D A A

A A

c 1 1

0

- -

é ù

= êë- - úû A I

M K M C

λ_ci = −ξ_iω_i+ iω_i 1−ξ_i² λ_i =α_i+ iβ_i

ω_i = α_i²+β_i² = λ_i ξ_i= − α_i

α_i²+β_i² =−α_i λ_i

(37)

2.7 ^穩態圖

SSI 識別法主要應用於無法量測輸入項之微振案例，將真實量測之輸出資料進行重建系統模型，已求得系統參數。由於應用於真實案例時，系統階數為未知參數，必須分別針對系統模態參數萃取，繪製識別模態識別圖，稱為穩態圖 (Stabilization diagram)。

本研究假設模態數為樓層總數，一個樓層為一個模態，因此系統階數為已知，一個模態頻率即為一個極點(pole)，每一次分析時，會變化 SSI 中的 Hankel 矩陣，得到一組固定階數之模態頻率，將其各次分析繪製在圖上，此圖即為本研究所用之穩態圖。

2.7.1 ^{穩態圖實例}

上節提到穩態圖的理論，當量測資料越多組時，其系統建構越完整，其分析之系統參數也越精準，本研究進行了一個簡單的數值模擬分析，模擬一個五層剪力構架模型，其參數如表2-1，模態頻率之解析解如表 2-2。模擬外力為白噪音，

並使用 Newmark-β 法得分別五層樓之加速度時間歷時，將這五層時間歷時進行 SSI 分析得之穩態圖如圖 2-1。圖中紅色的線條為 FFT 所繪製之頻率域圖，藍色點為上述章節SSI 所分析之極點，從圖 2-1 可以看出各群之藍色極點與 FFT 之頻率峰值一樣。若將各群之極點分群並求得該值，可得之該穩態圖各模態頻率值如表 2-3，1 至 5 層樓時間歷時分析之模態頻率與解析解之模態頻率平均誤差為 0.31%。

(38)

圖2-1. 所有樓層(1 層至 5 層)時間歷時分析之穩態圖(藍點為各模態的極點，紅線為FFT 圖)

(39)

表2-1. 各層數值模型參數 5 層 10 層樓層質量 (kg) 5000 5000

樓層勁度 (N/m)

表2-2. 5 層數值模型之頻率解析解(Hz)

1^st 2^nd 3^rd 4^th 5^th Modal freq. 3.005 8.771 13.827 17.762 20.259

表2-3. 1 層至 5 層分析(Hz)

在現地實驗時，往往只擺設一顆感測器於最頂層，量測最頂層之時間歷時並進行分析，本研究模擬只擺設一顆感測器情況，取第五層加速度時間歷時進行 SSI 分析，所繪製之穩態圖如圖 2-2，並分群求得各群之模態頻率如表 2-4，在與解析解進行比較，其平均誤差為2.23%，如表 2-5 所示。

圖2-2. 單一層樓(5 層)分析之穩態圖(藍點為各模態的極點，紅線為 FFT 圖) 2.2 10´ 7 2.2 10´ ⁷

(40)

表2-4. 單獨 5 層分析(Hz)

表2-5. 5 層數值模型之 5 個模態頻率平均誤差比較(Hz)

1^st 2^nd 3^rd 4^th 5^th Error(%) 解析解 3.005 8.771 13.827 17.762 20.259

所有樓層分析 3.006 8.765 13.787 17.657 20.371 0.31 第5 層樓分析 3.006 8.721 13.544 17.638 21.842 2.23

從上述兩個穩態圖中，可得之當輸入資料越齊全，也就是系統的資訊越完整時，得到的系統參數與理論越接近，反之較不精準。因此進行現地實測時，能建置越多感測器進行量測，使用SSI 分析得到之模態頻率越精準。

2.8 ^{三次樣條內插法}

上節提到當系統建置的越完整時，其系統參數也越精確，但往往在進行現地量測時，會受限於各種情況，使得系統不夠完整，為了彌補系統上的缺少，補足為設置感測器之區域，Limongelli (2003)提出利用三次樣條內插法(Cubic spline interpolation method 簡稱 CSI)得到未設置感測器的樓層之時間歷時，並且求得感測器最佳化位置，得到該棟結構物之模態形狀。因此本研究使用 CSI 法進行內插。

CSI 藉由一組資料來預測、修正進行預測得到另一組新的資料。為了達到最佳之預測結果，通常會採取劃分區域方式，對各自區間進行預測。主要目的在於避免資料的變動性過大，導致預測不精確，高階多項式之預測方式是現今常見且有效的預測方式。

建築物微振量測之最佳感測器配置

國立臺灣大學工學院土木工程學系 博士論文

Department of Civil Engineering College of Engineering National Taiwan University

Doctoral Dissertation

建築物微振量測之最佳感測器配置

Optimal Sensor Placements for Ambient-Vibration Monitoring of Building

楊晏瑜 Yang, Yen-Yu

指導教授：呂良正 教授 Advisor: Prof. Leu, Liang-Jenq

中華民國 112 年 1 月

January 2023

致 謝

摘要

ABSTRACT

目錄

圖目錄

表目錄

第一章 緒論

1.1 研究背景

1.2 研究動機與目的

1.3 文章架構

第二章 系統識別方法

2.1 微振量測

2.2 快速傅立業轉換 (Fast Fourier Transform)

2.3 狀態空間

( )

( )

å

2.3.1 連續時間狀態空間方程式

( ) ( ) ( )

( ) ( ) ( )

( )

( ) ( )

( )

( )

( )

( )

( )

( )

( )

( )

2.3.2 離散時間狀態空間方程式

( ) ( )

( )

( )

( )

( )

( )

( )

( )

[ ]

( ) ( ) ( )

( ) ( ) ( )

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )

( ) ( ) ( ) ( )

( ) ( ) ( )

( )

( ) ( ) ( ) ( ) ( ) ( ) ( )

2.3.3 隨機狀態空間方程式

( ) ( ) ( )

( )

( )

ò

( )

( )

( )

ò

( )

( ) ( )

ò

(

)

[ ]

ò

(

)

( )

( )

( )

( )

( )

國立臺灣大學工學院土木工程學系博士論文

指導教授：呂良正教授 Advisor: Prof. Leu, Liang-Jenq

致 ^謝

第一章緒論

1.1 ^研究背景

1.2 ^{研究動機與目的}

1.3 ^文章架構

第二章系統識別方法

2.1 ^微振量測

2.2 ^{快速傅立業轉換} (Fast Fourier Transform)

2.3 ^狀態空間

2.3.3 ^{隨機狀態空間方程式}

2.4 ^隨機過程

( ) ( ^{( )}

^{( )} (

2.5 ^{隨機子空間識別法}

( ) ( ^{( )} ^{( )} )

( ) ( ^{( ) ( )} )

( ) ( ^{( ) ( )} )

2.6 ^{系統模態參數}

2.7 ^穩態圖

2.7.1 ^{穩態圖實例}

2.8 ^{三次樣條內插法}