第5讲微分方程建模

数学建模算法与应用第6章

微分方程建模数学建模算法与应用6.4放射性废料的处理6.4.1问题的提出6.4.2模型的建立与求解subs是赋值函数，用数值替代符号变量替换函数

例如:输入subs(a+b,a,4)意思就是把a用4替换掉，返回4+b

也可以替换多个变量，例如：subs(cos(a)+sin(b),{a,b},[sym('alpha'),2])

分别用字符alpha替换a和2替换b，返回cos(alpha)+sin(2)vpa函数是什么意思？除了利用MATLAB 的help还可以借助于网络资源.simplify和simple是Matlab符号数学工具箱提供的两个简化函数，区别如下：

simplify的调用格式为：simplify(S)；对表达式S进行化简。

simple是通过对表达式尝试多种不同的方法（包括simplify）进行化简，以寻求符号表达式S的最简形式。6.5初值问题的Matlab数值解数学建模竞赛赛题讲评全国大学生数学建模竞赛1992年中国工业与应用数学学会(CSIAM)开始组织.1994年起教育部高教司和CSIAM共同举办(每年9月).

全国高校规模最大的学科竞赛内容

赛题：工程技术、管理科学中经过简化的实际问题

答卷：一篇包含模型假设、建立、求解、计算方法设计和计算机实现、结果分析和检验、模型改进等方面的论文形式3名大学生组队，在3天内完成的通讯比赛

可使用任何“死”材料（图书、计算机、软件、互联网等），但不得与队外任何人讨论宗旨创新意识团队精神重在参与公平竞争标准假设的合理性，建模的创造性，结果的正确性，表述的清晰性年份A题B题C题D题2004奥运会临时超市网点设计电力市场的输电阻塞管理饮酒驾车公务员招聘2005长江水质的评价和预测DVD在线租赁雨量预报方法的评价DVD在线租赁2006出版社的资源配置艾滋病疗法的评价和疗效的预测易拉罐形状和尺寸的最优设计煤矿瓦斯和煤尘的监测与控制2007中国人口增长预测

乘公交，看奥运手机“套餐”优惠几何

体能测试时间安排

2008数码相机定位高等教育收费标准探讨地面搜索NBA赛程的分析与评价2009制动器试验台的控制方法分析眼科病床的合理安排卫星和飞船的跟踪测控会议筹备2010储油罐的变位识别与罐容表标定上海世博会影响力定量评估输油管的布置学生宿舍设计方案评价

题目的特点题目来源:实际研究课题的简化、改编；有实际背景问题的编撰；合适的社会热点（或兴趣）问题题目背景尽量通俗易懂，涉及的专业知识不深题目需要的数学知识一般不超过本科的三门主干课（非数学专业）内容及统计、优化、计算等基本方法；专科题目力求少用大学数学内容解题所用的数学方法尽量多元化、综合化可以查阅到一些参考材料，但是无法照搬现成文献兼顾数据的处理与数据的收集基础性实践性综合性开放性挑战性A题：中国人口增长预测题目分析：题目属于那种类型：连续的、离散的？需要解决什么问题；最优化方案、预测模型、最短路径等等；可以用哪些相关模型、算法求解、需要什么数学工具；论文写作：（1）

问题提出和假设的合理性

论文中的假设要以严格、确切的数学语言来表达;所提出的假设确实是建立数学模型所必需的;假设应验证其合理性.（2）

模型的建立在作出假设后，我们就可以在论文中引进变量及其记号，抽象而确切地表达它们的关系，通过一定的数学方法，最后顺利地建立方程式或归纳为其他形式的数学问题，把得到数学模型的过程表达清楚，使读者获得判断模型科学性的一个依据.（3）模型的计算与分析把实际问题归结为一定的数学问题后，就要求解或进行分析.在数值求解时应对计算方法有所说明，并给出所使用软件的名称或者给出计算程序（通常以附录形式给出）.还可以用计算机软件绘制曲线和曲面示意图，来形象地表达数值计算结果。基于计算结果，可以用由分析方法得到一些对实践有所帮助的结论.（4）

模型的讨论

