古塔变形的数学模型

1、古塔变形的数学模型摘要 古塔被誉为中国古代杰出的高层建筑物，历史悠久，值得并需要我们的保护。本文研究了关于古塔变形的问题，古塔的变形与塔身的中心紧密相关，具体分析了古塔倾斜、弯曲、扭曲的变形情况及趋势。对于问题一，建立中心位置模型，采用多边形组合形心的算法，求的结果是表6，7，8，9中的数据；对于问题二，我们将古塔的倾斜、弯曲和扭曲等变形情况，分别给予合理的数学描述。对于倾斜变形，我们定义了倾斜角a ，即塔尖与底层中心的水平距离与塔高的比值；对于弯曲变形，我们定义了弯曲率K ，即用中心点所拟合出的空间曲线的曲率来描述古塔各处弯曲率；对于扭曲变形，我们定义了相对扭曲度，利用坐标的旋转变换角度描述

2、古塔的扭曲变形情况。利用空间曲线拟合、坐标变换等方法以及MATLAB 程序，分别求出了三个变形刻画量的量化指标。对于问题三，我们考虑通过古塔的倾斜、弯曲及扭曲程度来分析古塔的变形趋势，由于数据较少，我们利用灰色预测法来分析这三种变形，利用matlab得到的数据来分析古塔的变形趋势关键词：关键词： 组合图形的形心 拟合 中心变形曲线 灰色预测 一 问题重述古塔由于长时间承受自重、气温、风力等各种作用，偶然还要受地震、飓风的影响会产生各种变形，诸如倾斜、弯曲、扭曲等。为保护古塔，文物部门需适时对古塔进行观测，了解各种变形量，以制定必要的保护措施。现有一古塔已上千年历史，是我国重点保护文物。管理部门

3、委托测绘公司先后于1986年7月、1996年8月、2009年3月和2011年3月对该塔进行了4次观测。请根据题目提供的4次观测数据，讨论以下问题：1. 给出确定古塔各层中心位置的通用方法，并列表给出各次测量的古塔各层中心坐标。2. 分析该塔倾斜、弯曲、扭曲等变形情况。3. 分析该塔的变形趋势。二、 模型假设1，参照中国古塔的塔身样式假设塔身水平面是正八角形；2，假设古塔是均质物体，形心就是中心； 3，以地平面为X轴，建立空间直角坐标系X，Y，Z；4，假设塔尖面积极小，近视为一个点；5，假设模型4中古塔所在地不存在地震等这种偶然事件；6，假设古塔所在的地平面是水平面。三、符号说明Ci-第i个三角

4、形的形心（i=1,2,6）；Si-第i个三角形的面积（单位：m2)（i=1,2,6）；Ai-第i个三角形的形心坐标（单位：m）；Wt-第t年塔身的倾斜角(t=1986,1996,2009,2011)；Dt-第t年塔身在地面上的水平投影(t=1986,1996,2009,2011)；Ht-第t年塔上部与下部的相对位移(t=1986,1996,2009,2011)；Lt-第t年塔的拟合直线(t=1986,1996,2009,2011)；R-塔的倾斜率；mt-第t年中心直线x轴的法向量(t=1986,1996,2009,2011)；nt-第t年中心直线y轴的法向量(t=1986,1996,2009,

5、2011)；Jt-第t年的中心直线方向向量(t=1986,1996,2009,2011)；K-塔的曲率；四、模型准备（1） 给数据用图像形式表现出来题目所给的4次观测数据1986年、1996年塔13层的第五个点数据残缺。参照其10层、11层、12层的数据知道第五个点的数值在第四个点和第六个点的数值之间，另外分析后发现这两组数据对题目结果的影响不大。因此结合实际和理论，采取取相邻年份该层的第四个点和第六个点的平均值将两组数据补齐。（2） 用MATLAB软件将题中所给数据及补充的数据画图表现出来。五、模型的建立与求解问题一塔的变形趋势与塔身的中心有关，因此本题的主要目标是计算出塔中心的变化规律。问

