古塔的变形-数学建模论文-本科论文.doc

2013高教社杯全国大学生数学建模竞赛

承 诺 书

我们仔细阅读了《全国大学生数学建模竞赛章程》和《全国大学生数学建模竞赛参赛规则》（以下简称为“竞赛章程和参赛规则”，可从全国大学生数学建模竞赛网站下载）。
我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。
我们知道，抄袭别人的成果是违反竞赛章程和参赛规则的，如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。
我们郑重承诺，严格遵守竞赛章程和参赛规则，以保证竞赛的公正、公平性。如有违反竞赛章程和参赛规则的行为，我们将受到严肃处理。
我们授权全国大学生数学建模竞赛组委会，可将我们的论文以任何形式进行公开展示（包括进行网上公示，在书籍、期刊和其他媒体进行正式或非正式发表等）。

我们参赛选择的题号是（从A/B/C/D中选择一项填写）： C
我们的参赛报名号为（如果赛区设置报名号的话）：
所属学校（请填写完整的全名）： 山西建筑职业技术学院
参赛队员 (打印并签名) ：1. 王晋鹏
2. 李 帅
3. 徐 进
指导教师或指导教师组负责人 (打印并签名)： 原二保
（论文纸质版与电子版中的以上信息必须一致，只是电子版中无需签名。以上内容请仔细核对，提交后将不再允许做任何修改。如填写错误，论文可能被取消评奖资格。）

日期： 2013 年 9 月 日

赛区评阅编号（由赛区组委会评阅前进行编号）：

2013高教社杯全国大学生数学建模竞赛

编 号 专 用 页

赛区评阅编号（由赛区组委会评阅前进行编号）：

赛区评阅记录（可供赛区评阅时使用）：

评
阅
人

评
分

备
注

全国统一编号（由赛区组委会送交全国前编号）：

全国评阅编号（由全国组委会评阅前进行编号）:

C题 古塔的变形

摘要

文物保护是一件非常重要的事情，如何根据对古塔的测量数据了解其倾斜、弯曲、扭曲等情况并能够确定变形趋势是本文要解决的问题。针对问题一，本文将各层八个观测点看作同一平面上一个八边形的角点，利用公式，并编写C++程序计算其中心坐标，有较高的精确度。问题二要求确定古塔的倾斜、弯曲、扭曲等变形情况，结合测量专业知识，利用公式计算得到四年各层中心相对于下一层中心的偏心距（合位移）与倾斜角度以及四年古塔整体倾斜的合位移与倾斜角度：

1986年
1996年
2009年
2011年
合位移（m）
0.3779
0.3739
0.3778
0.3285
倾斜角（弧度）
0.0129
0.0121
0.0113
0.0103
利用公式计算反映各段的弯曲，利用；及 最终确定出各年古塔整体扭曲情况：

1986年
1996年
2009年
2011年
扭曲方向（度）
-34.3
-31.66
-36.52
-34.34
扭曲位移（m）
0.3752
0.3159
0.3035
0.2713
数据分析表明四年中底部1-3层及顶部11-13层扭曲现象严重，而中段4-10扭曲较平缓。应加强对底部及顶部的观测与维护。弯曲方面，四年相对变化不大，主要体现在各年不同层次的弯曲有一些变异，整体倾斜度不大且减小趋势。
关键词：重心（中心），倾斜、弯曲、扭曲 C++ MATLAB EXCEL

一 问题重述

由于长时间承受自重、气温、风力等各种作用，偶然还要受地震、飓风的影响，古塔会产生各种变形，诸如倾斜、弯曲、扭曲等。为保护古塔，文物部门需适时对古塔进行观测，了解各种变形量，以制定必要的保护措施。
某古塔已有上千年历史，是我国重点保护文物。管理部门委托测绘公司先后于1986年7月、1996年8月、2009年3月和2011年3月对该塔进行了4次观测。
请根据附件1提供的4次观测数据，讨论以下问题：
1.给出确定古塔各层中心位置的通用方法，并列表给出各次测量的古塔各层中心坐标。
2. 分析该塔倾斜、弯曲、扭曲等变形情况。
3. 分析该塔的变形趋势。
二 模型假设

