基于空間直線位置關系的古塔變形研究

摘? 要： 針對古塔變形問題，給出了一種基于極限思想的古塔各層中心的計算方法;利用空間直線的位置關系分析了古塔的各種變形，并利用MATLAB數值實驗對模型進行了計算從而得到了數值結果。

關鍵詞： 多邊形的中心;極限思想;空間直線的位置關系;古塔變形;MATLAB

中圖分類號： O29? ? 文獻標識碼： A? ? DOI：10.3969/j.issn.1003-6970.2020.10.062

本文著錄格式：游晉峰. 基于空間直線位置關系的古塔變形研究[J]. 軟件，2020，41（10）：246249

【Abstract】： In this paper， the deformation problem of ancient pagoda is studied. A calculation method of the center of each layer of ancient pagoda based on limit thought is given; the various deformations of ancient pagoda is analyzed by using the position relationship between spatial lines， and the model is calculated by using MATLAB numerical experiment， and the numerical results are obtained.

【Key words】： The center of polygon; Limit thought; The position relationship between spatial lines; The deformations of ancient pagoda; MATLAB

0? 引言

古塔長時間承受自重、氣溫、風力等作用和地震、颶風的影響，會產生傾斜、彎曲、扭曲等各種變形。2013年全國大學生數學建模競賽C題“古塔的變形”[1]，給出了1986年7月、1996年8月、2009年3月和2011年3月4次觀測數據，要求給出確定古塔各層中心位置的通用方法，并分析古塔的變形情況。

陳俊等人[2-4]利用投影多邊形模型來確定古塔各層的中心，采用最小二乘法建立了線性模型，并借助三維高次曲線方程建立古塔的曲率模型和撓率模型。但這種用中軸變形曲線刻畫物體的整體扭曲是有缺陷的[5]。蔡志杰等人[5]將古塔看作連續體，利用位移函數來刻畫古塔的變形。其中所用知識均涉及了較專業的高等數學知識。而數學建模應遵循一個原則：盡量采用簡單的數學工具[6-7]。

本文首先分析了題目所給數據，發現數據觀測起止點不同。于是，對數據進行了處理，同時繪制了古塔的三維網面圖和俯視圖，便于進一步觀察古塔。其次，給出了一種基于極限思想的古塔各層中心的計算方法。然后，利用空間直線的位置關系分析了古塔各種變形，利用Matlab軟件完成了對模型的求解，得到了數值結果。

1? 數據分析

通過分析題目所給的觀測數據發現，1986年和1996年觀測點的起止一致，2009年和2011年觀測點的起止一致。但前兩次觀測點的起止和后兩次的不同。

如圖1，對1986年第一層、1996年第二層、2009年第三層、2011年第四層的數據作圖可以看到，1986年和1996年的起止一致，2009年和2011年的起止一致。但前兩次觀測和后兩次觀測的起止點不同。

為方便后續的計算和變形分析，需要對數據進行處理，使四次觀測數據的起止點一致。

1.1? 數據準備

對2009年和2011年觀測點的起止順序進行修改，即，將每層的第7個觀測點、第8個觀測點的數據依次移動到第1個觀測點的數據之前。如圖2可以看到，數據修改后，4次觀測數據的起止保持一致。

1.2? 古塔觀測圖

利用題目觀測數據，繪制了古塔的三維網面圖和俯視圖，便于進一步觀察古塔。由于1996年、2009年和2011年三次觀測數據的作圖結果和1986年的類似，故此處只給出1986年的古塔三維網面圖（圖3）和俯視圖（圖4）。

從1986年古塔的俯視圖（如圖4）可以看出，古塔各層的結構近似于正八邊形，塔尖幾乎投影到了一點。

2? 基于極限思想的古塔各層中心的計算方法

2.1? 確定古塔在xoy平面上各層的中心

