古塔的变形问题研究

侯学慧

新疆警察学院信息安全工程系， 新疆 乌鲁木齐 830011

自东汉初年建成第一座塔约两千年以来，中华大地上出现的塔数逾万，虽然历经自然淘汰和人为损毁，据说尚余三千，仅此数也号称世界之最.古塔由于年代久远，长时间承受自重、气温、风力等各种作用，偶然还要受地震、飓风和洪涝的影响，古塔会产生各种变形，诸如倾斜、弯曲和扭曲等[1～4].本文探究古塔受外界环境影响而产生变形、弯曲、扭曲等情况并将其量化得到具体的数据以便预测未来许多年内，古塔的变形情况或者考虑古塔新建立未变形时的具体位置和高度等.

开封铁塔已有上千年历史，是我国重点保护文物.管理部门委托测绘公司先后于1986年7月、1996年8月、2009年3月和2011年3月对该塔进行了4次观测(图1).

根据图1提供的4次观测数据，本文主要讨论以下3个问题，并通过数据理论分析该古塔的变形情况：

(1)给出确定古塔各层中心位置的通用方法，并列表给出各次测量的古塔各层中心坐标；

(2)分析该塔的倾斜、弯曲和扭曲等变形情况；

(3)分析该塔的变形趋势.

1 问题分析

中国古塔成近千年，历尽沧桑，且不具备再生性.近千年来遭受到地震、冰雹、风灾和水患等自然灾害，尤其是在战争时期曾用飞机、大炮进行轰炸等一系列人为因素，历经沧桑，不少古塔身患重病、年久失修，导致古塔发生了各种形变，例如倾斜、弯曲、扭曲、下沉、构件断裂等，但这里我们目前只考虑倾斜、弯曲、扭曲这3种主要的形变.

1.1 确定古塔各层中心位置的通用方法

确定三种解决方案：

图1 用Mathmatic 计算每一层中心(圆心)坐标
Fig.1 Use Mathmatic to calculate the center coordinate of each layer

(1)我们利用Mathmatica作图可知，给出的每层8个点是出于顶点的位置，因此，我们可以用8个点的平均值点来代替中心位置(采用平均值法).

(2)中国古代建筑是对称结构，塔面多数采用偶数边的正多边形.因此采用求外接圆的圆心代替中心.

(3)用软件计算出圆心.

1.2 古塔的变形程度

对于变形问题我们只考虑3个主要因素：倾斜，扭曲和弯曲.解决时量化这3个量，求其变化角度.注意：给出的数据在2009年到2011年的坐标对应发生了变化，需要调整.

倾斜程度：构造每层中心点与古塔初建造时地面中心点的直角三角形，在直角三角形当中求解倾斜角.

扭曲程度：我们利用软件量化扭曲程度，将古塔同一标号的13个点在地面上的投影做函数拟合成直线(或者用一元线性回归方程)，从而分析拟合后的投影直线与变化后的投影直线的夹角，即为扭转角度.

1.3 分析古塔的变化趋势

我们采用数据直线拟合的方法，根据古塔的变化角度，分析其变化量，建立拟合函数方程，根据古塔的变形程度得到的数据，确定函数，从而预测n年以后，古塔的变化趋势.

2 古塔变形问题的解决

首先假设：

(1)古塔的倾斜程度，即倾斜角记为θ；

(2)古塔的扭曲程度，即倾斜角记为α；

(3)古塔同一标号位置13个点在地面上投影点集，拟合成的直线记为L.

2.1 确定古塔各层中心

2.1.1 确定古塔各层中心位置的通用方法

(1)方法一：求平均值法(用8个点的平衡点来代替中心点)

我们采用近似值算法，因为中国古代的古塔多采用对称结构，多数塔层都是由偶数边的正多边形构成的，这样的正多边形是相对对称均匀的.因为已知点的坐标，我们通过Mathmatice作图得到图2，其中给出的点是相对的顶点.因为已知的8个点是对称的，我们用这8个点坐标的平均值来代替每层塔面的中心位置.

(2)方法二： 将内切圆的圆心视为每一层塔面的中心点

用平均值求近似值，相对的误差较大，可改为采用内切圆求其圆心.

因为已有的8个点是在顶点位置，古代对于一个建筑物的建造总是以圆为内切圆或者外接圆建立的正多边形，因此可以认为这8个点一定是在某一个圆心为(a,b)，半径为r的圆上.这样问题就转化为已知圆上8个点的坐标，求出圆心即为中心.

图2 Mathmatice软件作出的古塔顶点在地面的投影Fig.2 To make the Graph of the vertex projection in the Ground

图3 近似八角古塔
Fig.3 Approximate Eight-Square

(3)利用C语言中的Mathmatica软件编写程序

求解问题转化为：已知图上的多个点的坐标，计算该圆的圆心与半径 ，从而确定正多边形的中心.

2.1.2 列表给出各层量古塔各层的中心坐标

(1)依据上述方法，利用EXCLE、Mathmatic软件，得到各层中心点的列表.如表 1所示.

(2)利用圆心去近似代替中心，我们给出1986年第一层的中心的计算方法，并得到中心坐标为(566.664 8,522.710 5,1.787 4)，其余各层方法相同.

2.2 分析古塔的倾斜、弯曲和扭曲等变形情况

2.2.1 古塔的倾斜程度分析

坐落在开封市铁塔是一座历史悠久的古塔，因此在1986年时，由于各种环境等因素，使古塔本身发生变形，造成倾斜.由于我们不了解古塔的背景与建造时间，我们只重点考虑1986年以后，并且假设在1986年时古塔的变形情况不严重，因此我们假设，将1986年古塔第一层中心的投影点，当作古塔建造时塔地面的中心O点(存在误差).