1.古塔变形测量的监测点选取符合变形观测的原则，测量数据有效。
2．楼层平面质地均匀。
3. 在古塔各层选定的8（或7个）个监测点在同一平面内。

三 符号约定
1．：第个观测点的坐标。
2. △x：各层中心相对于下一层中心在x轴方向的位移量。
2. △y：各层中心相对于下一层中心在y轴方向的位移量。
3． d：各层中心相对于底层中心的合位移。
4．：各层中心相对于底层中心的扭曲方向。
5．：各层中心扭曲方向的平均值，反映建筑物主体的扭曲方向。
6．h：各层中心相对于第一层中心的高度值。
四 模型分析建立与求解
问题一：确定古塔各层中心
分析:附件给出了测绘公司先后于1986年7月、1996年8月、2009年3月和2011年3月对古塔进行4次观测结果，即古塔十三层及塔尖各监测测点的坐标值。
在建筑物的变形测量中每次都选取固定监测点观测。按附件提供的测测点的坐标值作图（图一），显现出正八边形形状，所以可以假设古塔为八角形，各层所设8个（7个）监测点为各层八角形的塔角。各层所设8个（7个）监测点坐标值的z坐标值虽有所不同，但差值最大仅为0.194m，可以忽略高度变化，故各层中心计算时假设8个（7个）塔角在同一平面内。
查阅测量相关资料知道：对于外轮廓为正八边形的塔形垂直的建筑物，需对8条棱进行观测，每层观测8个角点（在同一平面内），从8个角点构成八边形得到的几何中心就是该层中心。
在附件给出的4次观测数据中，还发现1986年和1996年数据中十三层却监测点5。根据古塔已上千年历史，可推断古塔1986年和1996年进行测量时，13层塔角5处于损毁状态，到2009年第三次测量前已经修复。
图1
13层塔角5在两次测量无数据，故对1986年和1996年13层数据计算得到中心坐标误差较大，这对后面计算会产生影响。

模型建立：在高等数学中已知物理学中求平面薄片重心的公式为：
，
当薄片质量分布均匀，即密度为一常数时，分式上下可约分，此时分母转化为平面图形的面积，公式简化为：
，
此时平面薄片的重心完全由平面图形形状决定，即为平面图形的中心。利用这组公式进行计算，计算难度决定于图形形状，但即使是正八边形，要确定边界表达式困难也很大，使得积分计算难度增加。
工程计算中，对于非圆形建（构）筑物的几何中心，通常为通过测定其外廓的对称于几何中心的棱角点的坐标并取其平均值的办法来求得其几何中心的坐标，即使用公式
，
本公式使用简单，一般能满足计算需要，但文献表明，由于是测量数据，本身存在误差，而减小误差也是操作中需要注意的，故选取公式[1]，，约定，进行求解。

模型计算：利用上述公式计算各层古塔中心，通过编制C++程序实现。（附件一），计算结果为表1
表1
各年中 心坐标
层数
1986年中心坐标
1996年中心坐标
x
y
z
x
y
z
1
566.6649
522.7092
1.7874
566.6651
522.7089
1.7830
2
566.7218
522.6108
7.3203
566.7227
522.6698
7.3146
3
566.7775
522.6332
12.7553
566.7453
522.6794
12.7508
4
566.8092
522.6141
14.7896
566.8234
522.6012
14.7858
5
566.8684
522.5713
21.7205
566.8711
522.5685
21.7160
6
566.9168
522.5399
26.2351
566.9202
522.5365
26.2295
7
566.9511
522.5231
29.8369
566.9552
522.5181
29.8323
8
566.9845
522.5066
33.3509
566.9686
522.5023
33.3454
9
567.0182
522.4900
36.8549
567.0229
522.4850
36.8483
10
567.0491
522.4759
40.1721
567.0542
522.4706
40.1676
11
567.1021
522.4356
44.4409
567.1078
522.4298
44.4354
12
567.1550
522.3954
48.7119
567.1391
522.3991
48.7074
13
567.1496
522.5437
52.8343
567.1558
522.5375
52.8300
塔尖
567.2473
522.2438
55.1233
567.2544
522.2367
55.1198

