建模c题论文.docx
《建模c题论文.docx》由会员分享，可在线阅读，更多相关《建模c题论文.docx（20页珍藏版）》请在冰豆网上搜索。
建模c题论文
古塔的变形模型
摘要
古塔具有一定的历史文化价值，为更好地加强对国家重要文化遗产的保护，提出具体的保护措施和实施方案，因此分析古建筑的变形情况和变形趋势尤为重要。
本文通过对古塔的倾斜、弯曲、扭曲等各种形变进行分析研究，并建立数学模型，通过某古塔的实际测量数据计算验证，得出该古塔的各种形变数据，并对各种变形关于时间作出回归拟合，推测出未来几年的变化趋势。
对于问题一，通过查阅相关资料，结合给出数据，假设古塔的水平截面为近似正八边形，因为正多边形的中心与它的外接圆的圆心重合，已知八个点和七个点误差不大，所以运用拟合圆曲线的方法求出了中心点坐标。
对于问题二，本文全方位的分析了塔身的变形情况，主要包括了倾斜，弯曲，扭曲这三方面的变形。
以塔顶与塔底的倾斜角为指标来研究塔的倾斜情况，以曲率为指标来研究塔的弯曲情况，以偏转角度为指标来研究塔的扭曲情况，并借助MATLAB数学软件求解，综合分析出了塔身的变形情况，结果见文中表格。
对于问题三，本文主要从横向的方向进行研究，分别求时间（年份）与塔身的倾斜角度、塔身的每个层面上拐点的曲度、塔身的相对扭曲角度的趋势。
运用回归分析求解出塔身的倾斜与时间（年份）间的函数关系式，再用各个相对指标整体的分析出了该塔的整体变化趋势。
所以可根据此次预测，求出任意一年塔的整体变化趋势，便于日后对古塔进行检测和修复。
针对论文的实际情况，对论文的优缺点做了评价，在最后还给出了其他的改进方向和对论文中模型的推广思路。
关键词:
古塔变形；中心点坐标；回归分析；MATLAB；拟合；
一、问题重述
古塔具有一定的历史文化价值，但由于长时间承受自重、气温、风力等各种作用，偶然还会受地震、飓风的影响，古塔会产生各种变形，比如倾斜、弯曲、扭曲等。
为保护古塔，文物部门需适时对古塔进行观测，了解各种变形量，以制定必要的保护措施。
某古塔已有上千年历史，是我国重点保护文物。
管理部门委托测绘公司先后于1986年7月、1996年8月、2009年3月和2011年3月对该塔进行了4次观测。
依据附件1提供的4次观测数据，讨论以下三个问题：
问题一：
给出确定古塔各层中心位置的通用方法，并列表给出各次测量的古塔各层中心坐标；
问题二：
分析该塔倾斜、弯曲、扭曲等变形情况；
问题三：
分析该塔的变形趋势。
二、问题假设
1.假设四年的测量数据来自同一测绘地点；
2.假设给出的测量数据真实可靠，在允许的误差范围内；
3.假设古塔的倾斜、弯曲、扭曲三种变形相互之间影响较小；
4.假设该古塔的水平截面近似为正八边形，各层的测量点为八边形的顶点；
5.假设该古塔初建成时，塔身垂直于水平面；
6.假设缺失的数据均是随机产生的。
三、符号说明
：
古塔第
层的中心及其坐标；
：
古塔的每层第j个测量点及其坐标；
：
古塔的倾斜角；
：
古塔的扭曲角；
：
曲率；
空间曲线的挠率。
四、问题的分析
4.1问题的重要性分析
很多古塔都是重要的文物和景观，且具有一定的历史文化价值，但由于年代久远，古塔经受长期的自然损坏，所以对其进行有效的保护是一项十分必要而艰巨的任务。
基于古塔的重要性，需对其采取一定措施进行保护，进而延长古塔的寿命，这就要求对其进行观测，掌握其变形程度，这也是对古塔进行有效保护的前提和基础。
4.2问题的分析
对问题一的分析：
本题是关于几何图形中心的计算问题。
通过查阅资料了解到古塔通常为正多边形，结合附件1中给出的古塔各层观测点的个数，可把题目中所给古塔各层水平截面看做正八边形，则其中心点坐标即为重心坐标，但由于1986、1996两年测得的数据不完整，所以用这种方法求出的
坐标会有很大误差。
因为正多边形的中心与它的外接圆的圆心重合，已知八个点和七个点误差不大，所以采用MATLAB拟合圆曲线的方法找出中心点的坐标。
通过观察数据易知各点
坐标变化不大，比较稳定,即便数据不全也不影响数据的准确性，因此可采用取平均值的方法求出。
对问题二的分析：
本题分析该古塔的倾斜、弯曲、扭曲等变形情况。
在计算出中心坐标的基础上求出各中心的拟合直线，把这条直线与竖直方向的夹角作为古塔的倾斜角
，由次可计算出古塔的倾斜度；假设古塔各层中心和塔尖在同一平面内，其连线为平面内的一条曲线，用该曲线的弯曲程度描述古塔的弯曲程度；用各层相对于第1层的水平旋转角，来描述古塔的扭曲情况，找出其旋转角与高度的变化关系，用古塔的最大旋转角度作为反映古塔扭曲情况的扭曲角。
对问题三的分析：
本题求塔的变形趋势。
由于古塔的变形是长时间的累积作用，在问题二的基础上把各种形变与时间（年份）作回归拟合，并进行验证，在此基础上预测未来该古塔的变化趋势。
五、模型的建立与求解
5.1问题一的模型建立与求解
通过观察附件1提供的原始数据，利用MATLAB编程，做出三维图像，得到观察图像的俯视图（见图5.1.1）和立体图（见图5.1.2），易直观的发现每一层近似为正八边形。
图5.1.1古塔的俯视图
图5.1.2古塔的立体图
因为正多边形的中心与它的外接圆的圆心重合，已知八个点和七个点拟合出来的圆差别不会很大，所以采用MATLAB用函数nlinfit拟合圆曲线的方法找出圆心的坐标
此圆心坐标即为该层中心
坐标。
通过观察数据易知各点
坐标变化不大，比较稳定,所以
坐标用每层
坐标值的平均值作为中心的
坐标。
因此，可求出各层中心的三个分量
。
拟合后得到的古塔各层中心坐标如下表所示：
表5.1.11986年各层古塔中心坐标
层数
坐标
坐标
坐标
第一层
566.665
522.709
1.787375
第二层
566.7225
522.6713
7.32025
第三层
566.7786
522.6346
12.75525
第四层
566.8226
522.6056
17.07825
第五层
566.8698
522.5746
21.7205
第六层
566.9185
522.5443
26.23513
第七层
566.9521
522.5271
29.83688
第八层
566.9846
522.5102
33.35088
第九层
567.0176
522.4933
36.85488
第十层
567.0476
522.4789
40.17213
第十一层
567.1013
522.4385
44.44088
第十二层
567.1552
522.3981
48.71188
第十三层
567.2081
522.3578
52.83463
塔尖
567.2473
522.2438
55.12325
表5.1.21996年各层古塔中心坐标
层数
坐标
坐标
坐标
第一层
566.6652
522.7087
1.783
第二层
566.7234
522.6703
7.314625
第三层
566.7803
522.633
12.75075
第四层
566.8248
522.6034
17.07513
第五层
566.8726
522.5717
21.716
第六层
566.9219
522.5409
26.2295
第七层
566.9559
522.5232
29.83225
第八层
566.9888
522.506
33.34538
第九层
567.0223
522.4883
36.84825
第十层
567.0527
522.4737
40.16763
第十一层
567.107
522.4327
44.43538
第十二层
567.1612
522.3919
48.70738
第十三层
567.2142
522.3516
52.77438
塔尖
567.2544
522.2367
55.11975
表5.1.32009年各层古塔中心坐标
层数
坐标
坐标
坐标
第一层
566.7447
522.7004
1.7645
第二层
566.779
522.672
7.309
第三层
566.8118
522.6455
12.73225
第四层
566.8388
522.6232
17.06975
第五层
566.8671
522.6003
21.70938
第六层
566.9561
522.5512
26.211
第七层
566.9892
522.5289
29.82463
第八层
567.0417
522.4967
33.33988
第九层
567.094
522.464
36.84375
第十层
567.1484
522.41
40.16113
第十一层
567.1908
522.3695
44.43263
第十二层
567.2329
522.3298
48.69975
第十三层
567.2814
522.2844
52.81838
塔尖
567.336
522.2148
55.091
表5.1.42011年各层古塔中心坐标
层数
坐标
坐标
坐标
第一层
566.7448
522.7003
1.76325
第二层
566.7792
522.6718
7.2905
第三层
566.8122
522.6458
12.72688
第四层
566.8392
522.6228
17.052
第五层
566.8677
522.5997
21.70388
第六层
566.9567
522.5504
26.2045
第七层
566.99
522.5281
29.817
第八层
567.0425
522.4959
33.33663
第九层
567.095
522.463
36.82225
第十层
567.1495
522.4089
40.14413
第十一层
567.1919
522.3684
44.42488
第十二层
567.2342
522.3285
48.68388
第十三层
567.2827
522.2829
52.81313
塔尖
567.3375
522.2135
55.087
5.2问题二的模型建立与求解
5.2.1古塔的倾斜模型
利用问题1模型的解，用MATLAB画出塔在这四次观测中由各层中心拟合的直线的图像如下：
图5.2.1
从图中可以看出四条直线基本重合，这条拟合出来的直线与z轴的夹角即为该古塔的倾斜角
。
利用空间解析几何中的夹角公式可以求出古塔的倾斜角
，
假设两直线交点为
交点到地面的向量为
拟合出的直线与交点和地面间的向量为
则
然后带入求出各年份古塔的倾斜角
，如表5.2.1所示
表5.2.1各年份古塔的倾斜度
时间
1986年
1996年
2009年
2011年
倾斜度（
）
0.7517
0.7614
0.8531
0.8552
由上表可知，古塔的倾斜度有变大的趋势，且还可以看出，1986到1996年这十年间倾斜度变化较小。
之后，倾斜度变化很大，说明近年来古塔的变形较为严重。
5.2.2古塔的弯曲模型
基于一般情况假设各层中心点坐标
应该是一条连续曲线，本文研究塔身的弯曲情况就应该研究
关于
在每一点的弯曲度，则曲率公式
基于同样原理，中心点坐标是分散的点，创建差值曲率公式，来研究弯曲程度。
定义：
，
其中，
，
因此，当曲率
时，表示弯曲程度小，当曲率
时表示弯曲程度大，利用MATLAB计算得出1986年塔形的弯曲度为：
（0.0000，-0.0056，0.0448，-0.0031，0.0014，0.0000，-0.0001，-0.0043，0.0119，-0.0025，-0.0121，-0.0537）。
该指标数据显示塔身上方有轻微弯曲。
依次求出1996年、2009年、2011年的弯曲度，得出的规律和局部分析中的塔身的弯曲倾斜情况相一致。
所以，塔身的弯曲程度可以用曲率作为指标来衡量。
表5.2.2各年古塔的弯曲度
年份
1986年
1996年
2009年
2011年
弯曲度（曲率）
0.000014655
0.000014838
0.000015627
0.000015682
现将数据中每一年的中心坐标构成的主轴绘图描述，并根据已知中给出的点的坐标，
做出立体图像进行比较，可以明显的看出该塔的上部发生弯曲。
图5.2.1塔身整体立体图
5.2.3古塔的扭曲模型
古塔发生扭曲时，应该有一个旋转轴，我们假设古塔是绕着其中心线发生扭曲的，越靠近底面，扭曲转过的角度越小，离地面越高的地方扭曲时转过的角度越大。
除了曲率，空间曲线还存在挠率，由于挠率体现了密切平面的扭转状况，所以通常说它表示了曲线的扭曲程度，在这一问中，我们就用曲线的挠率来作为衡量古塔扭曲程度的标准。
以2011年数据为例，得到2011年中心点拟合的曲线为：
则
，
由挠率公式知：
根据公式，把方程
（1）带入进去，可以得到方程
（1）的挠率为0.001155764。
用同样的方法，我们可以拟合出其他三年的曲线，同样利用挠率公式，计算出其他三年的挠率，见下表：
表5.2.3古塔扭曲程度
年份
1986年
1996年
2009年
2011年
挠率
0.001304346
0.001304195
0.001155914
0.001155764
5.3问题三的模型建立与求解
该问题主要建立在问题二中对塔的倾斜，弯曲，扭曲等变形情况的分析基础上，预测未来塔的倾斜、弯曲、扭曲等变形趋势。
利用1986年到2011年共计4次测量的数据,结合建立的数学模型计算出古塔的倾斜角度、弯曲程度以及扭曲状况等结果，根据计算结果画出这三种情况的折线图。
5.3.1古塔倾斜的变形趋势的预测
根据问题二所求得的四个年份中古塔的倾斜角
的一系列数据，以时间
为自变量和以对应的倾斜角
为因变量，设
是
的线性函数，
它有4组观测值，使用MATLAB统计工具箱中的regress命令，计算线性回归系数结果整理如表5.2.4所示，
表5.2.4回归系数的计算结果
回归系数
回归系数的估计值
回归系数的置信区间
-9.0698
[-17.8322，-0.3073]
0.0049
[0.0005，0.0093]
=0.9208,F=23.2547,p=0.0404
的置信区间不含零点；p<0.05;
=18.5128较大。
这些充分检验了模型的有效性，可以对古塔的倾斜情况作出比较准确的预测。
拟合结果为
代入不同的年份，可以得到对应的倾斜角。
如：
2018年，古塔的倾斜角为0.8182.
5.3.2古塔弯曲的变形趋势的预测
根据问题二所求得的四个年份中古塔的弯曲度的一系列数据，以时间t为自变量和以对应的弯曲度K为因变量作出折线图，如图5.3.1所示。
图5.3.1古塔弯曲度随时间的变化图
从图7中可以看到，随着时间的推移，古塔的弯曲度有增大的趋势，但其增大的速率极小极小，使用线性回归的方法求出其1次项系数几乎为零。
仅从图中分析，曲线呈S形增长，2009年后，短期内其弯曲度没有明显增长，可以认为其弯曲度几乎没有变化，但无法对较长期的年份进行预测。
5.3.3古塔扭曲的变形趋势的预测
根据问题二所求得的四个年份中古塔的扭曲度的一系列数据，以时间t为自变量和以对应的扭曲度
为因变量作出散点图，如图5.3.2所示。
图5.3.2古塔扭曲度随时间的变化图
从图中可以看出，古塔的扭曲度随着时间变化起伏不定，类似周期变化，左右扭曲，使得其扭曲度时而增大，时而减小，可以预测2011年后的3、5年内其扭曲度依然有减小的趋势，更长时期后的情况，则因为数据太少，无法作出准确的预测。
六、模型的评价与推广
6.1模型的评价
优点：
通过充分利用题目中所给出的各次测量时所有的的观测数据，把误差降低到很小，使最后的结果更有说服力，并且可以很清楚看出古塔历年来的变化程度，也利于对古塔以后的变形趋势更加确定的推测。
缺点：
在讨论古塔的变形趋势时，没有研究数据的异常点，没有进行误差分析，一致的认为就是真实点，这对他的变形趋势会产生一定的误差，导致趋势预测不准确。
对题目中进行补充的缺失数据与实际测量的数据可能存在一定的误差，如果将模型运用于实践中，要考虑的因素有很多，例如大风，地震，海啸等都会对古塔的变形趋势有很大的影响。
6.2模型的推广
本模型可以应用于建筑物测量及变形的研究中，通过对建筑物变形程度的研究监控，及时了解建筑物的情况，避免危险的突然发生给社会带来较大的损失和危害；还可应用于预测长期在电脑前工作人员的脊椎的变形状况，通过记录，及时了解脊椎的变形趋势，有利于身体健康。
参考文献：
[1]梁海奎.古塔变形测量方法探讨.《城市勘探》2011年第03期
[2]姜启源、谢金星、叶俊，数学模型（第四版）[M]，北京：
高等教育出版社，
2011.
[3]胡志晓，古塔倾斜观测和数据分析，江苏建筑，第145期:
34-35,2011.
[4]李金泉、胡风学，高耸建筑物倾斜变形的观测方法和规律分析[J]，2008年全
国“三下”采煤学术会议论文集，52-55，2008.
[5]吕林良，解析几何，北京：
高等教育出版社，2006.5.
附录：
1.古塔八边形的轮廓图：
A1=xlsread（'c:
\附件1','sheet1','c4:
c108'）;
B1=xlsread（'c:
\附件1','sheet1','d4:
d108'）;
C1=xlsread（'c:
\附件1','sheet1','e4:
e108'）;
A2=xlsread（'c:
\附件1','sheet1','j4:
j108'）;
B2=xlsread（'c:
\附件1','sheet1','k4:
k108'）;
C2=xlsread（'c:
\附件1','sheet1','l4:
l108'）;
plot3（A1,B1,C1,'.-'）;
holdon
plot3（A2,B2,C2,'r.-'）;
2.1986年各层中心
clearall;
center=zeros（13,3）;
source=xlsread（'c:
\古塔1986','sheet1','c4:
e107'）;
fori=1:
13
data=source（（i-1）*8+1:
（i-1）*8+8,1:
2）;
f=@（p,data）（data（:
1）-p
（1））.^2+（data（:
2）-p
（2））.^2-p（3）^2;
p=nlinfit（data,zeros（size（data,1）,1）,f,[10010010]'）;
center（i,1）=p
（1）;
center（i,2）=p
（2）;
sum=0;
forj=1:
8
sum=sum+source（（i-1）*8+j,3）;
end
center（i,3）=sum/8;
end
3.求与Z轴夹角，倾斜度
clearall;
X=xlsread（'c:
\古塔1986年','中心','b2:
b15'）;
Y=xlsread（'c:
\古塔1986年','中心','c2:
c15'）;
Z=xlsread（'c:
\古塔1986年','中心','d2:
d15'）;
a=polyfit（X,Z,1）;
b=polyfit（Y,Z,1）;
XL1=[a
（1）,0,-1];
XL2=[0,b
（1）,-1];
p=XL1
（2）*XL2（3）-XL1（3）*XL2
（2）;
q=XL1（3）*XL2
（1）-XL1
（1）*XL2（3）;
r=XL1
（1）*XL2
（2）-XL1
（2）*XL2
（1）;
mol=sqrt（p*p+q*q+r*r）;
p=p/mol;
q=q/mol;
r=r/mol;
p=sign（r）*p;
q=sign（r）*q;
r=sign（r）*r;
arfa=acos（p）;
darfa=arfa*180/pi;
beit=acos（q）;
dbeit=beit*180/pi;
gama=acos（r）;
dgama=gama*180/pi;
XL=[X
（1）-X（14）,Y
（1）-Y（14）,Z
（1）-Z（14）];
moxl=sqrt（XL*XL'）;
XL=XL/moxl;
XL=XL*sign（XL（3））;
arfaline=acos（XL
（1））;
darfaline=arfaline*180/pi;
beitline=acos（XL
（2））;
dbeitline=beitline*180/pi;
gamaline=acos（XL（3））;
dgamaline=gamaline*180/pi;
tajianlianxian=[arfaline,beitline,gamaline];
nihezhixian=[arfa,beit,gama];
qingxiedu=zeros（13,3）;
fori=1:
13
XL=[X（i+1）-X（i）,Y（i+1）-Y（i）,Z（i+1）-Z（i）];
moxl=sqrt（XL*XL'）;
XL=XL/moxl;
XL=XL*sign（XL（3））;
ar=acos（XL
（1））;
dar=ar*180/pi;
bei=acos（XL
（2））;
dbei=bei*180/pi;
gam=acos（XL（3））;
dgam=gam*180/pi;
qingxiedu（i,1）=ar;
qingxiedu（i,2）=bei;
qingxiedu（i,3）=gam;
end
xlswrite（'c:
\古塔1986','倾斜度','b2:
d14'）
xlswrite（'c:
\古塔1986','倾斜度','b16:
d16'）
xlswrite（'c:
\古塔1986','倾斜度','b18:
d18'）
升级会员 