陳衛(wèi)忠，王 慶

(蘇州市職業(yè)大學(xué) 數(shù)理部，江蘇 蘇州 215104)

1 問題提出

在中國遼闊美麗的大地上，古塔形成了一種具有中國民族傳統(tǒng)特色的建筑類型.由于長時間承受自重、氣溫、風(fēng)力等各種作用，偶然還要受地震、臺風(fēng)的影響，古塔會產(chǎn)生各種變形，諸如傾斜、彎曲、扭曲等.為保護(hù)古塔，文物部門需適時對古塔進(jìn)行觀測，了解各種變形量，以制定必要的保護(hù)措施.

本文參照2013高教社杯全國大學(xué)生數(shù)學(xué)建模競賽D題，某管理部門委托測繪公司于1986年7月、1996年8月、2009年3月和2011年3月對我國重點保護(hù)文物某古塔進(jìn)行的4次數(shù)據(jù)觀測[1]，提出下列3個問題：

問題1，給出確定古塔各層中心位置的通用方法；

問題2，分析該塔傾斜、彎曲、扭曲等變形情況；

問題3，分析該塔的變形趨勢.

本文對這3個問題，力求通過對觀測數(shù)據(jù)進(jìn)行分析，建立數(shù)學(xué)模型，并通過模型預(yù)測古塔的變形趨勢.

2 問題分析和模型假設(shè)

為解決以上3個問題，通過對所列數(shù)據(jù)的分析，對3個問題分別作分析和模型假設(shè).

對問題1，由于題目提供的觀測數(shù)據(jù)有缺失，考慮用外推法或幾何方法補充1986年和1996年缺失的數(shù)據(jù)(古塔第13層第5點).確定古塔各層中心位置時，考慮各層基本處于同一平面內(nèi)，可先擬合出各層所在平面，將各測量點投影到擬合平面內(nèi)，然后再用均勻物體的重心公式[2]計算中心坐標(biāo).

對問題2，要分析該塔傾斜、彎曲、扭曲等變形情況，可對各層中心點作線性擬合[3]，中軸線與水平面法向的夾角可作為傾斜程度的度量.對各層中心點作曲線擬合，曲線上各點曲率的平均值可作為彎曲程度的度量.相鄰兩個平面的旋轉(zhuǎn)角度可作為扭曲程度的度量.

對問題3，可利用問題2中的各種變形，關(guān)于時間作擬合，推測未來幾年的變化情況.

3 模型的建立和求解

3.1 對問題1模型的建立和求解

首先對缺失數(shù)據(jù)進(jìn)行處理，將1986年和1996年古塔的第13層其他點的觀測數(shù)據(jù)的平均值，作為缺失的第5點的數(shù)據(jù).通過Matlab作出圖1，分析觀測數(shù)據(jù)的分布.由圖1可以直觀地看出古塔各層的觀測點不在同一平面上，因此要求出古塔各層的中心位置坐標(biāo)，需先擬合出各層觀測點所在平面.設(shè)平面π：Ax+By+Cz+D=0，則該層的第i個觀測點Qi(xi，yi，zi)到平面π的距離.利用最小二乘法，求出使該層上所有觀測點到平面π的距離之和最小時的系數(shù)A、B、C，擬合出平面π.

圖1 觀測數(shù)據(jù)的分布圖

再將各觀測點Qi(xi，yi，zi)投影到平面π內(nèi)，得到投影點，求出的坐標(biāo).平面π的法向量為n=(A,B,C)，則向量，即，又因為π∈ ，所以，則，從而得

最后，古塔各層中心坐標(biāo)用均勻物體的重心公式計算.如圖2所示，n邊形A1A2…An，頂點坐標(biāo)A1(x1，y1)，A2(x2，y2)，…，An(xn，yn)，連接A2An，A3An，…，An-2An，得到n-2個三角形：△A1A2An，△A2A3An，…，△An-2An-1An，它們的重心為其中，面積為σi，i=1，…，n-2，則n邊形A1A2…An的 重心為G(x,y)，其中：

圖2 重心公式計算圖

3.2 對問題2模型的建立和求解

利用Matlab畫出每年古塔各層的中心位置分布圖，如圖3所示，發(fā)現(xiàn)古塔發(fā)生傾斜和彎曲，各層的中心點不共線、也不共面.