各年中
心坐标
层数
2009年中心坐标
2011年中心坐标
x
y
z
x
y
z
1
566.7378
522.6979
1.7645
566.7413
522.7004
1.7633
2
566.7760
522.6712
7.3090
566.7762
522.6710
7.2905
3
566.8094
522.6437
12.7323
566.8098
522.6441
12.7269
4
566.8368
522.6208
14.7734
566.8561
522.6536
14.7622
5
566.8655
522.5970
21.7094
566.8706
522.6019
21.7039
6
566.9540
522.5475
26.2110
566.9547
522.5467
26.2045
7
566.9871
522.5254
29.8246
567.0282
522.6162
29.8170
8
567.0252
522.4756
33.3399
567.0402
522.4923
33.3366
9
567.0915
522.4604
36.8438
567.0925
522.4595
36.8223
10
567.1463
522.4068
40.1611
567.1474
522.4056
40.1441
11
567.1890
522.3668
44.4326
567.2131
522.3930
44.4249
12
567.2313
522.3273
48.6998
567.2326
522.3260
48.6839
13
567.2801
522.2814
52.8184
566.7762
522.6710
52.8131
塔尖
567.3360
522.2148
55.0910
567.3375
522.2135
55.0870

问题二：古塔的倾斜、弯曲、扭曲
分析：问题一对每次测量得到了古塔各层的中心坐标，分析古塔的倾斜、弯曲、扭曲等情况时，首先考虑视各层中心连线近似于一条直线，所以首先通过中心连线与铅垂直线夹角计算初步了解古塔的变形情况。
初解：利用最小二乘拟合方法能够确定已知数据点的平面拟合曲线。对于已知空间数据点，若数据点近似分布在一条直线附近，亦可进行空间直线拟合，因空间直线可视为两个平面相交所成的直线，故分别对两个方程进行平面数据拟合从而确定的取值[2]。
运行如下matlab程序，得到各年各层中心对应的空间拟合直线。
function nh
x1=input('input x1:');
x=x1';
y1=input('input y1:');
y=y1';
z1=input('input z1:');
z=z1';
F=[z;1 1 1 1 1 1 1 1 1 1 1 1 1 1];
M=F*F';
N=F*x';
O=F*y';
A=(M\N)'
B=(M\O)'
A
B
zz=1:10:60;xx=A(1)*zz+A(2);yy=B(1)*zz+B(2);subplot(2,2,3);plot3(xx,yy,zz,'r');hold on,plot3(x,y,z,'+');title('2009'),grid on
直线方程为：
1986年：x=0.01023z+566.6475，y=-0.0058z+522.6969
1996年：x=0.0103z+566.6435，y=-0.0064z+522.7193
2009年：x=0.0115z+566.6699，y=-0.0088z+522.7555
2011年：x=0.0085z+566.7331，y=-0.0066z+522.79
接下来应用空间两条直线夹角公式
，

运用matlab计算各条直线与铅垂线夹角如表2。
表2

1986年
1996年
2009年
2011年
角度（度）
0.6704
0.6948
0.8296
0.6166

初步结论：对表1、表2数据初步分析可知塔身既有倾斜又有扭曲。幅度并不大。
进一步求解
首先计算各层中心与底层中心相比得到的x偏移与y偏移，得到表4：

