湯 志 浩

(河南工業(yè)職業(yè)技術(shù)學(xué)院 基礎(chǔ)部， 河南 南陽 473009)

?

基于Matlab的古塔傾斜數(shù)據(jù)分析及監(jiān)測

湯 志 浩

(河南工業(yè)職業(yè)技術(shù)學(xué)院 基礎(chǔ)部， 河南 南陽 473009)

利用Matlab對某古塔1986年、1996年、2009年和2011年的觀測數(shù)據(jù)進(jìn)行了分析，依據(jù)觀測數(shù)值，找出了古塔每層的中心點(diǎn) ，對古塔的扭曲、彎曲、傾斜等進(jìn)行了分析，得出了隨著時(shí)間的推移，古塔的變形程度會(huì)進(jìn)一步增加的結(jié)論.

古塔； 傾斜； 彎曲； 扭曲

文物是國家不可再生的文化資源，尤其是有上千年歷史和文化的古塔，更是我們國家需要重點(diǎn)保護(hù)的古建筑.長期以來因各種因素的影響，古塔可能發(fā)生如傾斜、彎曲、扭曲等形變.對古塔的這些變形數(shù)據(jù)進(jìn)行分析和預(yù)測，對我國文化遺產(chǎn)的保護(hù)具有重要的意義.本文試用Matlab對這些數(shù)據(jù)進(jìn)行分析.

1 問題提出

長期以來，古塔因?yàn)樽陨碇亓?、寒暑交替、風(fēng)化等各種因素的影響，偶爾還可能遭受地震、颶風(fēng)等自然災(zāi)害的破壞，可能發(fā)生如傾斜、彎曲、扭曲等各種損壞性形變.為了保護(hù)沉淀了上千年厚重歷史、文化的古塔，文物管理部門會(huì)定期觀測古塔，記錄不同時(shí)期的變形數(shù)據(jù)，以根據(jù)數(shù)據(jù)擬定合適的保護(hù)方案.

某地區(qū)有一年代久遠(yuǎn)的古塔，屬于國家重點(diǎn)保護(hù)的古建筑.某測繪公司按照文物管理部門的要求對該塔進(jìn)行了4次觀測，分別是1986年7月、1996年8月、2009年3月和2011年3月.依據(jù)對古塔四年的觀測數(shù)據(jù)，本文將對三個(gè)問題進(jìn)行探討.

問題1：依據(jù)觀測數(shù)值，找出古塔每層的中心點(diǎn)，并給出通用的方法，且將古塔每層中心點(diǎn)的坐標(biāo)用表格列出來.

問題2：根據(jù)觀測值對古塔的扭曲、彎曲、傾斜等變形程度進(jìn)行分析.

問題3：結(jié)合問題1和問題2，分析和預(yù)測古塔的變形趨勢和方向.

2 問題分析

依據(jù)1986年觀測數(shù)據(jù)，在Matlab軟件中做出三維曲線圖[1]，見圖1.

圖1 由1986年數(shù)據(jù)得到的三維曲線

從圖1可以看出，該古塔為八角形，共有13層，外加一個(gè)塔尖層，各層的八個(gè)觀測點(diǎn)剛好取在古塔的八個(gè)角上.

2.1 模型的假設(shè)

(1)古塔初建成時(shí)，假設(shè)每層都與地面平行，八個(gè)測試點(diǎn)均在同一高度(z軸坐標(biāo)相同)，且都在每層的地面上.

(2)在條件(1)的前提下，假設(shè)古塔的每層中心均在同一z軸上，即每層中心的x,y坐標(biāo)相同(因?yàn)橥粚影藗€(gè)觀測點(diǎn)的z值變化很小).

(3)假設(shè)在建塔時(shí)塔身嚴(yán)格垂直于大地的水平面，則分析該塔時(shí)可把塔看作標(biāo)準(zhǔn)的八臺(tái)，所以在建造時(shí)塔的八條棱等長.

(4)假設(shè)各層的觀測點(diǎn)都取在古塔的地面上，而不是在中間的某個(gè)點(diǎn).

(5)假設(shè)古塔初建成時(shí)，沒有扭曲、彎曲、傾斜、下沉等現(xiàn)象.

(6)本文均以1986年的觀測數(shù)據(jù)作為基準(zhǔn)考慮.

2.2 問題1模型的建立與求解

該部分要求依據(jù)觀測數(shù)據(jù)，將古塔每層的中心點(diǎn)找出來，并給出通用的方法，且將古塔每層中心點(diǎn)的坐標(biāo)用表格列出來.本文取1986年第一層的八個(gè)觀測數(shù)據(jù)進(jìn)行具體計(jì)算和分析，數(shù)據(jù)見表1.

表1 1986年第一層八個(gè)點(diǎn)的橫縱豎坐標(biāo) m