圖3 中心位置分布圖

利用問題1求擬合平面及投影點.首先求出到各層中心點距離之和最小的平面π'，再將第i層中心點Qi(xi，yi，zi)，i=1，…，13，投影到平面π'：A'x+B'y+C'z+D'= 0 上，得到投影點，其中.再在平面π'內(nèi)選取O1'作為坐標(biāo)原點，O1'與O2'的連線O1'O2'作為X軸，過O1'作O1'O2'的垂線作為Y軸，過O1'作平面π'的垂線作為Z軸，得到空間坐標(biāo)系O1'XYZ，則各層中心點在平面π'的投影點在空間坐標(biāo)系O1'XYZ的坐標(biāo)為

分析該塔傾斜程度.在平面π'上對各層中心點的投影點作線性擬合，平面π'的法向量與Z軸的夾角θ可作為傾斜程度的度量.設(shè)各層中心點的投影點在平面π'上擬合的直線方向向量為{l，m，n}，則

分析該塔彎曲程度.在平面π'上對各層中心點的投影點作曲線擬合，曲線上各點曲率的平均值可作為彎曲程度的度量.設(shè)擬合的曲線方程為y=f(x)，曲線上各點曲率，則古塔的彎曲程度的度量為

分析該塔扭曲程度.相鄰兩個平面的旋轉(zhuǎn)角度可作為扭曲程度的度量.將古塔相鄰兩層的觀測點投影到坐標(biāo)面上.如圖4所示，某個觀測點的扭曲可用α表示，則

圖4 觀測點投影圖

其中O表示第j層上第i個觀測點投影到坐標(biāo)面上的向量.求出所有觀測點的扭曲，取平均值作為古塔相鄰兩層的扭曲.

3.3 對問題3模型的建立和求解

為了分析古塔的變形趨勢，在問題2的基礎(chǔ)上將古塔傾斜角度、彎曲率、扭曲度關(guān)于時間作擬合，推測未來幾年的變化情況.

1)古塔傾斜角度的時間擬合.由圖5可見，從1986年到1996年，古塔的傾斜角增加較大；2011之后古塔的傾斜角變化不大.

2)古塔彎曲率的時間擬合.由圖6可見，從1990年到1996年，古塔的彎曲較嚴(yán)重，而1996年到2009年，古塔的彎曲率逐年減小，2011年之后古塔的彎曲變化不大.

3)古塔扭曲度的時間擬合.由圖7可見，從1986年到2011年，古塔的扭曲度逐年增大，預(yù)計之后幾年內(nèi)古塔的扭曲度趨于平穩(wěn).

圖5 古塔傾斜角度的時間擬合圖

4 結(jié)論

圖6 古塔彎曲率的時間擬合圖

圖7 古塔扭曲度的時間擬合圖

通過文中給出確定古塔中心位置的方法，對古塔傾斜、彎曲、扭曲等變形情況進(jìn)行度量分析，古塔變形的原因和趨勢是：隨著時間的累積，古塔有變形跡象、沉降趨勢，這與古塔的高度有關(guān)系，同時受外界自然風(fēng)力、地震等的影響，也直接影響到古塔的扭曲、彎曲情況，所以建議有關(guān)部門定期對古塔進(jìn)行檢測維修.這一模型除了用于對古塔等古建筑的維護(hù)外，還可以對現(xiàn)代高樓大廈的維護(hù)提供指導(dǎo)性方法[4].

[1]教育部高等教育司，中國工業(yè)與應(yīng)用數(shù)學(xué)學(xué)會. 中國大學(xué)生數(shù)學(xué)建模競賽[EB/OL]. (2013-09-09)[2014-02-09].http：//www.mcm.edu.cn/html cn/node/a l ffc4c5587c8a6f96eacefb8dbcc34e.html.

[2]向愛琴. 質(zhì)面多邊形重心公式的再探討[J]. 新疆教育學(xué)院學(xué)報，1989(1)：58-60.

[3]龔楊. 空間直線擬合的一種方法[J]. 齊齊哈爾大學(xué)學(xué)報，2009(2)：64-67.

[4]梁海奎. 古塔變形測量方法探討[J]. 城市勘探，2011( 3)：113－118.