6、题1属于计算多边形组合形心的数学问题。因为1986年受自然、社会影响度较小，最接近古塔的原始图样，所以图形较为有参考价值，用MATLAB画出1986年塔身的大致图样为八边形，通过假设得出正多边形的形心就是中心，因此可运用组合图形形心公式计算出古塔各层中心坐标。问题一求古塔各层中心坐标，通过假设知道实质就是计算正多边形的中心，运用组合图形求形心的方法。将每层正八边形的塔面都分为六个三角形（图1-2）S1、S2、S3、S4、S5、S6，每年各层六个三角形的形心：C1(，)、C2(,)、C3(,)、C4(,)、C5(,)、C6(,)， 代入测量数据，用excel得出：19861986年各层分面积的形

7、心 塔层 坐标(m)A1A2A3A4A5A61565.6707，526.9707 564.3160，524.7823 565.1133，523.6693567.8320，525.6177 570.8977，523.7220570.3780，520.40572565.7107，526.7563564.4383，524.6417565.2177，523.5770567.8177，525.4710570.7933，523.6700570.3203，520.47733565.7503，526.5460564.5587，524.5040565.3203，523.4870567.8037，525.3273

8、567.6070，524.0603567.1797，520.98904565.7813，526.3787564.6540，524.3940565.4017，523.4147567.7920，525.2127570.6087，523.5783570.2180，520.60375565.8150，526.1987564.7567，524.2760565.4893，523.3370567.7800，525.0893570.5210，523.5347570.1693，520.66406565.8483，526.0230564.8583，524.1610565.5757，523.2620567.2693

9、，524.9697570.4320，523.4920570.1233，520.72177565.9323，525.8013565.0030，524.0360565.6087，522.1230568.3457，522.8950568.9073，520.9863567.9333，519.21108566.0143，525.5857565.1437，523.9140565.7230，522.1233568.2910，522.8597568.8363，521.0700567.9370，519.40309566.0963，525.3700565.2862，523.7907565.8408，522.120

10、7568.2398，522.8213568.7700，521.1480567.9455，519.771810566.1781，525.1621565.4221，523.6731565.9518，522.1205568.1901，522.7871568.7034，521.2268567.9501，519.771811566.2173，525.1621565.3411，522.7775565.9963，522.1686568.1775，522.8075568.6384，521.1750567.9245，519.905512566.4234，524.7179565.7663，523.4303566.

11、2073，522.0879568.1537，522.6799568.5723，521.3238567.8977，520.039513566.5415，524.5037565.9325，523.3134566.3304，522.0724568.1358，522.6281568.5088，521.3708568.5314，521.348219961996年各层分面积的形心 塔层 坐标(m)A1A2A3A4A5A61565.6709，526.9704564.3162，524.7820564.8909，522.3034568.4069，523.1387568.9919，520.6390567.6595

12、，518.45702565.7116，526.7554564.4392，524.6499564.9537，522.2315568.3784，523.0668567.0461，520.8848567.7212，518.58243565.7519，526.5444564.5603，524.5024565.1606，522.2207568.4056，523.0441 568.9876，520.7551567.7816，518.20574565.7835，526.3765564.6562，524.3918565.2669，522.1878568.4048，523.0065568.9855，520.80

13、05567.8299，518.80355565.8178，526.1958564.7594，524.2731567.2741，521.0509566.4548，524.1346568.5384，524.6946570.4308，523.52496563.6743，521.6294565.4214，520.6144567.3460，521.1294566.5451，523.9446568.4659，524.4789570.2377，523.47457563.8114，521.6232565.5112，520.6525565.6125，522.1191568.3495，522.8911568.91

14、11，520.9825567.9371，519.20718566.0185，525.5814565.1478，523.9097565.7272，522.1190568.2953，522.8554568.8405，521.0657567.9412，519.39879564.2483，521.6049565.7973，520.7773567.4596，521.2459566.6890，523.6426568.2713，522.7973567.9455，519.589310564.4428，521.5880565.9130，520.7977567.4925，521.2496566.7358，523.

15、5298568.327，524.0022569.2540，523.200911564.6655，521.5956566.0338，520.8531567.5067，521.2830566.6366，523.4134568.3226，523.8619569.7207，523.117012564.8876，521.6031566.1544，520.9094567.5207，521.3166566.9068，523.2978568.2923，523.7225569.5873，523.033713565.1091，521.6032566.2779，520.9556567.5416，521.341356