由表1可知，z坐標(biāo)的數(shù)據(jù)基本接近，因此，可將1986年第一層中心位置的豎(z)坐標(biāo)取這八個(gè)點(diǎn)的平均值，為1.78700，然后用Matlab軟件做出這八個(gè)觀測點(diǎn)的x,y(即橫縱坐標(biāo))的散點(diǎn)圖[1]，如圖2所示.

圖2 1986年觀測到的第一層數(shù)據(jù)的散點(diǎn)圖

圖1的散點(diǎn)圖，可用一個(gè)圓來擬合，如圖3所示.

圖3 第一層數(shù)據(jù)擬合圖

在Matlab軟件中調(diào)用circfit函數(shù)[2]可求出擬合圓的圓心的橫縱坐標(biāo)和半徑，由表1的數(shù)據(jù)得到的結(jié)果如下：

xc=566.65180，

yc=522.69640，

R=2.83090.

其中，xc表示1986年第一層數(shù)據(jù)擬合圓圓心的橫坐標(biāo)；yc表示1986年第一層數(shù)據(jù)擬合圓圓心的縱坐標(biāo)；R表示1986年第一層數(shù)據(jù)擬合圓的半徑.

因此，由表1可以得出，1986年古塔第一層中心點(diǎn)的橫縱豎坐標(biāo)為：(566.65180，522.69640,1.78700).1986年所有觀測點(diǎn)的平面俯視圖如圖4所示.

圖4 1986年觀測點(diǎn)的平面俯視圖

觀察圖4可以看出，每一層上的點(diǎn)都可以用一個(gè)圓來擬合，因此，對于1986年觀測的其他12層和塔尖，以及另外三年中的觀測數(shù)據(jù)全部可以采用上面的方法將古塔各層中心點(diǎn)的坐標(biāo)確定出來.表2-表5給出了四年中各層中心點(diǎn)的坐標(biāo)數(shù)據(jù).

表2 1986年每層中心點(diǎn)的橫縱豎坐標(biāo) m

表3 1996年每層中心點(diǎn)的橫縱豎坐標(biāo) m

表4 2009年每層中心點(diǎn)的橫縱豎坐標(biāo) m

表5 2011年每層中心點(diǎn)的橫縱豎坐標(biāo) m

2.3 問題2模型的建立與求解

2.3.1 傾斜

傾斜是指整體發(fā)生偏離豎向或者水平方向的偏斜.對于本題來說就是古塔偏離豎直方向的角度[3]，圖5為1986年中心點(diǎn)的散點(diǎn)圖.

圖5 1986年中心點(diǎn)的散點(diǎn)圖

由圖5可以看出，中心點(diǎn)基本在一條直線上，因此，本文以古塔中心線偏離豎直方向的角度作為古塔的傾斜角度.先取塔尖和第一層的中心點(diǎn)來進(jìn)行探討.

再次，把四年的測量數(shù)據(jù)分別代入公式可得古塔四次測量時(shí)的傾斜角，見表6.

表6 四次觀測時(shí)古塔的傾角

2.3.2 彎曲

構(gòu)件(經(jīng)常是線狀構(gòu)件)在外力作用下與中軸線偏轉(zhuǎn)的狀態(tài)稱為彎曲.

先將問題1求出的結(jié)果即四次觀測的各層中點(diǎn)在Matlab中擬合出曲線，并將其方程列出，然后用曲率公式算出四種情況下的曲率作為古塔的彎曲值(曲率是反映曲線偏離直線程度的量).

在Matlab中用三元線性回歸[4]擬合出四條曲線如圖6所示.

圖6 四年中八個(gè)中心點(diǎn)的擬合圖

曲線方程如下：

f1(x)=p1*x^5+p2*x^4+p3*x^3+p4*x^2+p5*x+p6.

p1=-4.476，p2=1.068e+004，

p3=-9.828e+006，

p4=4.279e+009，p5=-8.464e+011，p6=5.438e+013.

f2(x)=p1*x^5+p2*x^4+p3*x^3+p4*x^2+p5*x+p6.

p1=-6.88，p2=7003，p3=4.149e+006，

p4=-8.16e+009，p5=3.697e+012，p6=-5.5e+014.

f3(x)=p1*x^5+p2*x^4+p3*x^3+p4*x^2+p5*x+p6.

p1=-1.182，p2=3545，p3=-4.239e+006，

p4=2.528，p5=-7.518e+011，

p6=8.923e+013.

f4(x)=p1*x^5+p2*x^4+p3*x^3+p4*x^2+p5*x+p6.

p1=-4.476，p2=1.068e+004，p3=-9.828e+006，p4=4.279e+009,p5=-8.464e+011，p6=5.438e+013 .

曲率能體現(xiàn)出曲線對直線的偏離程度.在數(shù)學(xué)上通常用來反映曲線在某個(gè)點(diǎn)處的彎曲幅度.曲率取值越大，曲線彎曲越厲害.因此，本文可以通過計(jì)算四條曲線的曲率，使用Matlab軟件計(jì)算出各層的彎曲角.結(jié)果見表7.