表4
层数
1986年
1996年
2009年
2011年
△x
△y
△x
△y
△x
△y
△x
△y
2
0.0569
-0.0984
0.0576
-0.0391
0.0382
-0.0267
0.0349
-0.0294
3
0.1126
-0.0760
0.0802
-0.0295
0.0716
-0.0542
0.0685
-0.0563
4
0.1443
-0.0951
0.1583
-0.1077
0.0990
-0.0771
0.1148
-0.0468
5
0.2035
-0.1379
0.2060
-0.1404
0.1277
-0.1009
0.1293
-0.0985
6
0.2519
-0.1693
0.2551
-0.1724
0.2162
-0.1504
0.2134
-0.1537
7
0.2862
-0.1861
0.2901
-0.1908
0.2493
-0.1725
0.2869
-0.0842
8
0.3196
-0.2026
0.3035
-0.2066
0.2874
-0.2223
0.2989
-0.2081
9
0.3533
-0.2192
0.3578
-0.2239
0.3537
-0.2375
0.3512
-0.2409
10
0.3842
-0.2333
0.3891
-0.2383
0.4085
-0.2911
0.4061
-0.2948
11
0.4372
-0.2736
0.4427
-0.2791
0.4512
-0.3311
0.4718
-0.3074
12
0.4901
-0.3138
0.4740
-0.3098
0.4935
-0.3706
0.4913
-0.3744
13
0.4847
-0.1655
0.4907
-0.1714
0.5423
-0.4165
0.0349
-0.0294
塔尖
0.5824
-0.4655
0.5892
-0.4723
0.5982
-0.4831
0.5962
-0.4869
平均值
0.3159
-0.2028
0.3149
-0.1986
0.3028
-0.2257
0.2691
-0.1854
1．倾斜
层数　
1986
1996
2009
2011
合位移d
倾斜角α
合位移d
倾斜角α
合位移d
倾斜角α
合位移d
倾斜角α
1

2
0.1137
0.0205
0.0696
0.0126
0.0466
0.0084
0.0456
0.0083
3
0.1358
0.0124
0.0855
0.0078
0.0898
0.0082
0.0887
0.0081
4
0.1728
0.0133
0.1915
0.0147
0.1255
0.0097
0.1240
0.0096
5
0.2458
0.0123
0.2493
0.0125
0.1628
0.0082
0.1625
0.0082
6
0.3035
0.0124
0.3079
0.0126
0.2634
0.0108
0.2630
0.0108
7
0.3414
0.0122
0.3472
0.0124
0.3032
0.0108
0.2990
0.0107
8
0.3784
0.0120
0.3671
0.0116
0.3633
0.0115
0.3642
0.0115
9
0.4158
0.0119
0.4221
0.0120
0.4260
0.0122
0.4259
0.0122
10
0.4495
0.0117
0.4563
0.0119
0.5016
0.0131
0.5018
0.0131
11
0.5158
0.0121
0.5233
0.0123
0.5597
0.0131
0.5631
0.0132
12
0.5820
0.0124
0.5663
0.0121
0.6172
0.0132
0.6177
0.0132
13
0.5122
0.0100
0.5198
0.0102
0.6838
0.0134
0.0456
0.0009
塔尖
0.7455
0.0140
0.7551
0.0142
0.7689
0.0144
0.7698
0.0144
平均值
0.3779
0.0129
0.3739
0.0121
0.3778
0.0113
0.3285
0.0103
倾斜是建筑中心线或其墙、柱等，在不同高度的点对其相应底部点的偏移现象。高层建筑由于地基不均匀沉降而产生倾斜，可利用相互垂直的两个量来描述如图二所示，根据建筑物的设计，底层中心A点与顶层中心B点位于同一竖直直线上，当建筑物发生倾斜时，则B点相对A点移动了某一数值d，则建筑物的倾斜角为。计算得到其它年份各层中心与1986年各层中心的位移量及高度，求出倾斜角弧度数如表5，而各年份塔的合位移和倾斜度数如表6。
表5
表6

1986年
1996年
2009年
2011年
合位移（m）
0.3779
0.3739
0.3778
0.3285
倾斜角（弧度）
0.0129
0.0121
0.0113
0.0103

2.弯曲
观测古塔的弯曲变形就是考察古塔在外力作用下，其中心轴线由直线变成了曲线的情况。为此，计算各层中心相对于下一层中心的合位移以及各层中心相对于下一层中心的倾斜角（），通过比较合位移、倾斜角变化预测古塔的弯曲变形。见表7。
表7
1
9
8
6
年
层数
△xk
△y k
△z k
合位移d
角α（rad）
2
0.0569
-0.0984
5.5329
0.1137
1.5503
3
0.1126
-0.0760
5.4350
0.1358
1.5458
4
0.1443
-0.0951
2.0343
0.1728
1.4860
5
0.2035
-0.1379
6.9309
0.2458
1.5353
6
0.2519
-0.1693
4.5146
0.3035
1.5037
7
0.2862
-0.1861
3.6018
0.3414
1.4763
8
0.3196
-0.2026
3.5140
0.3784
1.4635
9
0.3533
-0.2192
3.5040
0.4158
1.4527
10
0.3842
-0.2333
3.3172
0.4495
1.4361
11
0.4372
-0.2736
4.2688
0.5158
1.4506
12
0.4901
-0.3138
4.2710
0.5820
1.4354
13
0.4847
-0.1655
4.1224
0.5122
1.4472
塔尖
0.5824
-0.4655
2.2890
0.7455
1.2559