表1 给出确定古塔各层中心坐标Tab.1 Give Tower the center coordinates of each layer

那么我们可以近似用坐标O(566.664 8, 522.710 5,0)(第一层的中心)来代替建筑物的塔底中心点.由于倾斜后，每层的中心点与地面的中心点O近似构成相似的直角三角形,那么我们可以在直角三角形当中，记O1(x1,y1,z1)，O0(x1,y1,z1),从而通过公式计算倾斜角公式

已知4年的塔顶坐标，从而近似确定这4年的塔顶中心坐标，因1986年和1996年的塔顶4个坐标都存在，因此直接求出平均中心点作为塔顶中心，而2009年和2011年的塔顶4个坐标不全，因此近似用第13层的中心点来作为塔顶的中心轴位置，我们根据算出的四年的塔顶中心点与塔底中心点的连线所在的直线与竖直方向存在倾斜角θ，利用EXCEL计算每一年的塔顶中心点坐标(表2)，根据倾斜角公式，我们近似得出1986年～2011年相应的倾斜角.表3给出倾斜角的角度值.

表2 每一年的塔顶中心点坐标Tab.2 Vertex coordinate of the top of tower in each year

表3 古塔倾斜角的度数Tab.3 Degree of inclination of tower

2.2.2 古塔的扭曲程度分析

因为我们对扭曲程度的变化，注重其点在地面上的投影，因此建立直角坐标系XOY，扭曲程度可分为2个方面：为了说明扭曲过程，我们都以1986年为基准.看1996年～2011年发生了多大程度的扭转.

(2)另一方面是1986年～2011年间古塔发生的扭转程度.我们选择计算出同一年当中同一个标号位置的不同层面的13个点在地面上投影所成的近似直线，考虑变化后直线的倾斜角.

这是标号均为1位置的1986年的13个顶点，画出在地面上的投影，用Mathmatice输入坐标，求出拟合的一次函数记为L1:y=k1x+b1；得到k1=-3.134 11.同上，我们选取标号为1的1996年的13个顶点在地面上的投影，拟合直线为L2:y=k1x+b1；得到k2x+b2,k2=-3.122 54，那么我们计算出L1,L2的夹角.

通过计算得出为1986年～1996年的第一个位置扭转的角度为0.061 5°.用相同的方法得到1986年～2011年对应的每个顶点的扭转度.在第一问题当中，同一年份的扭转角度算法相同，由于计算量较大，就不一一列出.

(3)弯曲程度

我们在考虑古塔弯曲时，中心轴会出现弯曲，可以采用根据已知点的坐标，拟合函数，计算中心轴的曲率，由Mathmatica作图可知，具体计算出较为明显弯曲的位置.也可以采用求解扭转角度的办法，计算拟合直线的夹角.在理想状态下，棱上的13个顶点，应在同一直线上，但是一旦出现明显弯曲，那么任意两点确定的直线应该会出现夹角，通过夹角判断其弯曲的具体位置.

2.3 近似用线性关系分析未来古塔的变形情况

倾斜程度的分析：1986年～1996年这十年，倾斜角增加了0.000 18°,平均每年增加0.000 018°.这十年的变化较为平稳.1996年～2009年这十三年的变化是比较迅速的，倾斜角增加了0.001 43°，平均每年增加0.000 143°，2009年～2011年这三年可能是由于时间较短，倾斜角仅增加了0.000 03°，预测往后的变化会越来越慢.

由此分析古塔的倾斜程度在1986年～1996年这十年间变化较慢；1996年～2009年这十三年间，变化较快，倾斜程度明显；2009年以后，基本趋于平稳，最后会达到一个上界：是古塔的最大倾斜角度.我们做一个短期预测，假设倾斜角随时间的函数是线性的，我们给出倾斜角的变化图.类似的分析扭曲角度的值，预测未来古塔的扭转程度.

图4 古塔发生扭转后顶点在地面投影Fig.4 Guta after the occurrence of a twisted vertex in the ground projection

图5 古塔发生扭转后每一层边在地面的投影Fig.5 The projection of each side of the ground in Guta after twisting

图6 倾斜角的变化图Fig.6 Variation of inclination angle

图7 扭转角的变化图
Fig.7 Variation of torsion angle

3 古塔模型的改进与评价

3.1 方案改进

(1)找到古塔开始建立的初始位置，计算出1986年的倾斜角；其余3个年份加入差额，即得出相对准确的倾斜角.(2)在求解古塔中心坐标时，可以利用C语言或Matlab 编程求出正八边形外接圆的圆心，中心的选取更为精确.(3)在扭曲程度上，可以有多种不同衡量的扭曲程度的量，可以进一步研究，准确到南北东西，方向各转动的度数，水平位移多少毫米.(4)弯曲上的测量，可以具体到量化，由于数据量过大，没有给出明确的数据解释.

3.2 对模型的评价

(1)优点：更为精确地给出了倾斜角，扭转角度的变化.(2)缺点：因为平衡位置是建塔的建造位置；但1986年已经开始变形，所以从1986年到2011年，1986年的倾斜角度存在偏差.

模型的评价：首先我们对数据进行了相关性分析、处理和量化，并且提出多种方法来解决并优化问题.同时利用Mathmatica软件进行计算，可信度高.论文中的图形与数据相结合更具有说服性.虽然此方法不是最优方法，但经过对模型拟合，数字量化已将误差减少到最少.