对所作的数学模型，可以作多方面的讨论.例如可以就不同的情景，探索模型将如何变化.或可以根据实际情况，改变文章一开始所作的某些假设，指出由此数学模型的变化.还可以用不同的数值方法进行计算，并比较所得的结果.有时不妨拓广思路，考虑由于建模方法的不同选择而引起的变化.通常，应该对所建立模型的优缺点加以讨论比较，并实事求是地指出模型的使用范围.2007A题近年来中国的人口发展出现了一些新的特点，例如，老龄化进程加速、出生人口性别比持续升高，以及乡村人口城镇化等因素，这些都影响着中国人口的增长。2007年初发布的《国家人口发展战略研究报告》(附录1)还做出了进一步的分析。关于中国人口问题已有多方面的研究，并积累了大量数据资料。附录2就是从《中国人口统计年鉴》上收集到的部分数据。试从中国的实际情况和人口增长的上述特点出发，参考附录2中的相关数据（也可以搜索相关文献和补充新的数据），建立中国人口增长的数学模型，并由此对中国人口增长的中短期和长期趋势做出预测；特别要指出你们模型中的优点与不足之处。实际或竞赛的建模问题课程学习与培训中的数学建模案例

引用案例

借鉴案例类比、引用和借鉴是很有效的方法：现状！数学建模方法=创造性的本质方法+模型类比方法数学建模=数学荐模+数学引模——一种初步实用的，基于类比与经验的方法但是：Logistic模型(阻滞增长模型)

基于Logistic模型，建立了含市、镇、乡人口相互流动关系的微分方程模型，求得全国总人口数在短期内将持续增长，到2010年、2020年分别为13.59亿和14.44亿，具有较好的中短期预测效果.1.中国人口增长的中短期预测从网上可查得2000-2005年的总人数如下表1表1年份/年总人数/万2000126583200112762720021284532003129227200412998820051307562006131448(1)数据预处理

题中所给5年我国人口1%调查数据是对人口的抽样调查数据，由于数据的不完备性，并不能由它来估计当时的全国总人口数.但基于抽样调查的等概率性，可以认为它所反应的市、镇、乡三个地区的人口比例及男女比例是与实际较为接近的.根据网上得到的具体的全国总人口数.进而可以得到各部分人口数.所得数据见表2。例：

表2

（2）模型建立考虑单一地区只受自然资源和环境条件限制Logistic模型：

其中为固有增长率.,分别表示第t年的人口数和此条件限制下的最大人口容量.因子体现了人口自身的增长趋势，因子则体现了资源和环境对人口增长的阻滞作用，人口增长是两个因子共同作用的结果.现在考虑市、镇、乡三个地区的总人口变化规律，假定在一段时期内，这三个地区各自均处于稳定状态，出生率和死亡率都没有太大变化。然而居民在三个地区间的迁移会对总人口产生一定的影响。对原Logistic模型进行修改，得到多地区条件下的改进模型：（3）模型求解未知参数为：，，：3.1求解r3.23.3表3迁移系数3.4表4模型对全国总人口的预测（4）模型分析结果分析由表4可以看出，我国人口在短时间内仍会继续增长，增长水平为每年1000万人左右，到2010年和2020年，全国总人口将达13.59亿人与14.44亿.2.灵敏性分析3.模型优点基于传统的Logistic模型，建立起市、镇、乡人口间的相互流动关系，综合考虑了出生率、死亡率与迁入迁出对人口增长的影响，对Logistic原模型做了合理的扩充，使其更符合实际情况，预测效果良好.同时根据实际情况，对难以求解的微分方程采用差分方程代替，降低了解题难度，使模型易于实现.4.模型缺点i.由于采用了差分方程的近似代替，整个预测结果均随着时间的增长而增大，这是不符合生物学规律的，故此模型对人口的中短期预测较好，长期预测则难以让人信服.ii.没有考虑性别比对年龄结构的影响，且只能对总人口数进行预测，无法对性别比和年龄结构做出相应的预测.5.改进方向：采用先进的工具软件对微分方程求解，得到模型预测的准确结果.ii.分段对模型进行调整和改进，减少累加误差.theend储油罐的变位识别与罐容表标定全国大学生数学建模竞赛2010年A题

题目题目的背景与分析问题（1）的通常做法问题（2）的通常做法一种建立储油量模型的新方法对学生论文的评述题目的启示