1
9
9
6
年
层数
△xk
△y k
△z k
合位移d
角α(rad)
2
0.0576
-0.0391
5.5316
0.0696
1.5582
3
0.0802
-0.0295
5.4362
0.0855
1.5551
4
0.1583
-0.1077
2.0350
0.1915
1.477
5
0.2060
-0.1404
6.9302
0.2493
1.5348
6
0.2551
-0.1724
4.5135
0.3079
1.5027
7
0.2901
-0.1908
3.6028
0.3472
1.4747
8
0.3035
-0.2066
3.5131
0.3671
1.4667
9
0.3578
-0.2239
3.5029
0.4221
1.4509
10
0.3891
-0.2383
3.3193
0.4563
1.4342
11
0.4427
-0.2791
4.2678
0.5233
1.4488
12
0.4740
-0.3098
4.2720
0.5663
1.439
13
0.4907
-0.1714
4.1226
0.5198
1.4454
塔尖
0.5892
-0.4723
2.2898
0.7551
1.2522

2
0
0
9
年
层数
△xk
△y k
△z k
合位移d
角α
2
0.0382
-0.0267
5.5445
0.0466
1.5624
3
0.0716
-0.0542
5.4233
0.0898
1.5542
4
0.0990
-0.0771
2.0411
0.1255
1.5094
5
0.1277
-0.1009
6.9360
0.1628
1.5473
6
0.2162
-0.1504
4.5016
0.2634
1.5124
7
0.2493
-0.1725
3.6136
0.0466
1.4871
8
0.2874
-0.2223
3.5153
0.0898
1.4678
9
0.3537
-0.2375
3.5039
0.1255
1.4498
10
0.4085
-0.2911
3.3173
0.1628
1.4207
11
0.4512
-0.3311
4.2715
0.2634
1.4405
12
0.4935
-0.3706
4.2672
0.3032
1.4272
13
0.5423
-0.4165
4.1186
0.3633
1.4063
塔尖
0.5982
-0.4831
2.2726
0.4260
1.2445

2
0
1
1
年
层数
△xk
△y k
△z k
合位移d
角α
2
0.0349
-0.0294
5.5272
0.0456
1.5625
3
0.0685
-0.0563
5.4364
0.0887
1.5545
4
0.1148
-0.0468
2.0353
0.1240
1.5100
5
0.1293
-0.0985
6.9417
0.1625
1.5474
6
0.2134
-0.1537
4.5006
0.2630
1.5124
7
0.2869
-0.0842
3.6125
0.2990
1.4882
8
0.2989
-0.2081
3.5196
0.3642
1.4677
9
0.3512
-0.2409
3.4857
0.4259
1.4492
10
0.4061
-0.2948
3.3218
0.5018
1.4209
11
0.4718
-0.3074
4.2808
0.5631
1.4400
12
0.4913
-0.3744
4.2590
0.6177
1.4268
13
0.0349
-0.0294
4.1292
0.0456
1.5597
塔尖
0.5962
-0.4869
2.2739
0.7698
1.2444
图3
用matlab绘制弯曲示意图，根据各层中心相对于下一层中心的合位移相沿各层高度变化值的分布情况，各年度数据分别拟合为二次多项式曲线图3，表现出各年份古塔的弯曲变化过程。
3.扭曲

研究古塔的扭曲变形是考察古塔在外力作用下的扭转，可以通过古塔各层中心相对于底层中心的合位移d和扭曲方向表示。
计算各层中心相对于底层中心的合位移d和扭曲方向的公式：
，
取各层中心相对于底层中心扭曲方向的平均值，作为该建筑物主体扭曲方向。
古塔各层中心对底层中心的合位移d和扭曲方向计算结果见表8。

