2014全国大学生数学建模竞赛A题论文.pdf

《2014全国大学生数学建模竞赛A题论文.pdf》由会员分享，可在线阅读，更多相关《2014全国大学生数学建模竞赛A题论文.pdf（18页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、1 2014高教社杯全国大学生数学建模竞赛承诺书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题我们知道，抄袭别人的成果是违反竞赛规则的,如果引用别人的成果或其他公开的资料（包括网上查到的资料)，必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性如有违反竞赛规则的行为，将受到严肃处理我们参赛选择的题号是（从A/B/C/D 中选择一项填写）赛区评阅编号（由赛区组委会评阅前进行编号）：2 201

2、4高教社杯全国大学生数学建模竞赛编 号 专 用 页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：1 嫦娥三号软着陆轨道设计与控制策略摘要本文针对嫦娥三号软着陆轨道设计与控制策略的实际问题，以理论力学（万有引力、开普勒定律、万能守恒定律等）和卫星力学知识为理论基础，结合微分方程和微元法，借助MATLAB 软件解决了题目所要求解的问题。针对问题（1），在合理的假设基础上，利用物理理论知识、解析几何知识和微元法，分析并求解出近月点和远月点的位置，即 139.1097

3、 。再 运用能 量守 恒定 律 和相 关数 据，计 算 出速度1v（近 月点 的 速度）=1750.78/m s,2v（远月点的速度）=1669.77/m s，最后利用曲线的切线方程，代入点（近月点与远月点）的坐标求值，计算出方向余弦即为相应的速度方向。针对问题（2）关键词：模糊评判，聚类分析，流体交通量，排队论，多元非线性回归一、问题重述嫦娥三号于 2013 年 12 月 2 日 1 时 30 分成功发射，12 月 6 日抵达月球轨道。嫦娥三号在着陆准备轨道上的运行质量为2.4t，其安装在下部的主减速发动机能够产生1500N到 7500N的可调节推力，其比冲（即单位质量的推进剂产生的推力）为

4、2940m/s，可以满足调整速度的控制要求。在四周安装有姿态调整发动机，在给定主减速发动机的推力方向后，能够自动通过多个发动机的脉冲组合实现各种姿态的调整控制。嫦娥三号的预定着陆点为 19.51W，44.12N，海拔为-2641m（见附件 1）。嫦娥三号在高速飞行的情况下，要保证准确地在月球预定区域内实现软着陆，关键问题是着陆轨道与控制策略的设计。其着陆轨道设计的基本要求：着陆准备轨道为近月点15km，远月点 100km的椭圆形轨道；着陆轨道为从近月点至着陆点，其软着陆过程共分为6 个阶段（见附件 2），要求满足每个阶段在关键点所处的状态；尽量减少软着陆过程的燃料消耗。根据上述的基本要求，请你

5、们建立数学模型解决下面的问题：（1）确定着陆准备轨道近月点和远月点的位置，以及嫦娥三号相应速度的大小与方向。（2）确定嫦娥三号的着陆轨道和在6 个阶段的最优控制策略。2（3）对于你们设计的着陆轨道和控制策略做相应的误差分析和敏感性分析。二、问题分析2.1 问题（1）的分析首先根据问题的假设、题目中所提供的数据及图片分析，可以知道嫦娥三号绕月球的轨道是由圆形轨道变为椭圆形轨道，借助开普勒定律、能量守恒定律求解出近月点的速度。为了确定近月点和元月点的精确位置及相应的速度方向，我们建立以赤道（月球的赤道）平面为xoy平面、月心为原点、月心与零度经线和零度纬线交线的交点的连线为坐标轴的坐标系和赤道（月

6、球的赤道）平面为xoy平面，为极轴（月球的极轴）为z轴建立空间直角坐标系，x轴与极坐标系的轴相重合。首先根据着陆点的经度、纬度及月球的半径求解出着陆点和近月点（带参数)的空间直角坐标。其次利用两点间的距离公式，并借助MATLAB 软件求解出近月点与着陆点最短距离。从而计算出（近月点的经度）=。最后利用卫星的轨迹是以月心为其中一个焦点，以近月点与远月点的距离为长轴的椭圆，从而求解出卫星的轨迹方程，再运用隐函数求导的应用的知识，求解出在近月点和远月点的方向导数，进而求解近月点和远月点方向余即为近月点和远月点的速度的方向。2.2 问题（2）的分析首先在根据题意，将嫦娥三号软着陆问题，分为 6 个阶段