通常加油站都有若干个储存燃油的地下储油罐，并且一般都有与之配套的“油位计量管理系统”，采用流量计和油位计来测量进/出油量与罐内油位高度等数据，通过预先标定的罐容表（即罐内油位高度与储油量的对应关系）进行实时计算，以得到罐内油位高度和储油量的变化情况。许多储油罐在使用一段时间后，由于地基变形等原因，使罐体的位置会发生纵向倾斜和横向偏转等变化（以下称为变位），从而导致罐容表发生改变。按照有关规定，需要定期对罐容表进行重新标定。图1是一种典型的储油罐尺寸及形状示意图，其主体为圆柱体，两端为球冠体。图2是其罐体纵向倾斜变位的示意图，图3是罐体横向偏转变位的截面示意图。题目油油浮子出油管油位探测装置注油口检查口地平线2m6m1m1m3m油位高度图1储油罐正面示意图油位探针油位探针α地平线图2储油罐纵向倾斜变位后示意图油油浮子出油管油位探测装置注油口检查口水平线图3储油罐截面示意图（b）横向偏转倾斜后正截面图β地平线垂直线油位探针（a）无偏转倾斜的正截面图油位探针油位探测装置3m(b)小椭圆油罐截面示意图α油油浮子出油管油位探针注油口水平线2.05mcm0.4m1.2m1.2m1.78m(a)小椭圆油罐正面示意图图4小椭圆型油罐形状及尺寸示意图

请你们用数学建模方法研究解决储油罐的变位识别与罐容表标定的问题。（1）为了掌握罐体变位后对罐容表的影响，利用如图4的小椭圆型储油罐（两端平头的椭圆柱体），分别对罐体无变位和倾斜角为=4.10的纵向变位两种情况做了实验，实验数据如附件1所示。请建立数学模型研究罐体变位后对罐容表的影响，并给出罐体变位后油位高度间隔为1cm的罐容表标定值。（2）对于图1所示的实际储油罐，试建立罐体变位后标定罐容表的数学模型，即罐内储油量与油位高度及变位参数（纵向倾斜角度和横向偏转角度

）之间的一般关系。请利用罐体变位后在进/出油过程中的实际检测数据（附件2），根据你们所建立的数学模型确定变位参数，并给出罐体变位后油位高度间隔为10cm的罐容表标定值。进一步利用附件2中的实际检测数据来分析检验你们模型的正确性与方法的可靠性。流水号C进油量/LD油位高度/mm采集时间说明1150159.022010-08-2010:32:18(1)罐体无变位进油，罐内油量初值262L；12100176.142010-08-2010:33:18(2)C列进油量是每次加入50L油后的累加值13150192.592010-08-2010:34:18(3)D列是原罐内初始油量加入相应油量后油位高度值。流水号C进油量/LD油位高度/mm采集时间说明11747.86411.292010-08-1915:10:2712797.86423.452010-08-1915:11:2713847.86438.332010-08-1915:12:3714897.86450.542010-08-1915:13:27(1)罐体纵向倾斜变位进油，纵向倾斜角4.10,罐内油量初值215L；(2)C列进油量是每次加入50L油后的累加值，个别的有小数零头；(3)D列是原罐内油量加入相应油量后油位高度值。附件1实验数据流水号C进油量/LD出油量/LE显示油高/mm显示油量容积/L采集时间备注说明201060.002632.2360448.882010-08-0108:00:49←开始采集2020149.092624.3060311.432010-08-0108:15:421)C列进油量,通常进油是一次性连续完成203068.452620.6760248.032010-08-0108:23:412)D列出油量是加油机不同时段的出油量vi2040199.272610.2960065.112010-08-0108:38:143)E列显示油高是进/出油vi后对应的罐内油位显示高度hi；205070.052606.6159999.692010-08-0108:53:084)F列显示油量容积是进/出油vi后,罐内油位5020显示高度为hi时的罐容表显示储油量Ｖi。

50360.81469.746086.74

511242486.2157704.792010-08-0811:43:402010-08-0812:52:52←一次性补充进油附件2检测数据来自于加油站设备研究与生产企业的一个实际课题.

（2）根据实际检测数据，正确识别罐体是如何变位的，具体变了多少？同时要给出罐容表的修正标定方法和结果，属于“反问题”。题目背景问题分析问题由两部分组成：（1）为了观察检验罐体变位对罐容表的影响，在已知变位参数的情况下，检测出油位高度和油量的对应数值，建模分析罐容表的变化规律，并给出修正的罐容表，属于“正问题”。问题（1）的通常做法