16、6.9817，523.1791568.2704，523.5812569.4667，522.946520092009年各层分面积的形心 塔层 坐标(m)A1A2A3A4A5A61562.4995，521.6767564.6728，520.3384567.1635，520.9321566.2995，524.4607568.7379，525.0660570.9705，523.73002562.7235，521.6673564.8138，520.3847567.1945，520.9648569.5556，521.5574569.5059，523.8796570.8206，523.67273562.94

17、24，521.6587564.9515，520.4308567.2241，520.9976566.3847，524.2815568.6415，524.8619570.6734，523.61734563.1175，521.6513565.0616，520.4670567.2484，521.0232566.4188，524.2100568.5920，524.7807570.5565，523.57275563.3050，521.6439565.1797，520.5062567.2741，521.0509566.4548，524.1346568.5384，524.6946570.4308，523.52

18、496563.6743，521.6294656.4214，520.6144567.3460，521.1294566.5451，523.9446568.4659，524.4789570.2377，523.43457563.8114，521.6232565.5112，520.6525567.3735，521.1561566.5799，523.8915568.4390，524.3948570.1656，523.39808564.0303，521.6142565.6546，520.7151567.4165，521.2012566.6344，523.7572568.3961，524.2642570.05

19、04，523.34339564.2483，521.6049565.7973，520.7773567.4596，521.2459566.6890，523.6426568.3533，524.1329569.9357，523.287610564.4428，521.5880565.9130，520.7977567.4925，521.2496566.7658，523.5298568.3527，524.0022569.8540，523.200711564.6655，521.5955566.0338，520.8535567.5067，521.2830566.8366，523.4134568.3226，523

20、.8619569.7207，523.117012564.8876，521.6031566.1544，520.9094567.5207，520.3166566.9068，523.2978568.2923，523.7225569.5873，523.033713565.1091，521.6032566.2779，520.9556567.5416，521.3413566.9817，523.1791568.2704，523.5812569.4667，522.946520112011年各层分面积的形心 塔层 坐标(m)A1A2A3A4A5A61562.49957，521.6766564.6856，5204

21、421564.6730，520.3383567.1637，520.9320567.1413，520.9436568.4800，523.13632562.7237，521.6671564.8140，520.3844567.1947，520.9645566.3432，524.3695568.7042，524.9621570.8208，523.67253562.9428，521.6591564.9518，520.4312567.2245，520.9979566.3851，524.2819568.6419，524.8622570.6737，523.61764563.1179，521.6509565.0

22、620，520.4666567.2488，521.0227566.4192，524.2096568.5924，524.7802570.5569，523.57225563.3056，521.6433565.1803，520.5056567.2747，521.0503566.4554，524.1339568.5391，524.6939570.4314，523.52436563.6750，521.6286565.4221，520.6136567.3467，521.1287566.5458，523.9439568.4666，524.4782570.2384，523.43387563.8122，521.

23、6224565.5120，520.6517567.3743，521.1553566.5807，523.8707568.4398，524.3940570.1664，523.39728564.0311，521.6133565.6554，520.7142567.4173，521.2003566.6352，523.7563568.3969，524.2633570.0512，523.34249564.2493，521.6040564.7983，520.7764567.4607，521.3884566.6899，523.1730568.3543，523.6633569.9366，522.818010564

24、.4439，521.5869565.9141，520.7966567.4936，521.2485566.7669，523.5287568.3538，524.0011569.8551，523.199611564.666，521.5944566.0349，520.8524567.5078，521.2819566.8377，523.4123568.3237，523.8607569.7218，523.115912564.8889，521.6018566.1557，520.9081567.5220，521.3153566.9081，523.2965568.2936，523.7212569.5886，52

25、3.032413565.1104，521.6018566.2792，520.9542567.5429，521.3398566.9830，523.1777568.2717，523.5797569.4680，522.9450每层六个三角形的面积S1、S2、S3、S4、S5、S6，代入测量数据，运行MATLAB软件得出：13层塔中各层划分后的三角形面积 塔层 面积(m2)S1S2S3S4S5S616.106714.639720.606821.002514.70746.256325.730713.640219.066219.390213.58645.839735.370612.696217.61211