表8
1986年古塔扭曲数据
层次
△x'(m)
△y'(m)
合位移d(m)
合位移方向θ（度）
2
0.0569
-0.0984
0.1137
-59.9917
3
0.1126
-0.0760
0.1358
-34.0349
4
0.1443
-0.0951
0.1728
-33.4036
5
0.2035
-0.1379
0.2458
-34.1405
6
0.2519
-0.1693
0.3035
-33.9219
7
0.2862
-0.1861
0.3414
-33.0505
8
0.3196
-0.2026
0.3784
-32.3877
9
0.3533
-0.2192
0.4158
-31.8331
10
0.3842
-0.2333
0.4495
-31.2835
11
0.4372
-0.2736
0.5158
-32.0546
12
0.4901
-0.3138
0.5820
-32.6471
13
0.4847
-0.1655
0.5122
-18.8619
塔尖
0.5824
-0.4655
0.7455
-38.6536
平均值
0.3159
-0.2028
0.3779
-34.3280

1996年古塔扭曲数据
层次　
△x'(m)
△y'(m)
合位移d(m)
合位移方向θ（度）
1
0.0576
-0.0391
0.0696
-34.1868
2
0.0802
-0.0295
0.0855
-20.2053
3
0.1583
-0.1077
0.1915
-34.2469
4
0.2060
-0.1404
0.2493
-34.2939
5
0.2551
-0.1724
0.3079
-34.0686
6
0.2901
-0.1908
0.3472
-33.3500
7
0.3035
-0.2066
0.3671
-34.2615
8
0.3578
-0.2239
0.4221
-32.0533
9
0.3891
-0.2383
0.4563
-31.5009
10
0.4427
-0.2791
0.5233
-32.2457
11
0.4740
-0.3098
0.5663
-33.1849
12
0.4907
-0.1714
0.5198
-19.2639
13
0.5892
-0.4723
0.7551
-38.7298
14
0.3149
-0.1986
0.3739
-31.6609

2009年古塔扭曲数据
层次
△x'(m)
△y'(m)
合位移d(m)
合位移△'k(m)
合位移方向θ（度）
2
0.0382
-0.0267
0.0466
0.0010
-34.9695
3
0.0716
-0.0542
0.0898
0.0763
-37.1439
4
0.0990
-0.0771
0.1255
0.1221
-37.9303
5
0.1277
-0.1009
0.1628
0.1311
-38.3329
6
0.2162
-0.1504
0.2634
-0.2528
-34.8422
7
0.2493
-0.1725
0.3032
-0.3002
-34.6984
8
0.2874
-0.2223
0.3633
0.3630
-37.7406
9
0.3537
-0.2375
0.4260
-0.3366
-33.8975
10
0.4085
-0.2911
0.5016
-0.2981
-35.4919
11
0.4512
-0.3311
0.5597
0.0903
-36.2904
12
0.4935
-0.3706
0.6172
0.4408
-36.9239
13
0.5423
-0.4165
0.6838
0.6756
-37.5442
塔尖
0.5982
-0.4831
0.7689
0.2464
-38.9437
平均值
0.3028
-0.2257
0.3778
0.0738
-36.5192

2011年古塔扭曲数据
层次
△x'(m)
△y'(m)
合位移d(m)
合位移△'k(m)
合位移方向θ（度）
2
0.0349
-0.0294
0.0456
0.0403
-40.1313
3
0.0685
-0.0563
0.0887
-0.0147
-39.4367
4
0.1148
-0.0468
0.1240
-0.1215
-22.1902
5
0.1293
-0.0985
0.1625
0.1509
-37.3188
6
0.2134
-0.1537
0.2630
-0.0895
-35.7812
7
0.2869
-0.0842
0.2990
-0.2369
-16.3643
8
0.2989
-0.2081
0.3642
-0.3472
-34.8640
9
0.3512
-0.2409
0.4259
-0.4241
-34.4650
10
0.4061
-0.2948
0.5018
-0.0666
-35.9952
11
0.4718
-0.3074
0.5631
-0.0653
-33.1029
12
0.4913
-0.3744
0.6177
0.5758
-37.3284
13
0.0349
-