根据附件1的实验数据画出罐体不变位（α=0）和罐体变位（α=4.10）时油位高度h和储油量V的曲线.实验数据表明，在同样的油位高度h下，罐体变位时储油量V减少.zx0z=HabS(H)z=H下的面积油位高度h时坐标y处z=H

油位高度h时坐标y处的截面面积罐体变位α时油位高度h和储油量(体积)V的数学模型化重积分为定积分zyxhyH0L1Lα油位计L2V(h,α)z=H高度h与储油量V的模型D3D2D1yh10L1L2αh2y1y2油位计高度h与储油量V的模型V(h,α)h/cmV/Lh/cmV/Lh/cmV/Lh/cmV/L1070.1340965.66702232.501003450.7220281.86501371.88802661.421103776.6430595.25601798.52903072.431204012.74=4.10的罐容表(部分)按照模型讨论变位对罐容表的影响与实验数据表示的影响一致.变位后储油量平均约小200L.模型结果与实验数据的比较=4.10=0模型结果与实验数据基本吻合，而在同样的油位高度h下，储油量V的实验数据均比模型结果小.可以有各种解释.是否要修正、怎样修正模型？问题（2）的通常做法1.罐体变位α,β时油位高度h和储油量V的数学模型hβD1D2D3整体考虑αhVxy0z

写出左右球面、圆柱面及油位高度h的油平面方程.

将罐体分为3个区域：D1,D2,D3

在3个区域中分别写出体积的重积分，并化为（对z的）定积分.油平面将罐体分为圆柱体和球缺，圆柱体按照问题（1）考虑z=H下的面积zx0rh'hβz=H截面面积S(y,h,α,β)

圆柱体体积V1

(h,α,β)

0xzyD1D2D3αh将罐体分为圆柱体和球缺，球缺部分单独考虑

精确计算球缺部分体积V2

将倾斜角α的油平面近似为水平面（α=0）计算球缺部分体积V2.直接将罐体的球缺部分按照体积不变条件折合成圆柱，一起计算.用各种方法得到罐体体积V(h,α,β),一般是积分表达式,给定α,β,可以数值计算油位高度h时的储油量V.问题（2）的通常做法2.根据储油量的数学模型V(h,α,β)和实测数据（附件2）辨识罐体变位参数α,β流水号C进油量/LD出油量/LE显示油高/mm显示油量容积/L采集时间备注说明201060.002632.2360448.882010-08-0108:00:49←开始采集（i=0）2020149.092624.3060311.432010-08-0108:15:42i=1hiViuiVi=Vi(hi,0,0)不能用！应该用数据ui,hi与模型V(h,α,β)

辨识参数α,β确定α,β应使储油量的改变ΔVi

与出油量ui相吻合！辨识准则确定α,β使Q最小

搜索法：按照α,β的可能范围划分区间，逐步搜索.

直接利用MATLAB的非线性曲线拟合程序lsqcurvefit辨识方法辨识结果ui,hi取进油前的300组数据h/cm20406080100V/L1065.803702.657371.3811756.6116664.62h/cm120140160180200V/L21941.1827450.7733066.9938667.2744128.48h/cm220240260280300V/L49322.4454109.9358329.2761768.9064026.17=2.110,β=4.310的修正罐容表(部分)罐体变位后显示储油量与实际储油量的比较模型及结果检验

取进油后的300组数据检验辨识结果α,β

α,β的灵敏性检验α的灵敏度β的灵敏度纵向倾斜角α比横向偏转角

β对罐容量的影响显著得多！一种建立模型V(h,α,β)的新方法xy0zhαhYZX0绕x轴转动α绕z轴转动βT2T1=T1T2?αhYZX0VVαxy0zhLdh-r油浮子坐标(0,h-r,L/2-d)过油浮子的油平面方程坐标变换xy0zhrsRVα油位高度h的储油区域为D油位高度h时的储油量V定义示性函数利用MATLAB的三重积分程序计算储油量triplequad(@(x,y,z)(-x*cos(a)*sin(b)+y*cos(a)*cos(b)+z*sin(a)<=cos(a