26、7.869412.535.446145.099511.970216.501616.705911.72115.140554.813711.214115.350915.505810.87774.822864.538710.502914.256414.389210.14.523974.059.267212.639112.73328.92594.002888.13911.158211.21717.85133.528293.17897.08179.77769.80146.84873.0811102.80146.16058.55718.54015.96732.6792112.47635.34497.531

27、27.4125.23612.3588122.16844.58736.56476.36134.558372.0604131.89313.90795.69685.42383.93922.6792多边形的形心（xc,yc,zc)公式： 将表1，2，3，4，5中的数据代入形心公式中，excel计算的结果如下1986 1986年塔的中心坐标 塔层 坐标(m)XcYcZc1567.1157524.37331.78742567.1355524.25837.32033566.38524.267212.75534567.6101523.857317.07835567.1858523.987921.7205656

28、7.2786523.784626.23517567.0062522.346129.83698566.9845522.506633.35099567.0212522.486836.854910567.0562522.470740.172111567.0363522.34644.440912567.1611522.389348.711913567.312522.402452.829319961996年塔的中心坐标 塔层 坐标(m)XcYcZc1566.6531522.76061.7832566.3558522.69047.31463566.7792522.631712.75084566.68885

29、23.089217.07515567.0364523.466621.7166567.04522.557926.22957566.9986522.132829.83238566.9059522.101133.34549566.8944521.920636.848310567.071522.304540.167611567.1643522.340144.435412567.2123522.299748.707413567.5044522.404552.825320092009年塔的中心坐标 塔层 坐标(m)XcYcZc 1566.7412522.70051.76452567.7195521.777

30、97.3093566.8095522.643712.73234566.8368522.623417.06985566.8656522.59721.70946566.9541522.547226.2117966.9871522.530529.82468567.0394522.493233.42949567.0915522.460436.853810567.1463522.589540.636811567.1897522.366844.432612567.2313522.069348.699813567.2799522.282252.818420112011年塔的中心坐标 塔层 坐标(m)XcYc

31、Zc 1565.8577520.91661.76332566.7785522.67067.29053566.8098522.644112.72694566.8374522.619717.0525566.8662522.596421.70396566.9548522.546826.20457566.9879522.524629.8178567.0402522.492333.33669567.0925522.260236.822310567.1474522.405740.269111567.1901522.365644.42512567.2326522.32648.683913567.241522

32、.32652.8131找出古塔各层中心位置的通用方法：将每层平面多边形划分为6个相互独立、互不相交的三角形，并求出每层的6个三角形的形心和面积，运用组合图形的形心坐标公式得出古塔各层中心坐标（xc，yc,zc），即中心位置。问题二类别塔身第一层的中心塔尖 坐标(m) 年份 xoyoZOXYZ1986567.1157524.37331.7874567.2473522.24455.12331996566.6531522.76061.783567.2544522.236755.11982009566.7412522.70051.7645567.336522.214855.0912011565.857