對于某次測量數據的第層（）在x軸和y軸的投影，求該多邊形的各邊的中點，并依次連線形成新的多邊形，再求新的多邊形的各邊的中點，……，如此循環多邊形會變小，直到該多邊形近似趨于一點，停止循環迭代。此時，得到的點為該多邊形的中心。

算法如下：

%賦初值

x_old=0;

y_old=0;

x_new=x（i，：）;

y_new=y（i，：）;

%迭代循環，直至兩次坐標差值比小

while （abs（x_old-x_new）>0.0000005） & （abs（y_ old-y_new）>0.0000005）

x_old=x_new;

y_old=y_new;

for j=1：size（x_new，2）-1

x_new（j）=（x_new（j）+x_new（j+1））/2;

y_new（j）=（y_new（j）+y_new（j+1））/2;

end

x_new（size（x_new，2））=x_new（1）;

y_new（size（y_new，2））=y_new（1）;

end

需要注意的是，由于頂層數據嚴重缺失，只考慮第1～13層。

本文中，考慮允許的觀測誤差為[8]，故設置當兩次迭代的值相差小于時，停止迭代。如圖5，對1986年的第一層的中心坐標利用上述算法進行計算，可以看到，多邊形逐漸趨于一點。

2.2? 確定古塔在z軸上各層的中心

2.1只考慮了古塔在xoy平面的中點，并未考慮z軸方向上的中心坐標。

利用求某次測量的第層（）的z軸方向上的坐標均值，來確定古塔第層（）的z軸方向上的中心坐標。

2.3? 確定的古塔各層中心坐標

確定古塔各層中心坐標的最終方案為：

（1）利用無窮逼近的思想，確定古塔第層（）x軸和y軸方向上的中心坐標;

（2）利用z軸方向上的坐標均值，確定古塔第層（）的z軸方向上的中心坐標。

（3）利用四次測量結果的均值確定1986年和1996年塔尖的中心坐標。

（4）直接利用題目所給的一次測量結果確定2009年和2011年塔尖的中心坐標。

最終確定的四次測量的古塔各層中心坐標，見表1。

3? 分析古塔的變形情況

針對古塔變形的情況，從古塔的墻體是否有裂縫或擠壓、古塔整體的傾斜量和傾斜角、古塔各層的沉降[9]、古塔的扭曲四個角度展開研究。

3.1? 古塔墻體的擠壓與裂縫

檢測墻體線段的長度，與前一年數據比較，根據發生的變化判定是否發生擠壓或裂縫使墻體變形。

各層各線段的長度與前一年數據比較，結果如下。

通過分析圖6～8發現，1996年到2009年測量數據的波動較大，即出現墻體有裂縫（圖中正數）和壓縮變形（圖中負數）的現象。其他兩個的比較數量級都在之間。

3.2? 傾斜量和傾斜角

利用2.3節得到的古塔各層中心坐標計算傾斜量和傾斜角。

古塔中心坐標的曲線如圖9。

3.2.1? 整體的傾斜量和傾斜角

利用第1層的中心坐標和塔尖的中心坐標連線形成的直線的傾斜量和傾斜角，估計整體的傾斜量和傾斜角。

整體傾斜量的計算公式如下：

整體傾斜角的計算公式如下：

利用Matlab軟件計算得到：

3.2.2? 彎曲，即各層的傾斜量和傾斜角

對于第層（）的傾斜量和傾斜角，利用第層的中心坐標和第層的中心坐標連線形成的直線進行計算。

第層（）傾斜量的計算公式如下：

第層（）傾斜角的計算公式如下：

通過計算得到：

第1層～第13層的傾斜量依次為：0.0493、0.0472、0.0392、0.0411、0.1015、0.0390、0.0608、0.0610、0.0768、0.0580、0.0578、0.0660、0.0862。

第1層～第13層的傾斜角依次為：0.5113、0.4977、0.5188、0.5062、1.2924、0.6178、0.9901、1.0028、1.3249、0.7764、0.7780、0.9163、2.1702。