7、依次为主减速、快速调整、粗避障、精避障、缓慢下降、自由下降，我们先将6 个阶段分为 4 个阶段，依次为第一阶段（主减速和快速调整）、第二阶段（粗避障）第三阶段(精避障),第四阶段（缓慢下降和自由下降）。其次在第一阶段粗避障阶段，嫦娥三号悬停在月球表面约2400米上方，对星下月表进行二维和三维成像，利用遗传算法的思想，从图像中先随机选取部分点，能直接从三维图像中得知该点的海拔高度，再分别扫描这些点附近的地貌，找出一些地势平坦的区域，我们用区域内所有点与中心点海拔的均方差作为地势判断依据之一，保留这些坐标，并进行重新组合，并改变某些坐标以便能获得其他新区域的坐标，再次搜索地势平坦的区域，重复进行多

8、次搜索，直到没有出现崎岖地势的时候，我们将此时地势最平坦的地方作为全局最优降落地点3 三、模型假设1、不考虑空间飞行器上各点因燃料消耗而产生的位移；2、在对卫星和空间飞行器进行轨道估计时，认为作用于其上的所有外力都通过其质心；3、卫星和空间飞行器的运动是在真空中进行的；4、卫星只受重力影响，空间飞行器除自身推力外只受重力影响；5、卫星的观测图片及数据精准；6、四、变量与符号说明0C一条车道的基本通行能力L连续车流的车头间距C n 条车道的基本通行能力y排队长度1x车流量2x横断面通行能力系数车流量3x持续时间五、模型建立与求解5.1 问题(1)的分析、模型建立与求解5.1.1 建模准备（1）开

9、普勒定律开普勒第一定律开普勒第一定律开普勒第一定律，也称椭圆定律：每一个行星都沿各自的椭圆轨道环绕太阳，而太阳则处在椭圆的一个焦点中。开普勒第二定律开普勒定律开普勒第二定律，也称面积定律：在相等时间内，太阳和运动着的行星的连线所扫过的面积都是相等的。这一定律实际揭示了行星绕太阳公转的角动量守恒。用公式表示为开普勒定律开普勒第4 三定律开普勒定律开普勒第三定律，也称调和定律：各个行星绕太阳公转周期的平方和它们的椭圆轨道的半长轴的立方成正比。由这一定律不难导出：行星与太阳之间的引力与半径的平方成反比。这是牛顿的万有引力定律 的 一 个 重 要 基 础。用 公 式 表 示 为32aKT开 普 勒 定