33、7520.91661.7633567.3375522.213555.087各年塔尖平面坐标(x,y）,各年塔身第一层的中心平面坐标(xo,yo)；塔身在地面上的水平投影D: 塔上部与下部的相对位移H: 在直角三角形OPQ中存在： tanWt即为倾斜度；倾斜角W: 倾斜率： r=arctanWt 代入表10的数据运算得出，W1986=1.5308,W1996=1.5558,W2009=1.5564,W2011=1.5339；r1986=2.1334,r1996=0.7975,r2009=0.7679,r2011=1.9677；即1986年的塔身倾斜了1.5308，1996年的塔身倾斜了1.555

34、8，2009年的塔身倾斜了1.5564，2011年的塔身倾斜了1.5339。5.3-2 问题2-1.2 相对倾斜模型对各年份各层的中心点数据（表6、7、8、9）进行空间直线拟合，利用MATLAB软件拟合各年份塔中心线，结果如下： L1986 x=0.0028z+567.0229 y=-0.04739z+524.57 L1996 x=0.015z+566.5093 y=-0.0474y+524.5711 L2009 x=0.0057z+566.9146 y=-0.0041y+522.5458 L2011 x=0.0181z+566.4005 y=0.0048y+522.1988由于两平面的交线与

35、这两平面的法线向量都垂直，所以可取 J1986= = J1996= = J2009= = J2011= 1 = 得出中心直线的方向向量：L1986：（0.0281，-0.0473,1）、L1996：(,-0.0041,1)L2009：(0.0057,-0.0041,1)、L2011：(0.181,-0.0048,1)两向量的夹角的余弦公式： 相对倾斜角： 代入1986年和1996年、1996年和2009年、2009年和2011年塔中心线的方向向量的数据求得W1=0.049,W2=0,W3=0.0245；即1996年相对于1986年倾斜了0.049，2009年相对于1996年没发生倾斜，2011

36、年相当于2009年倾斜了0.0245。5.3-2 问题2-1 弯曲模型运用表6、7、8中竖坐标的数据在MATLAB中进行三次多项式的拟合，得到：y1986=0.015x3-0.3853x2y1996=0.0152x3-0.3855x2y2009=0.0152x3-0.385x2y2011=0.0149x3-0.3799x2对拟合后的直线进行一阶求导、二阶求导，1986=0.045x219861996=0.0456x219962009=0.0456x22009 2011=0.0447x22011曲率计算公式为： MATLAB软件运行后的结果：1986年每层曲率：0.0028、0.0033、0.0

37、038、0.0041、0.0041、0.0036、0.0025、0.0010、0.0008、0.0023、0.0035、0.0041、0.00421996年每层曲率：0.0028、0.0033、0.0038、0.0041、0.0040、0.0034、0.0023、0.0008、0.0009、0.0024、0.0035、0.0040、0.00412009年每层曲率：0.0028、0.0033、0.0038、0.0040、0.0040、0.0034、0.0023、0.0007、0.0009、0.0024、0.0035、0.0040、0.00402011年每层曲率：0.0005、0.0028、0.0

38、030、0.0021、0.0013、0.0007、0.0004、0.0003、0.0002、0.0001、0.0001、0.0000、0.0000对于古塔的扭曲，我们以古塔每层的中心点作为研究对象，研究其横纵坐标的变化，通过比较相邻两年的中心坐标，算出对应的每层中心点偏移的平面距离作为扭曲距离，设每层扭曲距离的权重为wi，每层中心点的扭曲距离为mi，则建立古塔总扭曲距离型为：总扭曲距离：直接用计算器算得扭曲距离如下表所示： 1986年-1996年1996年-2009年2009年-2011年m11.67770.10661.9907m21.75111.64081.2971m31.68350.032

39、60.0005m41.19950.48870.0037m50.54230.88620.0008m61.24970.08660.0008m71.01511.08690.0059m81.00670.95110.0012m90.58020.57470.2002m100.16690.29480.1838m110.12810.03690.0013m120.10320.23120.2567m13 0.21890.25570.0586 结合总扭曲距离公式和matlab软件算得，相应年份的总扭曲距离为1996年总扭曲距离11.9399, 2009年总扭曲距离8.4095，2011年总扭曲距离为3.8508问题

40、三本题是分析古塔的变形情况。本文中，我们认为建筑物变形由建筑物的倾斜、弯曲、扭曲等因素共同造成。由于附录只给出了四次统计的数据，而我们的目标是分析古塔未来多年的变化趋势，因此我们采用信息不完全、不充分的预测系统灰色预测对古塔未来的变形趋势进行预测。我们建立灰色预测模型GM（2，1）模型分别从古塔的倾斜，弯曲，扭曲三个方面来研究古塔变形的趋势。 由第二问所得到的的数据结果可知三种变形的原始数据序列 带入灰色预测模型，用matlab进行求解， 倾斜角预测函数x = 1.55725*t + 0.0257409*exp(-9.68306*t) + 1.28226epsilon = 0 0.0024 -

41、0.0014 -0.0009delta = 0 0.0016 0.0009 0.0005同理曲率的预测函数及残差、相对误差值楼层10.0028*t - 4.33681e-19*exp(t) + 0.0028epsilon = 0.0000 0.0000 0.0000 -0.0023delta =0.0000 0.0000 0.0000 4.6000楼层2x =0.0033*t - 8.67362e-19*exp(t) + 0.0033epsilon = 1.0e-03 * 0.0000 0.0000 0.0000 -0.5000delta = 0.0000 0.0000 0.0000 0.17

42、86楼层3x =0.0038*t + 8.67362e-19*exp(t) + 0.0038epsilon = 1.0e-03 * -0.0000 -0.0000 -0.0000 -0.8000delta =0.0000 0.0000 0.0000 0.2667楼层4x =0.0041*t - 8.67362e-19*exp(t) + 0.0041epsilon = 0.0000 0.0000 0.0000 -0.0020delta =0.0000 0.0000 0.0000 0.9524楼层5x = 0.00405094*t + 0.0000499822*exp(0.981481*t) +

43、0.00405002epsilon = -0.0000 -0.0001 -0.0003 -0.0033delta =0.0000 0.0336 0.0684 2.5729楼层6x =0.00350385*t + 0.0000998521*exp(0.962963*t) + 0.00350015epsilon = 0.0000 -0.0003 -0.0005 -0.0039delta =0.0000 0.0781 0.1551 5.5906楼层7x = 0.00240556*t + 0.0000996914*exp(0.947368*t) + 0.00240031epsilon = -0.000

44、0 -0.0003 -0.0005 -0.0031delta =0.0000 0.1143 0.2224 7.6310楼层8x =0.00106667*t - 0.000133333*exp(0.5*t) + 0.00113333epsilon = 1.0e-03 * 0 -0.1802 -0.2241 -0.5315delta = 0 0.2252 0.3201 1.7718楼层9x =0.000853333*t - 0.0000497778*exp(1.07143*t) + 0.000849778epsilon = 1.0e-03 * 0.0000 0.1422 0.3256 0.1611

45、delta =0.0000 0.1580 0.3618 0.8056楼层10x = 0.00235106*t - 0.0000499774*exp(1.02174*t) + 0.00234998epsilon = -0.0000 0.0001 0.0003 -0.0016delta =0.0000 0.0574 0.1232 15.6529楼层11x = 0.0035*t + 0.0035epsilon = 0 0 -0.0000 -0.0034delta = 0 0 0.0000 34.0000楼层12x = 0.00405063*t + 0.000049992*exp(0.9875*t)

46、+ 0.00405001epsilon = -0.0000 -0.0001 -0.0003 -0.0047delta =0.0000 0.0337 0.0692 Inf楼层13x = 0.00415451*t + 0.0000472603*exp(0.962523*t) + 0.00415274epsilon = 0 -0.0001 -0.0004 -0.0047delta = 0 0.0319 0.0887 Inf扭曲距离的预测函数及残差、相对误差楼层1x =0.963311*t + 0.389552*exp(1.83387*t) + 1.28815epsilon = -0.0000 -2.

47、9050 -11.7910delta =0.0000 27.2515 5.9230楼层二x = 1.80323*t - 0.0767589*exp(0.679081*t) + 1.82786epsilon = 0 -0.0878 -0.3590delta = 0 0.0535 0.2768楼层三x = 1.71689 - 0.0333857*exp(-50.4299*t) - 0.000136527*tepsilon = 1.0e-03 * 0.0000 -0.6491 0.6365delta =0.0000 0.0199 1.2731楼层四x = 6.00555 - 4.80605*exp(

48、-0.465567*t) - 1.03804*tepsilon = 0.0000 -0.2622 -0.0813delta =0.0000 0.5365 21.9736楼层5x = 0.638507*t - 0.0692928*exp(1.38841*t) + 0.611593epsilon = 0 0.4562 0.1979delta = 0 0.5147 247.4058楼层6x = 1.34971 - 0.100011*exp(-12.5559*t) - 0.00603342*tepsilon = 0 -0.0074 0.0068delta = 0 0.0852 8.5413楼层7x = 1.01957*t - 0.0041934*exp(1.06642*t) + 1.01929epsilon = 0.0000 0.0753 -0.9905delta = 0.0000 0.0693 167.8756楼层