各層的傾斜量和傾斜角折線圖，如圖10和圖11。

從圖10、11可以看到，第5、9、13層的傾斜量和傾斜角相對比較明顯，因此要對第5、9、13層加強防護。

3.3? 古塔各層的沉降

利用2.3節得到的古塔各層中心的z軸方向的坐標值、、、的差值

分析古塔的各層沉降。

例如：比較1986年和1996年z軸方向的差值，計算公式如下：

1986年和1996年之間、1996年和2009年之間、2009年和2011年之間各層的沉降為見表2（負數表示下沉）。

觀察表2發現，古塔各層均有不同程度的下沉，下沉范圍在0.0013 m～0.0288 m之間。

3.4? 古塔的扭曲變形

利用題目所給的觀測數據，通過作圖觀察四次測量的第層（）的扭曲變形程度，見圖12。

通過觀察圖像發現，第1、2、3層四年之間幾乎沒發生扭曲，第4、5、7～13層四年之間發生了輕微的的扭曲，第6層1996年和2009年之間的扭曲較嚴重。

下面以第六層1996年和2009年的數據為例，計算該層此兩年間的扭曲變形情況，如圖13。

對于第條線段的扭曲度，其計算公式如下：

計算得到，對應的八條線段的扭曲度依次為：1.1823、1.6083、0.3730、0.3603、0.8891、1.8980、2.1883、0.1329。

4? 結語

數學建模講求盡量采用簡單的模型，不盲目追求所謂的“高大上”的方法和“全能新”的模型[10]。本文利用初等方法解決古塔變形問題，提出了利用取多邊形各邊中點的方法無限逼近古塔各層中心，從而得到古塔各層的中心坐標;利用空間直線的位置關系，分析了古塔可能的各種變形情況。對于古塔變形，不僅分析了題目要求的傾斜、彎曲、扭曲等變形，而且分析了題目古塔墻體是否有裂縫、古塔各層的沉降。該數據結果與陳俊的數據結果雖略有差距，但對于初接觸數學建模的學生較容易理解。

參考文獻

[1]2013年高教社杯全國大學生數學建模競賽賽題. [DB/OL]. （2013-09-13）[2020-07-18]. http：//www.mcm.edu.cn/problem/ 2013/2013.html.

[2]陳俊. 古塔的形變研究[J]. 數學建模及其應用， 2014年2月第3卷第1期： 67-76.

[3]劉世杰. 基于多層中軸點擬合的古塔變形檢測[J]. 同濟大學學報（自然科學版）， 2018年3月， 第46卷第3期： 401-405.

[4]侯學慧. 古塔的變形問題研究[J]. 山西師范大學（自然科學版）， 2019年3月， 第33卷第1期： 47-53.

[5]蔡志杰， 譚永基. 用測量數據分析古塔的變形[J]. 數學建模及其應用， 2013年11月， 第2卷第5-6期： 14-19.

[6]姜啟源， 謝金星， 葉俊. 數學模型（第四版）[M]. 北京： 高等教育出版社， 2011. 1， 第14頁.

[7]徐仁旭， 孔亞仙. 數學實驗與建模[M]. 長沙： 湖南師范大學出版社， 2011. 6.

[8]百度文庫. 建筑施工測量技術要求及允許偏差. [DB/OL]. （2010-11-15）[2020-07-18]. http：//wenku.baidu.com/view/ 7c5c5f1aa8114431b90dd854.html.

[9]搜狐新聞. 地震加速西安大雁塔變形沒有造成破壞性影響. [DB/OL]. （2008-05-26）[2020-07-18]. http：//news.sohu.com/ 20080526/n257078702.shtml.

[10]韓中庚. 數學建模競賽論文的寫作方法[J]. 數學建模及其應用， 2017年6月， 第6卷第2期： 42-48.