10、 律这里，是行星公转轨道半长轴，是行星公转周期，是常数。（2）万有引力万有引力：任意两个质点有通过连心线方向上的力相互吸引。该引力大小与它们质量的乘积成正比与它们距离的平方成反比，与两物体的化学组成和其间介质种类无关。即：122M MFGr，其中12MM，为两物体的质量，11226.6710.GN m kg（牛顿每平方米二次方千克）5.1.2 模型的建立根据以上的分析，建立以月球赤道平面为xOy平面,月心为原点O、Ox为月心与零度经线和零度纬线交线的交点的连线，Oz为极轴（月球的极轴），Oy与Ox和Oz满足右手标架，建立空间直角坐标系(如图 5-1 所示）。图 5-1 卫星绕月轨迹及软着陆轨迹

11、由于着陆点在球面上且近月点与远月点是由月球的经度、纬度及高度唯一确定，在此为了便于计算将极坐标转化为空间直角坐标，并代数题中相关数据，反解出经度。极坐标转化为空间直角坐标即：sincossinsincosxryrzr（5.1.1）5 sin(90-sin(90-cos(90-xryrzr）cos(-)）cos(-)）（5.1.2）距离公式：222000()()()dxxyyzz（5.1.3）其中：为纬度；为经度；r为嫦娥三号距月心的距离；d为嫦娥三号距着陆点的距离；根据能量守恒、开普勒第二定律（面积定律），建立以下模型即：112222121122r vr vmvmghmvmgH（5.1.4）则

12、近月点的速度，近月点的速度：221222121222212()2()g Hh rvrrg Hh rvrr（5.1.5）其中：m为卫星的质量，1h为海拔高度，h近月点距月球表面的距离；101rhrh，201rHrh，0r月球半径，H远月点距月球表面的距离，g月球重力加速度，1v近月点的速度，2v近月点的速度。5.1.3 模型的求解5.1.3.1 近月点与远月点的位置根据题目所给数据以上分析，可知：010,15000,1737013,2641hm rm hm将以上数据代入（5.1.1）式可得，着陆点及近月点的空间直角坐标分别为：000000000sin(90)cossin(9019.51)cos4

13、4.12sin(90)sinsin(9019.51)sin44.12cos(90)cos(9019.51)xrryrrzrr（5.1.6）00sin(90-)cossin(90-)sincos(90-xrhyrhzr）cos(-)=(r）sin(-)=-(r）=0（5.1.7）6 再将（5.1.6）式和（5.1.7）式代入（5.1.3）式可得关于与d（近月点和 着陆点 距 离）的 函 数，？利 用 Mathematica 5.0编 程求 解 可 得：-139.1075.1.3.2近月点与远月点的速度大小及方向近月点与远月点的速度方向，即为相应速度在x轴与y轴方向上的投影（如图 5-2 所示）图

14、 5-2 近月点与远月点的速度方向示意图由图易知：5.2 模型二的建立5.2.1 模型准备5.2.1.1系统模型 1、着陆器的动力下降段一般从15km左右的轨道高度开始，下降到月球表面的时间比较短，在几百秒范围内，所以可以不考虑月球引力摄动。月球自转速度比较小，也可忽略。因此，可以利用二体模型描述系统的运动。建立图 5-2 所示的着陆坐标系，并假设着陆轨道在纵向平面内，令月心为坐标原点，Oy指向动力下降段的开始制动点，Ox指向着陆器的开始运动方向。则着陆器的质心动力学方程可描述如下:22(/)sin/(/)cos2/SPrvvFmrrFmvrmFI式中:,r和m分别为着陆器的月心距、极角、角速

15、度和质量；v为着陆器沿 r 方向上的速度；F 为制动发动机的推力(固定的常值或 0）；SPI为其比冲；为月球引力常数；为发动机推力与当地水平线的夹角即推力方向角。7 图 5-3 月球软着陆坐标系动力下降的初始条件由霍曼变轨后的椭圆轨道近月点确定，终端条件为着陆器在月面实现软着陆。令初始时刻00t，终端时刻ft不定，则相应的初始条件为0000,0,Lorrh v终端约束为,0,0fLffrrv式中:Lr为月球半径；0h为初始轨道高度；o为轨道角速度。月球软着陆的最优轨道设计就是要在满足上述初始条件和终端约束的前提下，调整推力大小和方向9使得着陆器实现燃料最优软着陆，即要求以下性能指标达最大。0f

16、tJmdt5.2.1.2 模型归一化在轨道优化过程中，由于各状态变量的量级相差较大，寻优过程中可能会导致有效位数的丢失。通过归一化处理可以克服这一缺9点，提高。计算精度。令00,reftefrrmm，则/,/,/,/refrefrefrefSpSprefrrrvv vvrIIr8 23/,/,/,/,/refrefrefrefrefrefrefrefFFFFmvrmm mrttt,/,refreftrv。那么，着陆器的动力学方程可改为：22(/)sin/(/)cos2/SPrvvFmrrFmvrmFI相应的初始条件和终端约束变为：300001,0,/orvr10/,0,0fffrrr v性能指

17、标改写为：0ftJmdt5.2.4 模型评判由以上计算可知，多车道道路通行能力从中心至边缘车道依次递减.视频一中撞车位置在距道路中心一、二条车道上，因而可行车道为第三条车道；视频二中撞车位置在距道路中心二、三条车道上，因而可行车道为第一条车道.而从计算中可得，可上述结论，即视频二的事故所处断面实际通行能力要比视频一要强.这与实际情况比较吻合。5.3 问题（3）的模型建立与求解当上游交通需求量大于事发路段现有的通行量，到达车流在事故地点陆续减慢速度甚至停车而集结成密度较高的队列，事故接触后，由于数段通过能力的恢复，排队车辆又陆续加速而疏散成一列具有适当密度的车队，排队消散完毕后，车流就会恢复顺畅

18、的交通状态.9 5.3.1 模型推导由车流波动理论可知，波速公式为：,XYXY(QQ)/(KK).X YW（1）式中：xyW，为集散波的波速，Km/h；XQ、YQ为前后两种车流状态的流量，辆/h；XK、YK为前后两种车流状态的密度，辆/Km.根据交通流模型可知，交通量Q、行车速度 v、车流密度 K 三者的关系为：QK.v（2）速度-密度线性关系模型：fj(1 K/K.vv）（3）式中：fv为畅行速度，即车流密度为零时,车辆的最大速度；jK为阻塞密度，即车流密集到所有车辆无法移动时的密度；由以上（1）、（2）、（3）式可以推导出波速与密度的关系：xy,fjK/K(1.KX YWv）（4）5.3.

19、2 模型的建立与求解事故发生后排队长及消散时间的计算图为事故发生后累计车辆-时间图，实线表示交通需求流量，点划线表示通过能力.为叙述简便，对所有符号说明如下：事故发生时堵塞了部分车道，该路段通行能力下降1S；相应密度上升1sK；交通事故处理所需时间为；事故解除后到车队消散前通行能力回升为2S；车流密度相应地下降为2sK.其中路段的通行能力由图中点划线的斜率来表示.路段上游交通需求流量为12QQ、.由图中实线斜率表示；持续时间为12TT、；相应车流密度为12KK、.在图中，由车流波动理论可知，波速公式为：,XYXY(QQ)/(KK).X YW首先假设两波相遇之前该路段需求量始终未变，OA与 CB

20、 相交处表示排队向上游的延伸达到的最远处，设两波相遇时的时间为T，集结波波速为12W，消散波波速为2,3W，10 则根据两波相遇时波传动的距离相等这一关系可知：12230(TT).WTW (5）其中11111211=|(1)|;sfsfKKQSWVKKK12122312=|(1)|.ssfssfKKSSWVKKK则:T=122f02V(1)(T1(T).SSsfKKk xxKSQ若 T 1T，则说明在车队消散之前该路段上游要求流量发生了变化，需求流量变为2Q，相应的密度为2K.所以（5）式改写为121341230(TT)W(TT).WTW (6)其中212134f21WV(1).SSfKKQS

21、KKK则12021112()T.SISSKKKTKKTKK）根据公式；121f022(V(1)(TT)K/SSSSfKKTTK可解出本次事故引起的排队长.由资料可知车队消散时间为:13/.TTx V其中3V为路段通行能力为2S时的行车速度，22/;SvmSK,XYXY(QQ)/(KK);X YWQK.v所以，根据以上推导可以得到排队长度（x）关于横断面实际通行能力（2S）、拥堵持续时间（消散时间1T）、路段上游车流量(Q)的关系式如下：122f02V(1)(T1(T).SSsfKKk xxKsQ5.4 问题（4）的模型建立与求解5.4.1 数值分析确立各个参数值在此问中我们采用非线性回归法：非

22、线性回归是指因变量y 对于1.m（不是自变量）是非线性.MATLAB 统计工具箱的命令 nlinfit等不仅可以给出拟合的回归系数及其置信区间，而且可以给出预测值及置信区间等.多元线性回归分析模型:11 其中都是与无关的未知参数，称为线性回归系数.现在得到 n个独立观测数据，其中为 y的观测值，分别为的观测值，由上式得:在进行非线性回归的时候，道路通行能力以道路通行系数表示.我们选取了视频 1情况下的六次排队情况（由于每次持续时间较短，因此可以假设排队长度为120，）道路通行系数的选取见问题一的图像.排队长度 y车流量1x横断面通行能力系数2x持续时间3x120 10 0.6 10 120 8

23、5 0.5 300 120 55 0.6 172 120 24 0.5 65 120 26 0.8 60 120 30 0.6 90 120 90 0.7 560 第 一 步：以 回 归 系 数 和 自 变量 为 输 入 变 量，将 要 拟 合 的 模 型 写 成 匿 名 函数:1121233()/)()xxxy=(；第二步：用nlinfit计算回归系数，用nlparci计算回归系数的置信区间，用nlpredci计算预测值及其置信区间；第三步：确立各参数，得出123、的参考值为0.0001,1.0243,-0.3253所以123.(0.0001 1.0243)/)(-0.3253).xxxy=

24、(000015.4.2 模型建立与求解基于第三问的模型建立求解过程：由第三问我们可以得出：123.(0.0001 1.0243)/)(-0.3253).xxxy=(00001在该方程中，有三个变量，分别为道路上游车流量、通行系数以及堵车时间，假设从开始排队到排队位置到达上游路口的时间内上游车流量保持不变，即为1500pcu/h，假设该车流为连续流，于是设定2x为 0.5.解该等式：可得出，堵车时间间隔大约为8 分钟，即从事故发生开始，经过大约8 分钟后排队长度将达到上游六、模型的分析与评价在模型一中，综合考虑了驾驶时间、车身长度、车道折算系数等多种因素，使得计算出来的畅通率更加的正确.但是我们

25、为了简化模型，我们没有考虑到模型的岔道对车流的影响.在问题二的模型中，用模糊评判法对附件1 和附件 2 两种状态下的通行能力进行定性的研究；用聚类分析模型对附件1和附件 2两种状态下的通行能力进行定量上的分析，理论严谨.在问题三中，通过建立基于车流波动理论的交通流模型，推导过程思路清晰.但是12 我们是把车流运输当做连续的流体处理，这个因素使得结果有些误差.在问题四中的模型中，通过问题三所得关系式，进行参数分析，通过多元非线性回归，得出了数值表达式，进而求出最短时间，思路明朗.但在模拟的多元非线性自身也存在一些误差.七、模型的改进由于交通问题是一个复杂的系统，可以考虑进行基于元胞自动机的交通流

26、计算机模拟研究.在问题一中，可以进一步把统计的区间段进行细分，进行差值拟合，是图像更加接近于真实情况；在问题二中，同理可以进行细分；在问题三中，引入相位变换等，可以考虑添加修正系数，获取相关数据，并用BP神经网络进行训练.八、模型的应用与推广为了更好的反映车道被占用对城市道路通行能力的影响，问题三的模型给出了车辆排队长度与车流量、间隔时间、通行能力的关系.然而在实际问题中，还有许多非确定性因素.且车辆排队长度与车流量、间隔时间、通行能力以及其他非确定性因素都存在相互影响的关系，经过相关数据提取技术，在获取相关数据后，可以结合片最小二乘回归分析进行处理：第一步：提取两变量组的第一对成分，使之相关

27、性达到最大；第二步：建立回归方程；第三步：残差替换；第四步：进行偏最小二陈回归分析法；第五步：交叉有效性检验.参考文献：1 姜启源数学模型（M）第三版北京：高等教育出版社，20038.2 张德丰 MATLAB 数值计算方法（M）北京：机械工业出版社，20101.3 汪晓银，周保平数学建模与数学实验（M）北京：科学出版社，20102.4孔惠惠，秦超，李新波，李引珍，交通事故引起的排队长度及消散时间的估算，第二十七卷，第 5期，2004-12-23 5 石磊等，数学建模优秀论文，北京：清华大学出版社，2011 年 9 月第 1 版6 姜启源，数学模型（第二版），北京：高等教育出版社，1993 年

28、8 月第 2 版7 王炜，过秀成等编著交通工程学南京：东南大学出版社，2000 8 杨肇夏计算机模拟及其应用北京：中国铁道出版社，1999 9 谌红模糊数学在国民经济中的应用武汉：华中理工大学出版社，1994 1 附录附录一：程序1问题二中模型一：模糊评判模型程序一（源代码：p1.m）clc,clear a=0.1 0.3 0.4 0.1 0.1 0.2 0.4 0.4 0 0 0.1 0.3 0.5 0.1 0 0 0.1 0.1 0.4 0.4 0 0 0.1 0.4 0.5 0.2 0.4 0.3 0.1 0;w=0.6 0.4;w1=0.5 0.25 0.25;w2=0.5 0.25

29、0.25;b(1,:)=w1*a(1:3,:);b(2,:)=w2*a(4:6,:);c=w*b;程序二（源代码：p2.m）clc,clear a=0.1 0.3 0.4 0.2 0 0 0.4 0.2 0.2 0.2 0.2 0.3 0.4 0.1 0 0.1 0.1 0.1 0.4 0.3 0.1 0.1 0.1 0.4 0.3 0.2 0.4 0.2 0.1 0.1;w=0.6 0.4;w1=0.25 0.25 0.5;w2=0.25 0.25 0.5;b(1,:)=w1*a(1:3,:);b(2,:)=w2*a(4:6,:);c=w*b 2问题二中模型二：聚类分析视频 1 程序代码（源

30、代码：moxinger.m）clc,clear load hua.txt gj=zscore(hua);%数据标准化y=pdist(hua);%求对象间的欧氏距离,每行是一个对象2 z=linkage(y,average);%按类平均法聚类h=dendrogram(z);%画聚类图set(h,Color,k,LineWidth,1.3)%把聚类图线的颜色改成黑色，线宽加粗for k=3:5 fprintf(划分成%d类的结果如下：n,k)T=cluster(z,maxclust,k);%把样本点划分成 k 类 for i=1:k tm=find(T=i);%求第 i 类的对象 tm=resha

31、pe(tm,1,length(tm);%变成行向量 fprintf(第%d类的有%sn,i,int2str(tm);%显示分类结果 end if k=5 break end fprintf(*n);end 视频二程序代码（源代码：p4.m）：clc,clear load hu.txt gj=zscore(hu);%数据标准化y=pdist(hu);%求对象间的欧氏距离,每行是一个对象z=linkage(y,average);%按类平均法聚类h=dendrogram(z);%画聚类图set(h,Color,k,LineWidth,1.3)%把聚类图线的颜色改成黑色，线宽加粗for k=3:5 f

32、printf(划分成%d类的结果如下：n,k)T=cluster(z,maxclust,k);%把样本点划分成 k 类 for i=1:k tm=find(T=i);%求第 i 类的对象 tm=reshape(tm,1,length(tm);%变成行向量 fprintf(第%d类的有%sn,i,int2str(tm);%显示分类结果3 end if k=5 break end fprintf(*n);end 3问题四中模型：非线性回归法（源代码：p5.m）xy0=120 1200 0.6 10120 930 0.5 300120 840 0.4 172120 760 0.5 65120 110

33、0 0.3 60120 900 0.6 90120 890 0.7 560;x=xy0(:,2:4);y=xy0(:,1);huaxue=(beta,x)(beta(1)-x(:,1).*beta(1).*beta(2)./x(:,2).*(beta(3).*x(:,1)-x(:,3);beta0=0.05,0.005,0.02;beta,r,j=nlinfit(x,y,huaxue,beta0)betaci=nlparci(beta,r,jacobian,j)yhat,delta=nlpredci(huaxue,x,beta,r,jacobian,j)运行结果：beta=0.0001 1.0

34、243-0.3253 r=28.4765-8.0898 13.1246 65.7630-55.1525 54.3847-3.4178 j=1.0e+006*0.8056 0.0001-0.0003 1.1274 0.0001-0.0002 0.9407 0.0001-0.0002 0.4774 0.0001-0.0001 1.5416 0.0002-0.0005 4 0.5775 0.0001-0.0002 1.0863 0.0001-0.0001 betaci=1.0e+003*-0.0005 0.0005-4.2144 4.2165-0.0015 0.0008 yhat=91.5235 128.0898 106.8754 54.2370 175.1525 65.6153 123.4178 delta=1.0e+004*1.3405 1.8908 1.5770 0.7937 2.5870 0.9604 1.8150