表7 各塔層三次測量的轉(zhuǎn)角 弧度

2.3.3 扭曲

幾何體因?yàn)橥饬ψ饔冒l(fā)生橫截面方向的順(逆)時(shí)針扭轉(zhuǎn)的變形狀態(tài)稱為扭曲，從第一層開始幾乎每一層和其上一層都會(huì)有微小的傾角，這些小的角度加起來就是扭轉(zhuǎn)角度，圖7為1986年觀測點(diǎn)對角線連線的俯視圖，如果沒有扭轉(zhuǎn)，這些對角線的俯視圖應(yīng)該重合在一條線上，而這些對角線連線的俯視圖呈分散狀態(tài)，說明古塔的扭曲可以用圖中八個(gè)分支的角度中最大的那個(gè)角度作為古塔的扭曲度.而且，隨著年份的增加，古塔的扭曲程度會(huì)越來越大.

圖7 1986年觀測點(diǎn)對角連線俯視圖

通過在Matlab中對數(shù)據(jù)的處理，給出了四年中對角線在水平面上投影產(chǎn)生的傾斜度的最大值，即古塔的扭曲度,見表8.

表8 古塔扭曲角度 弧度

2.4 問題3求解

問題3為分析該古塔的變形趨勢.從問題2的求解可以看出，古塔在四次測量傾角時(shí)，與z軸的夾角都達(dá)到了4.2度，這個(gè)角度隨著各種因素的作用會(huì)進(jìn)一步變大.由圖7可以看出，整個(gè)塔身呈現(xiàn)順時(shí)針旋轉(zhuǎn)的扭曲狀態(tài)，四年中古塔在順時(shí)針方向上轉(zhuǎn)動(dòng)的角度由0.0135變?yōu)?.0139，扭轉(zhuǎn)角度增加.而且由圖4可知，觀測點(diǎn)在左上角分布稀疏，右下角分布密集，這表明古塔在東南角有下沉的變化趨勢.從圖6對四年古塔的各層中心點(diǎn)的描繪可以看出，古塔的彎曲、扭曲、傾斜在各層、各個(gè)方向都存在.因?yàn)殚L期受自身重力、寒暑交替、風(fēng)化等各種因素的影響，偶爾還要遭受地震、颶風(fēng)的破壞，古塔的彎曲角度、扭曲角度、傾斜角度、下沉趨勢，如果不人為保護(hù)的話，都會(huì)進(jìn)一步變大.

圖8是以1986年的數(shù)據(jù)為例做出的古塔第三層底面的三維圖，它進(jìn)一步證明了對古塔的扭曲、彎曲、旋轉(zhuǎn)、傾斜、下沉的變化趨勢分析是正確的.

3 模型的評價(jià)與改進(jìn)

(1)在問題1中用圓來擬合，取圓的中心作為古塔各層的中心位置，這在將來的研究中可以采用橢圓擬合，這樣會(huì)更具有一般性，因?yàn)閳A是橢圓的特例.

(2)運(yùn)用Matlab和Excel軟件畫圖，進(jìn)行函數(shù)調(diào)用和編程求解，產(chǎn)生的誤差較小，結(jié)果準(zhǔn)確、合理.

圖8 1986年古塔第三層底面三維圖

(3)采用函數(shù)擬合法建立的模型原理簡單，符合實(shí)際，簡化了算法，計(jì)算方便.

(4)所建模型和使用的方法通用性強(qiáng)，對類似的問題均有較強(qiáng)的參考意義.

(5)在計(jì)算中心坐標(biāo)和扭轉(zhuǎn)角度時(shí)充分運(yùn)用了垂直投影的數(shù)據(jù)處理方法.

[1] 顏文勇.數(shù)學(xué)建模[M].北京：高等教育出版社，2011：160-163.

[2] 張志涌,楊祖櫻.Matlab教程[M].北京：北京航空航天大學(xué)出版社，2010：130-135.

[3] 胡志曉.古塔傾斜觀測和數(shù)據(jù)分析[J].江蘇建筑,2011(6):34-44.

[4] 韓中庚．?dāng)?shù)學(xué)建模方法及其應(yīng)用[M].北京：高等教育出版社，2009:122-123.

責(zé)任編輯：金 欣

The Matlab-based data analysis of an ancient tower

TANG Zhi-hao

(Department of Basic Teaching, Henan Polytechnic Institute, Nanyang, Henan 473009, China)

This study is based on the Matlab-based data analyses of a tower in 1986, 1996, 2009 and 2011 in relation to the center point of each of the ancient tower, distortion of the ancient tower, bending, and tilting.

ancient tower; tilting; bending; distortion

2015-05-27

湯志浩 (1980-)，女，湖北十堰人，河南工業(yè)職業(yè)技術(shù)學(xué)院講師，碩士，研究方向：應(yīng)用數(shù)學(xué).

1009-4873(2015)04-0015-06

O242.1

A