《数学电工数模》PPT课件.ppt
120页
电工数学建模竞赛 一．总则 n全国大学生电工数学建模竞赛（以下简称竞赛）是中国电机工程学会电工数学专委会主办的面向全国大学生的科技活动，目的是提高学生的综合素质、增强创新意识、培养学生应用数学知识解决实际工程问题的能力，激发学生学习数学的积极性，同时也将推动高校的教学改革与教育创新的进程 二．竞赛内容 n竞赛题目一般来源于电工、近代数学及经济管理等方面，经过适当的简化、加工的实际问题，主要包括：1．信息处理问题；2．控制理论及应用问题；3．运筹与决策问题；4．电路与电磁场理论相关问题 二．竞赛内容 n参赛学生应学过普通高校的工科数学课程及相关专业的专业基础知识，不要求参赛者预先掌握深入的专门知识竞赛题目比较灵活，能够使参赛学生充分发挥其创造能力参赛者要根据题目要求，完成一篇包括模型的假设、建立和求解、算法的设计和计算机实现、结果的分析和检验、模型的改进等方面的论文（即答卷）答卷的评定以假设的合理性、建模的创造性、结果的可行性和文字的清晰度为主要标准 三．竞赛形式、规则和纪律 n1．全国统一竞赛题目，采取网上竞赛方式，以学校为单位进行2．竞赛一般在每年11月末的三天内举行。
3．本科生以队为单位参赛，每队3人，专业不限每队可设一名指导教师（或教师组），对参赛学生进行赛前的及赛前准备工作竞赛期间指导教师要回避参赛队员，禁止进行指导或参与讨论 三．竞赛形式、规则和纪律 n4．竞赛期间参赛学生可以使用各种图书资料、计算机和软件以及在网上浏览，不能与队外任何人讨论5．竞赛题目将按照规定时间准时在指定的网站公布，参赛队员在规定的时间内完成答卷，并准时在网上交卷6．各参赛院校应责成有关职能部门负责竞赛的组织和纪律监督工作，保证竞赛的规范性和公正性 数学建模 Ø 数学模型数学模型（（Mathematical Model）） 是用数学符号、数学式子、程序、图形等对实际课是用数学符号、数学式子、程序、图形等对实际课题本质属性的抽象而又简洁的刻划，它或解释某些题本质属性的抽象而又简洁的刻划，它或解释某些客观现象，或能预测未来的发展规律，或能为控制客观现象，或能预测未来的发展规律，或能为控制某一现象的发展提供某种意义下的最优策略或较好某一现象的发展提供某种意义下的最优策略或较好策略 Ø 数学建模数学建模（（Mathematical Modeling）） 应用知识从实际课题中抽象、提炼出数学模型的过应用知识从实际课题中抽象、提炼出数学模型的过程。
程§1.1 数学模型与数学建模数学模型与数学建模 Ø 1.了解问题的实际背景，明确建模目的，收集掌握了解问题的实际背景，明确建模目的，收集掌握必要的数据资料必要的数据资料Ø 2. 通过对资料的分析计通过对资料的分析计 算，算， 找出起主要作用的因找出起主要作用的因素，经必要的精炼、简化，提出若干符合客观实际的素，经必要的精炼、简化，提出若干符合客观实际的假设Ø 3.在所作假设的基础上，利用适当的数学工具去刻在所作假设的基础上，利用适当的数学工具去刻划各变量之间的关系，建立相应的数学结构划各变量之间的关系，建立相应的数学结构 ——即即建立数学模型建立数学模型 Ø 4.模型求解模型求解 Ø 5.模型的分析与检验模型的分析与检验 在难以得出解析解时，也在难以得出解析解时，也应当借助应当借助 计算机计算机 求出数值求出数值解 §§1.21.2 数学建模的一般步骤数学建模的一般步骤实体信实体信息息(数据数据)假设假设建模建模求解求解验证验证应用应用 §1.3 数学模型的分类数学模型的分类分类标准分类标准分类标准分类标准具体类别具体类别具体类别具体类别对某个实际问题对某个实际问题了解的深入程度了解的深入程度白箱模型、灰箱模型、黑箱模型白箱模型、灰箱模型、黑箱模型模型中变量的特模型中变量的特征征连续型模型、离散型模型或确定性连续型模型、离散型模型或确定性模型、随机型模型等模型、随机型模型等建模中所用的数建模中所用的数学方法学方法初等模型、微分方程模型、差分方初等模型、微分方程模型、差分方程模型、优化模型等程模型、优化模型等研究课题的实际研究课题的实际范畴范畴人口模型、生态系统模型人口模型、生态系统模型 、交通、交通流模型、经济模型、流模型、经济模型、 基因模型等基因模型等 ①①数学建模实践的数学建模实践的 每一步中都每一步中都 蕴含着能力上的蕴含着能力上的 锻炼，在锻炼，在调查研究阶段，需调查研究阶段，需 要用到要用到观察能力观察能力、、分析能力分析能力和和数据处理数据处理能力能力等。
在提出假设等在提出假设 时，又需要用到时，又需要用到 想象力和归纳想象力和归纳 简化简化能力 ②②在真正开始自己的研究之前，还应当尽可能先了解一下在真正开始自己的研究之前，还应当尽可能先了解一下前人或别人的工作，使自己的工前人或别人的工作，使自己的工 作成为别人研究工作作成为别人研究工作 的的继续而不是别人工作的重复，你可以把某些已知的研究结继续而不是别人工作的重复，你可以把某些已知的研究结果用作你的假设，去探索新的奥秘因此我们还应当学会果用作你的假设，去探索新的奥秘因此我们还应当学会在尽可能短的时间在尽可能短的时间 内内查到并学会查到并学会我想应用的知识的本领我想应用的知识的本领 ③③还需要你多少要有点还需要你多少要有点 创新的能力创新的能力这种能力不是生来就这种能力不是生来就有的，建模实践就为你提供了一个培养创新能力的机会有的，建模实践就为你提供了一个培养创新能力的机会 §1.4 数学建模与能力的培养数学建模与能力的培养 开设数学建模课的主要目的为了提高学开设数学建模课的主要目的为了提高学 生的生的综合素质综合素质，增强，增强 应用数学知识应用数学知识 解决实际问解决实际问 题的本领。
题的本领 •例例1 某人平时下班总是按预定时间到达某处，然某人平时下班总是按预定时间到达某处，然然后他妻子开车接他回家有一天，他比平时提早然后他妻子开车接他回家有一天，他比平时提早了三十分钟到达该处，于是此人就沿着妻子来接他了三十分钟到达该处，于是此人就沿着妻子来接他的方向步行回去并在途中遇到了妻子，这一天，他的方向步行回去并在途中遇到了妻子，这一天，他比平时提前了十分钟到家，问此人共步行了多长时比平时提前了十分钟到家，问此人共步行了多长时间？间？ §§1.5 1.5 一些简单实例一些简单实例 似乎条件不够哦似乎条件不够哦 换一种想法，问题就迎刃而换一种想法，问题就迎刃而解了假如他的妻子遇到他后仍解了假如他的妻子遇到他后仍载着他开往会合地点，那么这一载着他开往会合地点，那么这一天他就不会提前回家了提前的天他就不会提前回家了提前的十分钟时间从何而来？十分钟时间从何而来？ 显然是由于节省了从相遇点到显然是由于节省了从相遇点到会合点，又从会合点返回相遇点这一会合点，又从会合点返回相遇点这一段路的缘故，故由相遇点到会合点需段路的缘故，故由相遇点到会合点需开开5分钟。
而此人提前了三十分钟到分钟而此人提前了三十分钟到达会合点，故相遇时他已步行了二十达会合点，故相遇时他已步行了二十五分钟 请思考一下，本题解答中隐含了哪些假设请思考一下，本题解答中隐含了哪些假设请思考一下，本题解答中隐含了哪些假设请思考一下，本题解答中隐含了哪些假设 ？？ •例例2 2 交通灯在绿灯转换成红灯时，有一交通灯在绿灯转换成红灯时，有一个过渡状态个过渡状态————亮一段时间的黄灯请分亮一段时间的黄灯请分析黄灯应当亮多久析黄灯应当亮多久设想一下黄灯的作用是什么，不难看设想一下黄灯的作用是什么，不难看出，黄灯起的是警告的作用，意思是出，黄灯起的是警告的作用，意思是马上要转红灯了，假如你能停住，请马上要转红灯了，假如你能停住，请立即停车停车是需要时间的，在这立即停车停车是需要时间的，在这段时间内，车辆仍将向前行驶一段距段时间内，车辆仍将向前行驶一段距离离 L这就是说，在离街口距离为这就是说，在离街口距离为 L处处存在着一条停车线（尽管它没被画在存在着一条停车线（尽管它没被画在地上），见图对于那些黄灯亮时已地上），见图对于那些黄灯亮时已过线的车辆，则应当保证它们仍能穿过线的车辆，则应当保证它们仍能穿过马路。
过马路 马路的宽度马路的宽度 D是容易测得是容易测得 的，问题的关键在的，问题的关键在 于于L的确定为确定的确定为确定 L，还应当将，还应当将 L划分为两段：划分为两段：L1和和L2，，其中其中 L1是司机在发现黄灯亮及判断应当刹是司机在发现黄灯亮及判断应当刹车的反应时间内驶过的路程车的反应时间内驶过的路程 ，，L2为刹车制动后为刹车制动后车辆驶过的路程车辆驶过的路程L1较容易计算，交通部门对司较容易计算，交通部门对司机的平均反应时间机的平均反应时间 t1早有测算，反应时间过长早有测算，反应时间过长将考不出驾照），而此街道的行驶速度将考不出驾照），而此街道的行驶速度 v 也是也是交管部门早已定好的，目的是使交通流量最大，交管部门早已定好的，目的是使交通流量最大，可另建模型研究，从而可另建模型研究，从而 L1=v*t1刹车距离刹车距离 L2既可用曲线拟合方法得出，也可利用牛顿第二定既可用曲线拟合方法得出，也可利用牛顿第二定律计算出来律计算出来 黄灯究竟应当亮多久现在已经变得清楚多了第黄灯究竟应当亮多久现在已经变得清楚多了第一步，先计算出一步，先计算出 L应多大才能使看见黄灯的司机应多大才能使看见黄灯的司机停得住车。
第二步，黄灯亮的时间应当让已过线停得住车第二步，黄灯亮的时间应当让已过线的车顺利穿过马路，即的车顺利穿过马路，即T 至少应当达到至少应当达到 （（L+D））/v DL 初等模型初等模型初等模型初等模型 某航空母舰派其护卫舰去搜寻其跳伞的飞行某航空母舰派其护卫舰去搜寻其跳伞的飞行员，护卫舰找到飞行员后，航母通知它尽快员，护卫舰找到飞行员后，航母通知它尽快 返回与其汇合并通报了航母当前的航速与方返回与其汇合并通报了航母当前的航速与方 向，问护卫舰应怎样航行，才能与航母汇合向，问护卫舰应怎样航行，才能与航母汇合§2.1 舰艇舰艇的会合的会合 令：令：则上式可简记成则上式可简记成 ：：A(0,b)XYB(0,-b)P(x,y)O航母航母 护卫舰护卫舰 θ1 θ2 即：即：可化为：可化为：记记v2/ v1=a通常通常a>1 则则汇合点汇合点 p必位于此圆上必位于此圆上 （护卫舰的路线方程）（护卫舰的路线方程）（航母的路线方程（航母的路线方程 ））即可求出即可求出P点的坐标和点的坐标和θ2 的值本模型虽简单，但分析本模型虽简单，但分析极极清晰清晰且且易于实际应用易于实际应用 §2.2 双层玻璃的功效双层玻璃的功效在寒冷的北方，在寒冷的北方， 许多住房的玻璃窗都是双层许多住房的玻璃窗都是双层玻璃的，现在我们来建立一个简单玻璃的，现在我们来建立一个简单 的数学模的数学模型，研究一下双层玻璃到底有多大的功效。
型，研究一下双层玻璃到底有多大的功效比较两座其他条件完全相同的房屋，它们的比较两座其他条件完全相同的房屋，它们的差异仅仅在窗户不同差异仅仅在窗户不同 不妨可以提出以下不妨可以提出以下 假设假设：：1、设室内热量的流失是热传导、设室内热量的流失是热传导引起的，不存在户内外的空气对引起的，不存在户内外的空气对流2、室内温、室内温 度度T1与户外温与户外温 度度T2均均为常数3、玻璃是均匀的，热传导系数、玻璃是均匀的，热传导系数为常数 设玻璃的热传导系数设玻璃的热传导系数 为为k1，空气的，空气的热传导系数热传导系数 为为k2，单位时间通过单，单位时间通过单位面积由温度高的一侧流向温度低位面积由温度高的一侧流向温度低的一侧的热量为的一侧的热量为θ ddl室室外外T2室室内内T1TaTb由热传导公式由热传导公式 θ=kΔT/d 解得：解得： 此函数的图形为此函数的图形为dd室室外外T2室室内内T1类似有类似有 一般一般故故记记h=l/d并令并令f(h)= 01234567891000.10.20.30.40.50.60.70.80.91hf(h)考虑到考虑到美观美观和使用上和使用上 的的方便方便，，h不必取得过大，例如，可不必取得过大，例如，可 取取h=3，即，即l=3d，此时房屋热量的损失不超过单层玻璃窗，此时房屋热量的损失不超过单层玻璃窗时的时的 3% 。
§2.3 崖高的估算崖高的估算假如你站在崖顶且身上带着一只具有跑表功假如你站在崖顶且身上带着一只具有跑表功 能的计算器，你也许会出于好奇心想用扔下能的计算器，你也许会出于好奇心想用扔下 一块石头听回声的方法来估计山崖的高度，一块石头听回声的方法来估计山崖的高度， 假定你能准确地测定时间，你又怎样来推算假定你能准确地测定时间，你又怎样来推算 山崖的高度呢，请你分析一下这一问题山崖的高度呢，请你分析一下这一问题我有一只具有跑我有一只具有跑 表功能的计算器表功能的计算器 方法一方法一假定空气阻力不计，可以直接利用自由落体运动的公式假定空气阻力不计，可以直接利用自由落体运动的公式来计算例如，来计算例如， 设设t=4秒，秒，g=9.81米米/秒秒2，则可求得，则可求得h≈78.5米 我学过微积分，我可以做我学过微积分，我可以做 得更好，呵呵得更好，呵呵 除去地球吸引力外，对石块下落影响最大的当除去地球吸引力外，对石块下落影响最大的当 属属空气阻力空气阻力根据流体力学知识，此时可设空气阻力正比于石块下落的根据流体力学知识，此时可设空气阻力正比于石块下落的速度，阻力系速度，阻力系 数数K为常数，因而，由牛顿第二定律可得：为常数，因而，由牛顿第二定律可得： 令令k=K/m，，解得解得 代入初始条件代入初始条件 v(0)=0，得，得c=－－g/k，故有，故有 再积分一次，得：再积分一次，得： 若设若设k=0.05并仍设并仍设 t=4秒，则可求秒，则可求 得得h≈73.6米。
米 听到回声再按跑表，计算得到的时间中包含了听到回声再按跑表，计算得到的时间中包含了 反应时间反应时间 进一步深入考虑进一步深入考虑进一步深入考虑进一步深入考虑不妨设不妨设平均反应时间平均反应时间 为为0.1秒秒 ，假如仍，假如仍 设设t=4秒，扣除反秒，扣除反应时间后应应时间后应 为为3.9秒，代入秒，代入 式式①①，求得，求得h≈69.9米 ①①多测几次，取平均值多测几次，取平均值再一步深入考虑再一步深入考虑再一步深入考虑再一步深入考虑代入初始条代入初始条 件件h(0)=0，得到计算山崖高度的公式：，得到计算山崖高度的公式： 将将e-kt用泰勒公式展开并用泰勒公式展开并 令令k→ 0+ ，即可，即可得出前面不考虑空气阻力时的结果得出前面不考虑空气阻力时的结果 还应考虑还应考虑回声回声传回来所需要的时间为此，令石块下落传回来所需要的时间为此，令石块下落 的真正时间的真正时间 为为t1，声音传回来的时间记，声音传回来的时间记 为为t2，还得解一个，还得解一个方程组：方程组： 这一方程组是这一方程组是非线性非线性的，求的，求解不太容易，解不太容易，为了估算崖高为了估算崖高竟要去解一个竟要去解一个非线性主程组非线性主程组似乎不合情理似乎不合情理 相对于石块速度，声音速度要快得多，我们可相对于石块速度，声音速度要快得多，我们可 用方法二先求一次用方法二先求一次 h，令，令t2=h/340，校正，校正t，求石，求石块下落时间块下落时间 t1≈t-t2将将t1代入式代入式①①再算一次，得出再算一次，得出崖高的近似值。
例如，崖高的近似值例如， 若若h=69.9米，则米，则 t2≈0.21秒，故秒，故 t1≈3.69秒，求得秒，求得 h≈62.3米 最小二乘法最小二乘法 插值方法插值方法 当问题的机理非常不清楚难以直接利用其他知当问题的机理非常不清楚难以直接利用其他知识来建模时，一个较为自然的方法是利用数据识来建模时，一个较为自然的方法是利用数据进行曲线拟合，找出变量之间的近似依赖关系进行曲线拟合，找出变量之间的近似依赖关系即函数关系即函数关系§2.4 经验模型经验模型 最小二乘法最小二乘法设经实际测量已得设经实际测量已得 到到n组数据（组数据（xi , yi），），i=1,…, n将数据画在平面直角坐标系中，见画在平面直角坐标系中，见 图如果建模者判断图如果建模者判断 这这n个点很个点很象是分布在某条直线附近，令象是分布在某条直线附近，令 该直线方程该直线方程 为为y=ax+b，进而，进而利用数据来求参利用数据来求参 数数a和和b由于该直线只是数据近似满足的关由于该直线只是数据近似满足的关系式，故系式，故 yi-(axi+b)=0一般不成立，但我们希望一般不成立，但我们希望 最小最小此式对此式对a和和b的偏导数均的偏导数均 为为0，，解相应方程组，求得：解相应方程组，求得： y=ax+byO(xi ,yi)x其中其中 和和 分别为分别为xi和和yi的平均值的平均值 如果建模者判断变量间的关系并非线性关系而是其他类型的函数，如果建模者判断变量间的关系并非线性关系而是其他类型的函数，则可作则可作 变量替换变量替换使之转化为线性关系或用类似方法使之转化为线性关系或用类似方法拟合拟合。
显然，运动员体重越大，他能举起的重量也越大，但举重显然，运动员体重越大，他能举起的重量也越大，但举重成绩和运动员体重到底是怎样关系的，不同量级运动员的成绩和运动员体重到底是怎样关系的，不同量级运动员的成绩又如何比较优劣呢？运动成绩是包括生理条件、心理成绩又如何比较优劣呢？运动成绩是包括生理条件、心理因素等等众多相关因素共同作用的结果，要建立精确的模因素等等众多相关因素共同作用的结果，要建立精确的模型至少现在还无法办到但我们拥有大量的比赛成绩纪录，型至少现在还无法办到但我们拥有大量的比赛成绩纪录，根据这些数据不妨可以建立一些经验模型为简单起见，根据这些数据不妨可以建立一些经验模型为简单起见，我们不妨取表中的数据为例我们不妨取表中的数据为例例例1（举重成绩的比较）（举重成绩的比较）举重举重是一种一般人都能看懂的运动，它共分是一种一般人都能看懂的运动，它共分九个重量级，有两种主要的比赛方法：抓举九个重量级，有两种主要的比赛方法：抓举和挺举 表中给出了到表中给出了到1977年底为止九个年底为止九个重量级的世界纪录重量级的世界纪录255200110以上以上237.518511022118090207.517082.5195157.575180141.567.5161.513060151120.55614110952挺挺举举（公斤）（公斤）抓抓举举（公斤）（公斤）成成绩绩重量重量级级（上限体（上限体重）重） 模型模型1（线性模型）（线性模型） 将数据画在直角坐标系中可以发现，运动成绩与体将数据画在直角坐标系中可以发现，运动成绩与体量近似满足线性关系，只有量近似满足线性关系，只有110公斤级有点例外，两公斤级有点例外，两项成绩都显得较低。
应用前面叙述的方法可求出近项成绩都显得较低应用前面叙述的方法可求出近似关似关 系式系式L=kB+C，其中，其中B为体重，为体重，L为举重成绩为举重成绩你在作图你在作图 时时L轴可以放轴可以放 在在50公斤或公斤或52公斤处，因为公斤处，因为没有更轻级别的比赛，具体计算留给同学自己去完没有更轻级别的比赛，具体计算留给同学自己去完成 模型模型2（幂函数模型）（幂函数模型） 线性模型并未得到广泛的接受，要改进结果，能够线性模型并未得到广泛的接受，要改进结果，能够想到的自然首先是幂函数模型，即令想到的自然首先是幂函数模型，即令L=kBa，对此式，对此式取对数，得取对数，得 到到lnL=lnk+a lnB将原始数据也取对数，将原始数据也取对数，问题即转化了线性模型，可用最小二乘法求出参数问题即转化了线性模型，可用最小二乘法求出参数几十年前英国和爱尔兰采用的比较举重成绩优劣几十年前英国和爱尔兰采用的比较举重成绩优劣 的的Austin公式公式：：L′=L/B3/4就是用这一方法求得的就是用这一方法求得的 模型模型3（经典模型）（经典模型） 经典模型是根据生理学中的已知结果和比例关系推导出来的经典模型是根据生理学中的已知结果和比例关系推导出来的公式，应当说，它并不属于经验公式。
为建立数学模型，先公式，应当说，它并不属于经验公式为建立数学模型，先提出如下一些假设：提出如下一些假设： (1)举重成绩正比于选手肌肉的平均横截举重成绩正比于选手肌肉的平均横截 面积面积A，即，即L=k1A(2)A正比于身高正比于身高 L的平方，即的平方，即 A=k2L2(3)体重正比于身高体重正比于身高 L的三次方，的三次方， 即即B=k3L3根据上述假设，可得根据上述假设，可得 显然，显然，K越大则成绩越好，故可用越大则成绩越好，故可用 来比较选手来比较选手比赛成绩的优劣比赛成绩的优劣 模型模型4（（O’ Carroll公式）公式） 经验公式的主要依据是比例关系，其假设条件非常粗糙，可经验公式的主要依据是比例关系，其假设条件非常粗糙，可信度不大，因而大多数人认为它不能令人信服信度不大，因而大多数人认为它不能令人信服1967年，年，O’ Carroll基于动物学和统计分析得出了一个现在被广泛使用的基于动物学和统计分析得出了一个现在被广泛使用的公式O’ Carroll模型的假设条件是：模型的假设条件是： (1) L=k1Aa, a<1 (2) A=k2Lb, b<2 (3) B-Bo =k3L3 假设假设(1)、、(2)是解剖学中的统计规律，在假设是解剖学中的统计规律，在假设 （（3）中）中O’ Carroll将体重划分成两部分：将体重划分成两部分：B=B0+B1，B0为非肌肉重量。
为非肌肉重量 故有：故有： 根据三条假设可根据三条假设可 得得L=k(B-B0)β,k和和β为两个常数，为两个常数， 此外，根据统计结果，他此外，根据统计结果，他 得出得出B0≈35公斤，公斤， k越大成绩越好因而越大成绩越好因而 建议根据建议根据 来比来比 较选手成绩的优劣较选手成绩的优劣 模型模型5（（Vorobyev公式）公式） 这是一个前苏联使用的公式建模者认为举重选手举起的不这是一个前苏联使用的公式建模者认为举重选手举起的不光是重物，也提高了自己的重心，故其举起的总重量为光是重物，也提高了自己的重心，故其举起的总重量为L+B，可以看出，他们更重视的是腿部肌肉的爆发力应用与模，可以看出，他们更重视的是腿部肌肉的爆发力应用与模型型4类似的方法，得出了按类似的方法，得出了按 的大小比较成绩优劣的建议的大小比较成绩优劣的建议 上述公式具有各不相同的基准，无法相互比较为了使公式具上述公式具有各不相同的基准，无法相互比较为了使公式具有可比性，需要对公式稍作处理。
例如，我们可以要求各公式有可比性，需要对公式稍作处理例如，我们可以要求各公式均满足均满足在在 B=75公斤时有公斤时有 L’=L，则上述各公式化为：，则上述各公式化为：（（1））Austin公式：公式：（（2））经典经典公式：公式：（（3））O’ Carroll公式：公式：（（4））Vorobyev公式：公式： 将公式（将公式（1））—（（4）用来比较）用来比较1976年奥运会的抓举成绩，各年奥运会的抓举成绩，各公式对九个级别冠军成绩的优劣排序如表公式对九个级别冠军成绩的优劣排序如表 所示，比较结果较所示，比较结果较为一致，例如，对前三名的取法是完全一致的，其他排序的为一致，例如，对前三名的取法是完全一致的，其他排序的差异也较为微小差异也较为微小 138.5(8)141.9(7)135.6(7)131.8(8)175110150.3(2)152.9(2)150.5(2)148.3(2)17090152.1(1)153.5(1)152.2(1)151.3(1)162.542.5145.0(6)145.0(5)145.0(3)145.0(6)14575145.8(5)144.7(6)144.8(5)146.1(5)13567.5147.7(3)146.2(3)145.0(3)147.8(3)12560146.6(4)145.7(4)142.8(6)146.3(4)117.556138.8(7)139.7(8)134.0(8)138.2(7)10552VorobyevO’ Carroll经典公式经典公式Austin抓举成绩抓举成绩(公斤公斤)体重体重(公斤公斤) 我们希望建立一个我们希望建立一个 体重体重与与身高身高之间的关系式，不难看出两者之间的关系式，不难看出两者之间的关系不易通过机理的分析得出，不妨可以采取之间的关系不易通过机理的分析得出，不妨可以采取 统计统计方法方法，用数据来拟合出与实际情况较为相符的经验公式。
用数据来拟合出与实际情况较为相符的经验公式 为为此，我们先作一番抽样调查，测量了十五个不同高度的人的此，我们先作一番抽样调查，测量了十五个不同高度的人的体重，列成了体重，列成了 下表，在抽样时，各高度的人都需经适当挑选，下表，在抽样时，各高度的人都需经适当挑选，既不要太胖也不要太瘦既不要太胖也不要太瘦例例2 体重与身高的体重与身高的 关系关系将表中的数画将表中的数画 到到h-w平面上，你会发现这些数据分布很接近某平面上，你会发现这些数据分布很接近某一指数曲线为此，一指数曲线为此， 对对h和和w均取对数，令均取对数，令x=lnh，，y=lnw，将，将（（xi，yi）再画到）再画到x-y平面中去（平面中去（i=1,…,15），这次你会发现这），这次你会发现这些点几乎就分布在一条直线附近，令此直线的些点几乎就分布在一条直线附近，令此直线的 方程为方程为y=ax+b，用最小二乘法求，用最小二乘法求 得得a≈2.3,b≈2.82，故可取，故可取y=2.32x+2.84，即，即lnw=2.32lnh+2.84，故有，故有w=17.1h2.327566595451体重体重 w（公斤）（公斤）1.851.781.711.671.63身高身高 h（米）（米）5048413527体重体重 w（公斤）（公斤）1.601.551.511.351.26身高身高 h（米）（米）2017151210体重体重 w（公斤）（公斤）1.121.080.960.860.75身高身高 h（米）（米） 在使用在使用 最小二乘法最小二乘法 时，我们并未要求得到的拟合曲线一定时，我们并未要求得到的拟合曲线一定要经过所有的样本点，而只是要求要经过所有的样本点，而只是要求 了了总偏差最小总偏差最小。
当实际当实际问题要求拟合曲线必须问题要求拟合曲线必须 经过样本点经过样本点 时，我们可以应用数值时，我们可以应用数值逼近中的逼近中的 插值法插值法根据实际问题的不同要求，存在多种不同的插值方法，有根据实际问题的不同要求，存在多种不同的插值方法，有只要求过样本点的只要求过样本点的 拉格朗日插值拉格朗日插值 法法、、牛顿插值法牛顿插值法 等，有既等，有既要求过插值点（即样本点）又对插值点处的导数有所要求要求过插值点（即样本点）又对插值点处的导数有所要求的的样条（样条（Spline）插值）插值，甚至还有对插值曲线的凹凸也有，甚至还有对插值曲线的凹凸也有要求的要求的B样条插值法样条插值法 本课不准备详细介绍这些细致的插本课不准备详细介绍这些细致的插值方法，只是提请读者注意，在建立经验模型时，插值法值方法，只是提请读者注意，在建立经验模型时，插值法也是可以使用的数学工具之一也是可以使用的数学工具之一 插值方法插值方法 对插值法感兴趣的对插值法感兴趣的 同学可以查同学可以查阅相关书籍，例如由阅相关书籍，例如由 李岳生李岳生编著上编著上海科学技术出版社出版的海科学技术出版社出版的《《样条与样条与插值插值》》（（1983年出版）等。
年出版）等 在建立数学模型时常常需要确定一些参数，选什么量为参数，在建立数学模型时常常需要确定一些参数，选什么量为参数，怎样选取参数，其中也有一些技巧，参数选得不好，会使问怎样选取参数，其中也有一些技巧，参数选得不好，会使问题变得复杂难解，给自己增添许多不必要的麻烦确定参数题变得复杂难解，给自己增添许多不必要的麻烦确定参数以后，一般需要利用数据来获得这些参数的具体取值，例如以后，一般需要利用数据来获得这些参数的具体取值，例如在使用经验方法建模时，假如你准备用线性函在使用经验方法建模时，假如你准备用线性函 数数ax+b来表来表达变量间的关系，你还要用最小二乘法去求出参达变量间的关系，你还要用最小二乘法去求出参 数数a、、b的的值，这一过程被称值，这一过程被称 为为“参数识参数识 别别”总之，参数的选取应使总之，参数的选取应使其后的识别尽可能简便，让我们来考察一个实例其后的识别尽可能简便，让我们来考察一个实例§2.5 参数识别参数识别 例例3 录像带还能录多长时间录像带还能录多长时间录像机上有一个四位计数器，一盘录像机上有一个四位计数器，一盘 180分钟分钟的录像带在开始计数时为的录像带在开始计数时为 0000，到结束时计，到结束时计数为数为1849，实际走时为，实际走时为185分分20秒。
我们从秒我们从0084观察到观察到0147共用时间共用时间3分分21秒若录像秒若录像机目前的计数为机目前的计数为1428，问是否还能录下一个，问是否还能录下一个 60分钟的节目？分钟的节目？ rθRl由由得到得到又又 因和因和 得得 积分得到积分得到即即从而有从而有我们希望建立一个录像带已录像时我们希望建立一个录像带已录像时 间间t与计数器计与计数器计 数数n之间的函数关系为建立一个正确的模型，首之间的函数关系为建立一个正确的模型，首 先必先必须搞清哪些量是常量，哪些量是变量首先，录像须搞清哪些量是常量，哪些量是变量首先，录像 带带的厚的厚 度度W是常量，它被绕在一个半径是常量，它被绕在一个半径 为为r的园盘上，的园盘上，见图磁带转动中线速见图磁带转动中线速 度度v显然也是常数，否则图象显然也是常数，否则图象声音必然会失真此外，计数器的读声音必然会失真此外，计数器的读 数数n与转过的圈与转过的圈数有关，从而与转过的角数有关，从而与转过的角 度度θ成正比 rθRl 此式中的三个参数此式中的三个参数ω、、v和和r均不易精确测得，均不易精确测得，虽然我们可以从上式解出虽然我们可以从上式解出t与与n的函数关系，的函数关系，但效果不佳，故令但效果不佳，故令 则可将上式简化为：则可将上式简化为： 故故令令上式又可化简记成上式又可化简记成 t= an2+bn t= an2+bn rθRl上式以上式以a、b为参数显然是一个十分明智的为参数显然是一个十分明智的做法，它为公式的最终确立即参数求解提做法，它为公式的最终确立即参数求解提供了方便。
将已知条件代入，得方程组：供了方便将已知条件代入，得方程组： 从后两式中消从后两式中消 去去t1，解得，解得a=0.0000291, b=0.04646,故故t=0.0000291 n2+0.04646n，令，令n=1428，得到，得到t=125.69（分）由于一盒录像带实际可录像时间为（分）由于一盒录像带实际可录像时间为185.33分，分，故尚可录像时间故尚可录像时间 为为59.64分，已不能再录下一个分，已不能再录下一个60分钟分钟的节目了的节目了 物理量大都带有量纲，其中基本量纲通常是质量（用物理量大都带有量纲，其中基本量纲通常是质量（用M表示）表示）、长度（、长度（ 用用L表示）、时间（表示）、时间（ 用用T表示），有时还有温度表示），有时还有温度（用（用Θ表示）其他物理量的量纲可以用这些基本量纲来表表示）其他物理量的量纲可以用这些基本量纲来表示，如速度的量纲为示，如速度的量纲为LT-1，加速度的量纲为，加速度的量纲为 LT-2，力的量纲，力的量纲为为 MLT-2，功的量纲为，功的量纲为 ML2T-2等 §2.6 量纲分析法建模量纲分析法建模量纲分析的原理量纲分析的原理 是：当度量量纲的基本单位改变时，公式是：当度量量纲的基本单位改变时，公式本身并不改变，例如，无论长度取什么单位，矩形的面积总本身并不改变，例如，无论长度取什么单位，矩形的面积总等于长乘宽，即公式等于长乘宽，即公式 S=ab并不改变。
此外，在公式中只有并不改变此外，在公式中只有量纲相同的量才能进行加减运算，例如面积与长度是不允许量纲相同的量才能进行加减运算，例如面积与长度是不允许作加减运算的，这些限制在一定程度上限定了公式的可取范作加减运算的，这些限制在一定程度上限定了公式的可取范围，即一切公式都要求其所有的项具有相同的量纲，具有这围，即一切公式都要求其所有的项具有相同的量纲，具有这种性质的公式被称为种性质的公式被称为 是是“量纲齐次量纲齐次”的 例例3 在万有引力公式中，引力常数在万有引力公式中，引力常数G是有量纲的，根据量是有量纲的，根据量纲齐次性，纲齐次性，G的量纲为的量纲为M-1L3T-2，其实，在一量纲齐次的公，其实，在一量纲齐次的公式中除以其任何一项，即可使其任何一项化为无量纲，因式中除以其任何一项，即可使其任何一项化为无量纲，因此任一公式均可改写成其相关量的无量纲常数或无量纲变此任一公式均可改写成其相关量的无量纲常数或无量纲变量的函数例如，与万有引力公式量的函数例如，与万有引力公式 相关的物理量有：相关的物理量有：G、m1、m2、r和和F。
现考察这些量的无量纲乘积现考察这些量的无量纲乘积 的量纲为的量纲为由于由于 是无量纲的量，故应有：是无量纲的量，故应有： 此方程组中存在两个自由变量，其解构成一个二维线性空此方程组中存在两个自由变量，其解构成一个二维线性空间取（（a,b））=（（1,0））和和（（a,b））=（（0,1），），得到方程组解得到方程组解空间的一组基空间的一组基 （（1,0,2,-2,-1））和和（（0,1,-1,0,0），），所有由这些所有由这些量组成的无量纲乘积均可用这两个解的线性组合表示两量组成的无量纲乘积均可用这两个解的线性组合表示两个基向量对应的无量纲乘积分别为：个基向量对应的无量纲乘积分别为：而万有引力定律则可写而万有引力定律则可写 成成f(π1,π2)=0，其对应的显函数为：，其对应的显函数为：π1=g(π2)，即即 万有引万有引力定律力定律 定理定理2.1 （（Backinghamπ定理）方程当且仅当可以表定理）方程当且仅当可以表 示为示为 f（π1，π2…）=0时时才是量纲齐次的，其中才是量纲齐次的，其中 f是某是某一函数，一函数，π1，π2…为问题所包含的变量与常数的无量为问题所包含的变量与常数的无量 纲乘积。
纲乘积设此变换的零空间为设此变换的零空间为 m维的，取此零空间的一组基维的，取此零空间的一组基e1,……,em，并将其扩充，并将其扩充 为为k维欧氏空间的一组基维欧氏空间的一组基e1,……,em, em+1,……ek 令令πi=g-1(ei), i=1,…,k,显然，显然，π1，…, πm是无量纲的，而是无量纲的，而πm+1,…, πk是有量纲的（若是有量纲的（若k>m）由于公式量纲齐次当且仅当它可用无量纲的量表示，）由于公式量纲齐次当且仅当它可用无量纲的量表示，故方程当且仅当可写故方程当且仅当可写 成成f(π1，…, πm)=0时才是量纲齐次时才是量纲齐次的，定理证毕的，定理证毕 证证 设设x1,…,xk为方程中出现的变量与常数为方程中出现的变量与常数, ,对这些变量对这些变量与常数的任一乘积与常数的任一乘积 ,令令 函数函数g建立了建立了xi(i=1,…,k)的乘积所组成的空间的乘积所组成的空间 与与k维欧维欧氏空间之间的一个一一对应现设涉及到的基本量纲有氏空间之间的一个一一对应。
现设涉及到的基本量纲有n个个,它们它们 为为y1,…,yn.用这些基本量纲来表达用这些基本量纲来表达 该该xi的乘幂的乘幂,设此乘幂的量纲为设此乘幂的量纲为 令令易见易见dg-1是是k维欧氏空间维欧氏空间 到到n维欧氏空间的一个变换，这维欧氏空间的一个变换，这里的里的g-1为为g的逆变换的逆变换 例例4（理想单摆的摆动周期）（理想单摆的摆动周期）考察质量集中于距支点为考察质量集中于距支点为 l 的质点上的无阻的质点上的无阻尼尼 单摆，（如图），其运动为某周单摆，（如图），其运动为某周 期期 t 的的左右摆动，现希望得到周期左右摆动，现希望得到周期 t 与其他量之间与其他量之间的的 关系θlmg考考察察，，的的量量纲纲为为MaLb+dTc-2b若若无无量量纲纲，，则则有有 量纲分析法虽然简单，但使用时在技巧方面的要求较高，稍量纲分析法虽然简单，但使用时在技巧方面的要求较高，稍一疏忽就会导出荒谬的结果或根本得不出任何有用的结果一疏忽就会导出荒谬的结果或根本得不出任何有用的结果首先，它要求建模者对研究的问题有正确而充分的了解，能首先，它要求建模者对研究的问题有正确而充分的了解，能正确列出与该问题相关的量及相关的基本量纲，容易看出，正确列出与该问题相关的量及相关的基本量纲，容易看出，其后的分析正是通过对这些量的量纲研究而得出的，列多或其后的分析正是通过对这些量的量纲研究而得出的，列多或列少均不可能得出有用的结果。
其次，在为寻找无量纲量而列少均不可能得出有用的结果其次，在为寻找无量纲量而求解齐次线性方程组时，基向量组有无穷多种取法，如何选求解齐次线性方程组时，基向量组有无穷多种取法，如何选取也很重要，此时需依靠经验，并非任取一组基都能得出有取也很重要，此时需依靠经验，并非任取一组基都能得出有用的结果此外，建模者在使用量纲分析法时对结果也不应用的结果此外，建模者在使用量纲分析法时对结果也不应抱有不切实际的过高要求，量纲分析法的基础是公式的量纲抱有不切实际的过高要求，量纲分析法的基础是公式的量纲齐次性，仅凭这一点又怎么可能得出十分深刻的结果，例如，齐次性，仅凭这一点又怎么可能得出十分深刻的结果，例如，公式可能包含某些无量纲常数或无量纲变量，对它们之间的公式可能包含某些无量纲常数或无量纲变量，对它们之间的关系，量纲分析法根本无法加以研究关系，量纲分析法根本无法加以研究 §2.7 赛艇成绩的比较赛艇成绩的比较(比例模型比例模型)八人赛艇比赛和举重比赛一样，分成八人赛艇比赛和举重比赛一样，分成86公斤公斤的重量级和的重量级和 73公斤的轻量级公斤的轻量级1971年，年，T.A.McMahon比较了比较了1964-1970年期间两次年期间两次奥运会和两次世锦赛成绩，发现奥运会和两次世锦赛成绩，发现 86公斤级比公斤级比73公斤级的成绩大约好公斤级的成绩大约好5%，产生这一差异的，产生这一差异的原因何在呢？原因何在呢？ 我们将以我们将以L表示轻量级、以表示轻量级、以H表示重表示重量级，用量级，用S表示赛艇的浸水面积，表示赛艇的浸水面积，v表示赛艇速度，表示赛艇速度，W表示选手体重，表示选手体重，P表示选手的输出功率，表示选手的输出功率，I表示赛程，表示赛程，T表示比赛成绩（时间）。
表示比赛成绩（时间） 考察优秀赛艇选手在比赛中的实际表现可以发现，整个赛程大考察优秀赛艇选手在比赛中的实际表现可以发现，整个赛程大致可以分三个阶段，致可以分三个阶段， 即即初始时刻的加速阶段初始时刻的加速阶段、、中途的匀速阶中途的匀速阶段段和和到达终点的冲刺阶段到达终点的冲刺阶段 由于赛程较长，可以略去前后两由于赛程较长，可以略去前后两段而段而只考虑中间一段只考虑中间一段 ，为此，提出以下建模假设为此，提出以下建模假设1）设赛艇浸水部分的摩擦力是唯一阻力，摩擦力）设赛艇浸水部分的摩擦力是唯一阻力，摩擦力f正比正比 于于Sv2,（见流体力学），空气阻力等其他因素不计见流体力学），空气阻力等其他因素不计2）同一量级的选手有相同的体重）同一量级的选手有相同的体重W，选手的输出功，选手的输出功 率率P正比于正比于W，且效率大体相同且效率大体相同由由假设假设1，，，故，故 竞赛成绩竞赛成绩记比例系数记比例系数 为为k，则有，则有： 故故由由假设假设2，， 故故令令WH=86,WL=73,则有则有由于由于SL略小于略小于SH，故轻量级所化时间比重量级所化时间，故轻量级所化时间比重量级所化时间约约 多多5%左右。
左右 §2.8 方桌问题方桌问题将一张四条腿的方桌放在不平的地面上，不将一张四条腿的方桌放在不平的地面上，不 允许将桌子移到别处，但允许其绕中心旋转允许将桌子移到别处，但允许其绕中心旋转 ，是否总能设法使其四条腿同时落地？，是否总能设法使其四条腿同时落地？ 不附加任何条件，答案不附加任何条件，答案 显然显然 是否定的，是否定的， 因此我们因此我们假设假设 （（1））地面为连续曲面地面为连续曲面 （（2））方桌的四条腿长度相方桌的四条腿长度相同同 （（3））相对于地面的弯曲程相对于地面的弯曲程度而言，方桌的腿是足够长度而言，方桌的腿是足够长的的 （（4））方桌的腿只要有一点方桌的腿只要有一点接触地面就算着地接触地面就算着地总可以使三条腿总可以使三条腿同时着地同时着地 现在，我们来证明：如果上述假设条件成立，那么答案是肯定现在，我们来证明：如果上述假设条件成立，那么答案是肯定的以方桌的中心为坐标原点作直角坐标系如的以方桌的中心为坐标原点作直角坐标系如 图所示，方桌图所示，方桌的四条腿分别在的四条腿分别在A、、B、、C、、D处，处，A、、C的初始位置在的初始位置在x轴上，轴上，而而B、、D则在则在y轴上，当方桌绕中轴上，当方桌绕中 心心0旋转时，对角线旋转时，对角线 AC与与x轴轴的夹角记为的夹角记为θ。
容易看出，当四条腿尚未全部着地时，腿到地面的距离是不确容易看出，当四条腿尚未全部着地时，腿到地面的距离是不确定的为消除这一不确定性，令定的为消除这一不确定性，令 f(θ)为为A、、C离地距离之和，离地距离之和，g(θ)为为B、、D离地距离之和，它们的值离地距离之和，它们的值 由由θ唯一确定由唯一确定由假设假设（（1），），f(θ)、g(θ)均为均为θ的连续函数又的连续函数又 由由假设（假设（3），），三条腿三条腿总能同时着地，总能同时着地， 故故f(θ)g(θ)=0必成立（必成立（ ）不妨设）不妨设f(0)=0,g(0)>0（若（若g(0)也为也为0，则初始时刻已四条腿着地，不必，则初始时刻已四条腿着地，不必再旋转），于是问题归结为：再旋转），于是问题归结为：yxθCDABo已知已知f(θ)、g(θ)均为均为θ的连续函数，的连续函数，f(0)=0,g(0)>0且对任意且对任意θ有有f(θ)g(θ)=0，求证存在某一，求证存在某一θ0，使，使f(θ0)=g(θ0)=0 （证法一）（证法一）当当θ=π/2时，时，AC与与BD互换位置，故互换位置，故f(π/2)>0 , g(π/2)=0。
作作h(θ)=f(θ)-g(θ)，显然，显然，h(θ)也是也是θ的连续函数，的连续函数，h(0)=f(0)-g(0)<0而而h(π/2)=f(π/2)-g(π/2)>0，由连续函数的取，由连续函数的取零值定理，存在零值定理，存在 θo，0<θo <π/2，h(θ0)=0，即，即f(θo)=g(θo)又由于又由于f(θo)g(θo)=0，故必有，故必有f(θo)=g(θo)=0，证毕 （证法二）（证法二）同证一可得同证一可得 f(π/2)>0,g(π/2)=0令θo =sup {ζ |f (ζ)=0，0≤ζ<θ},显然显然θ0 <π/2因为f 连续，由上确界定义必连续，由上确界定义必 有有f(θ0)=0，且对任意小，且对任意小 的的ε>0，总有，总有δ>0且且δ<ε，使，使f(θ0+δ)>0因为f(θ0+δ)g (θo+δ)=0，故必有，故必有g (θ0+δ)=0，由，由δ可任意小且可任意小且g连续，可知必连续，可知必 有有 g (θ0)=0，证毕证法，证毕证法二除用二除用 到到f、g的连续性外，还用到了上确界的性质的连续性外，还用到了上确界的性质 在解决实际问题时，注意观察和善于想象是十分重要的，在解决实际问题时，注意观察和善于想象是十分重要的，观察与想象不仅能发现问题隐含的某些属性，有时还能顺观察与想象不仅能发现问题隐含的某些属性，有时还能顺理成章地找到解决实际问题的钥匙。
本节的几个例子说明，理成章地找到解决实际问题的钥匙本节的几个例子说明，猜测也是一种想象力没有合理而又大胆的猜测，很难做猜测也是一种想象力没有合理而又大胆的猜测，很难做出具有创新性的结果开普勒的三大定律（尤其是后两条）出具有创新性的结果开普勒的三大定律（尤其是后两条）并非一眼就能看出的，它们隐含在行星运动的轨迹之中，并非一眼就能看出的，它们隐含在行星运动的轨迹之中，隐含在第谷记录下来的一大堆数据之中历史上这样的例隐含在第谷记录下来的一大堆数据之中历史上这样的例子实在太多了在获得了一定数量的资料数据后，人们常子实在太多了在获得了一定数量的资料数据后，人们常常会先去猜测某些结果，然后试图去证明它猜测一经证常会先去猜测某些结果，然后试图去证明它猜测一经证明就成了定理，而定理一旦插上想象的翅膀，又常常会被明就成了定理，而定理一旦插上想象的翅膀，又常常会被推广出许多更为广泛的结果即使猜测被证明是错误的，推广出许多更为广泛的结果即使猜测被证明是错误的，结果也决不是一无所获的失败而常常是对问题的更为深入结果也决不是一无所获的失败而常常是对问题的更为深入的了解 §2.9最短路径与最速方案问题最短路径与最速方案问题 例例5（最短路径问题）（最短路径问题） 设有一个半径为设有一个半径为 r 的圆形湖，圆心为的圆形湖，圆心为 O。
A、、B 位于湖的两侧，位于湖的两侧，AB连线过连线过O，见图现拟从现拟从A点步行到点步行到B点，在不得进入湖中的限点，在不得进入湖中的限 制下，问怎样的路径最近制下，问怎样的路径最近 ABOr将湖想象成凸出地面的木桩，将湖想象成凸出地面的木桩， 在在AB间拉一根软线，当间拉一根软线，当线被拉紧时将得到最短路径根据这样的想象，猜测线被拉紧时将得到最短路径根据这样的想象，猜测 可以如下得到最短路径：可以如下得到最短路径： 过过A作圆的切线切圆于作圆的切线切圆于E，过，过B作圆的切线切圆作圆的切线切圆 于于F最短路径为由线最短路径为由线 段段AE、弧、弧EF和线段和线段FB连接而成的连续曲线（根据对称性，连接而成的连续曲线（根据对称性，AE′，，弧弧E′F′，，F′B连接而成的连续曲线也是）连接而成的连续曲线也是）EFE′F′ 以上只是一种猜测，现在来证明这一猜测是正确的为此，以上只是一种猜测，现在来证明这一猜测是正确的为此，先介绍一下凸集与凸集的性质先介绍一下凸集与凸集的性质定义定义2.1（（凸集凸集）称集合）称集合 R为凸集，若为凸集，若x1、x2∈R及及λ∈[0，1]，总有总有λx1+（1-λ）x2∈R。
即若即若x1、x2∈R，则，则x1、x2的连线必整个地落的连线必整个地落 在在R中定理定理2.2（（分离定理分离定理）对平面中的凸）对平面中的凸 集集R与与R外的一点外的一点K，，存在直线存在直线 l , l 分离分离R与与K，即，即R与与K分别位于分别位于 l 的两侧（注：的两侧（注：对一般的凸对一般的凸 集集R与与R外的一点外的一点K，则存在超平面分，则存在超平面分 离离R与与K），见图klR下面证明猜想下面证明猜想 猜测证明如下：猜测证明如下：（方法一）（方法一）显然，显然， 由由AE、、EF、、FB及及AE′，，E′F′，，F′B围成围成的区域的区域 R是一凸集利用是一凸集利用分离定理分离定理易证最短径不可能经过易证最短径不可能经过R外的点，若不然，设外的点，若不然，设 Γ为最短路径，为最短路径，Γ过过R外的一点外的一点M，则，则必存在直必存在直 线线l分离分离M与与R，由于路径，由于路径Γ是连续曲线，由是连续曲线，由A沿沿Γ到到M，必交，必交l于于M1，由，由M沿沿Γ到到B又必交又必交l于于M2这样，直线这样，直线 段段M1M2的长度必小于路的长度必小于路 径径M1MM2的长度，与的长度，与Γ是是A到到B的的最短路径矛盾，至此，我们已证明最短路径必在凸集最短路径矛盾，至此，我们已证明最短路径必在凸集R内。
内不妨设路径经湖的上方到达不妨设路径经湖的上方到达B点，则弧点，则弧EF必在路径必在路径F上，又上，又直线段直线段AE是由是由A至至E的最短路径，直线的最短路径，直线FB是由是由F到到B的最短的最短路径，猜测得证路径，猜测得证ABOrEFE′F′M1M2MΓl 还可用还可用微积分微积分方法求弧长，根据计算证方法求弧长，根据计算证明满足限止条件的其他连续曲线必具有明满足限止条件的其他连续曲线必具有更大的长度；此外，本猜测也可用更大的长度；此外，本猜测也可用平面平面几何几何知识加以证明等知识加以证明等 根据猜测不难看出，根据猜测不难看出， 例例5中的条件可以大大放中的条件可以大大放松，可以不必松，可以不必 设设AB过圆心，甚至可不必设湖过圆心，甚至可不必设湖是圆形的例如对是圆形的例如对 下图，我们可断定由下图，我们可断定由A至至B的最短路径必的最短路径必 为为l1与与l2之一，其证明也不难类之一，其证明也不难类似给出 ABl1l2D到此为止，我们的研讨还只局限于平面之中，到此为止，我们的研讨还只局限于平面之中，其实上述猜测可十分自然地推广到一般空间其实上述猜测可十分自然地推广到一般空间中去。
中去1973年，年，J.W.Craggs证明了以上结果：证明了以上结果：若可行区域的边界是光滑曲面则最短路径必由下列弧组若可行区域的边界是光滑曲面则最短路径必由下列弧组成，它们或者是空间中的自然最短曲线，或者是可行区域成，它们或者是空间中的自然最短曲线，或者是可行区域的边界弧而且，组成最短路径的各段弧在连接点处必定的边界弧而且，组成最短路径的各段弧在连接点处必定相切 例例6 6 一辆汽车停于一辆汽车停于 A A处并垂直于处并垂直于ABAB方向，此方向，此汽车可转的最小圆半径为汽车可转的最小圆半径为 R，求不倒车而由，求不倒车而由 A A到到B B的最短路径的最短路径解解（情况（情况1））若若|AB|>2R，最短路径由，最短路径由 弧弧AC与切线与切线BC组组成（见成（见图图①① ）情况（情况2））若若|AB|<2R，则最短路径必居，则最短路径必居 于于图图②②（（a）、）、（（b））两曲线之中可以证明，两曲线之中可以证明， （（b））中的曲中的曲 线线ABC更短AR2RBRC①①②②ABoC（（a））CABo1o2（（b）） 例例7 驾驶一辆停于驾驶一辆停于A处与处与AB成成θ1角度的汽角度的汽车到车到B处去，已知处去，已知B处要求的停车方向必须处要求的停车方向必须 与与 AB成成θ2角，试找出最短路径（除可转角，试找出最短路径（除可转 的最小圆半径为的最小圆半径为R外，不受其他限止）。
外，不受其他限止）解解 根据根据Craggs定理并稍加分析可知，最短路径应定理并稍加分析可知，最短路径应在 l1与与 l2中，见图，比较中，见图，比较 l1与与 l2的长度，即可得到的长度，即可得到最短路径最短路径Al1l2Bθ2θ1 最速方案问题最速方案问题例例8 将一辆急待修理的汽车由静止开始沿一将一辆急待修理的汽车由静止开始沿一 直线方向推至相隔直线方向推至相隔 S米的修车处，设阻力不米的修车处，设阻力不 计计 ，推车人能使车得到的推力，推车人能使车得到的推力 f 满足满足:-B≤f≤A , f>0为推力，为推力，f<0为拉力问怎样为拉力问怎样推车可使车最快停于修车处推车可使车最快停于修车处 设该车的运动速度设该车的运动速度 为为υ=υ(t)，根据根据题意，题意，υ(0)= υ(T)=0，其中，其中T为推车所为推车所花的全部时间由花的全部时间由 于于-B≤f≤A，且，且f=mυ′，，可知可知-b≤υ′≤a（其中（其中m为汽车质为汽车质量，量，a=A/m,b=B/m） 据此不难将本例归纳为据此不难将本例归纳为 如如下的数学模型：下的数学模型： min T υ(0)=υ(T)=0 此问题为一泛函极值问题，求解十分困难，为得出一个最速方此问题为一泛函极值问题，求解十分困难，为得出一个最速方案。
我们作如下猜测：案我们作如下猜测：猜测猜测 最速方案为以最大推力将车推到某处，然后以最大拉最速方案为以最大推力将车推到某处，然后以最大拉力拉之，使之恰好停于修车处，其中转换点应计算求出力拉之，使之恰好停于修车处，其中转换点应计算求出证明证明 设设υ=υ(t)为在最速推车方案下汽车的速度，则有为在最速推车方案下汽车的速度，则有 设此方案不同于我们的猜测设此方案不同于我们的猜测现 从从O点出发，作射线点出发，作射线 y=at；从；从（（t,0））出发，作直出发，作直 线线y=-b(t-T)交交y=at于于A，由于，由于 ，曲线，曲线υ=υ(t)必位于三角形区必位于三角形区 域域DAT的内部，从而有的内部，从而有ΔOAT的面积大于的面积大于S在O到到T之间任取一之间任取一点点 T′ ，过，过T′作作AT的平行线交的平行线交OA于于A′ 显然ΔOA′T′的面积的面积S（（T′））是是T′的连续函数，当的连续函数，当T′=0时时S（0）=0，当，当T′=T时，时，S（T）>S，故由连续函数的性质存在，故由连续函数的性质存在 某某T′ oυtAT′TA′Sy=aty=-b(t-T) 圆周率是人类获得的最古老的数学概念圆周率是人类获得的最古老的数学概念之一，早在大约之一，早在大约37003700年前（即公元前年前（即公元前17001700年左右）的古埃及人就已经在年左右）的古埃及人就已经在 用用256/81256/81（约（约3.16053.1605）作为）作为ππ的近似值了几千年的近似值了几千年来，人们一直没有停止过求来，人们一直没有停止过求ππ的努力§2.10 π的计算的计算 古古 典典 方方 法法 分分 析析 方方 法法 其其 它它 方方 法法 概率方法概率方法 数值积分方法数值积分方法 古典方法古典方法用什么方法来计用什么方法来计 算算ππ的近似值呢？显然，不可能仅根据的近似值呢？显然，不可能仅根据圆周率的定义，用圆的周长去除以直径起先，人们采圆周率的定义，用圆的周长去除以直径起先，人们采用的都是用圆内接正多边形和圆外切正多边形来逼近的用的都是用圆内接正多边形和圆外切正多边形来逼近的古典方法古典方法6边形边形12边形边形24边形边形圆圆 阿基米德曾用圆内接阿基米德曾用圆内接 9696边形和圆外切边形和圆外切9696边形夹逼的方法证明了边形夹逼的方法证明了由由和和 导出导出 公元公元5世纪，祖冲之指出世纪，祖冲之指出比西方得到同样结比西方得到同样结果几乎早了果几乎早了1000年年 十五世纪中叶，阿尔十五世纪中叶，阿尔··卡西给出卡西给出ππ的的16位小数，打破了祖冲之的纪录位小数，打破了祖冲之的纪录 1579年年, ,韦达证明韦达证明 1630年年, ,最后一位用古典方法求最后一位用古典方法求ππ的人的人格林伯格也只求到了格林伯格也只求到了π的第的第39位小数位小数 分析方法分析方法从十七世纪中叶起，人们开始用更先进的从十七世纪中叶起，人们开始用更先进的分析方法来求分析方法来求ππ的近似值，其中应用的主的近似值，其中应用的主要工具是收敛的无穷乘积和无穷级数，在要工具是收敛的无穷乘积和无穷级数，在本节中我们将介绍一些用此类方法求本节中我们将介绍一些用此类方法求ππ近近似值的实例。 似值的实例 取取取取 1656年，沃里斯年，沃里斯( (Wallis) )证明证明 在微积分中我们学过泰勒级数，其中有在微积分中我们学过泰勒级数，其中有当当 取取取取 在中学数学中证明过下面的等式在中学数学中证明过下面的等式左边三个正方形左边三个正方形组成的矩形中，组成的矩形中， 由由 和和 可得可得 和和 的展开式的收敛速度的展开式的收敛速度都比都比 快得多快得多ACBD 麦琴麦琴( (Machin) )给出给出(Machin公式公式) )记记 ，， ，得，得此式求得了此式求得了π的第的第100位小数且全部正确位小数且全部正确 其它方法其它方法除用古典方法与分析方法求除用古典方法与分析方法求ππ的近似值以的近似值以外，还有人用其他方法来求外，还有人用其他方法来求ππ的近似值的近似值这里我们将介绍两种方法：这里我们将介绍两种方法： 概率方法概率方法 数值积分方法数值积分方法 概率方法概率方法取一个二维数组（取一个二维数组（x,yx,y），取一个充分大的），取一个充分大的正整正整 数数n n，重复，重复n n次，每次独立地从次，每次独立地从 （（0 0，，1 1）中随机地取一对）中随机地取一对 数数x x和和y y ，分别检验，分别检验x x2 2+y+y2 2≤1≤1是否成立。 是否成立 设设n n次试验中等式成立次试验中等式成立的共有的共有m m次，令次，令π≈4m/nπ≈4m/n但这种方法很难得到但这种方法很难得到π的较好的近似值的较好的近似值 Ø 数值积分方法数值积分方法还可用其它数值积还可用其它数值积分公式来求，但用分公式来求，但用此类方法此类方法效果也很效果也很难做得比用幂级数难做得比用幂级数展开更好展开更好 微分方程建模 §3.13.1 为什么要用三级火箭来发射人造卫星为什么要用三级火箭来发射人造卫星构造数学模型，以说明为什么不能用一级火箭而必须用多构造数学模型，以说明为什么不能用一级火箭而必须用多级火箭来发射人造卫星？为什么一般都采用三级火箭系统级火箭来发射人造卫星？为什么一般都采用三级火箭系统？？ 1 1、为什么不能用一级火箭发射人造卫星、为什么不能用一级火箭发射人造卫星? ? （（1 1）卫星能在轨道上运动的最低速度）卫星能在轨道上运动的最低速度 假设：假设：（（i i）） 卫星轨道为过地球中心的某一平面上的圆，卫星卫星轨道为过地球中心的某一平面上的圆，卫星 在此轨道上作匀速圆周运动在此轨道上作匀速圆周运动。 （（iiii））地球是固定于空间中的均匀球体，其它星球对卫地球是固定于空间中的均匀球体，其它星球对卫 星的引力忽略不计星的引力忽略不计 分析：分析：根据牛顿第三定律，地球对卫星的引力为根据牛顿第三定律，地球对卫星的引力为: 在地面有在地面有: :得得: : k=gR2 R R为地球半径，为地球半径，约为约为64006400公里公里 故引力故引力: : 假设(ii) dmm-dmvu-v假设(i)卫星所受到的引力也就是它作匀速圆周运动的向心力卫星所受到的引力也就是它作匀速圆周运动的向心力故又有故又有: :从而从而: :设设g=9.81g=9.81米米/ /秒秒2 2，得，得: : 卫星离地面高度卫星离地面高度 ( (公里公里) )卫星速度卫星速度 ( (公里公里/ /秒秒) )100100200200400400600600800800100010007.807.807.697.697.587.587.477.477.377.377.867.86（（2 2）火箭推进力及速度的分析）火箭推进力及速度的分析 假设：假设：火箭重力及空气阻力均不计火箭重力及空气阻力均不计 分析：分析：记火箭在时刻记火箭在时刻t的质量和速度分别为的质量和速度分别为m(t)和和υ(t) 有：有：记火箭喷出的气体相对于火箭的速度为记火箭喷出的气体相对于火箭的速度为u（（常数），常数）， 由动量守恒定理：由动量守恒定理： υυ0 0和和m m0 0一定的情况下，一定的情况下，火箭速度火箭速度υυ(t)(t)由喷发由喷发速度速度u u及质量比决定。 及质量比决定 故：故：由此解得：由此解得：( (3.11) ) （（2 2）火箭推进力及速度的分析）火箭推进力及速度的分析 现将火箭现将火箭——卫星系统的质量分成三部分：卫星系统的质量分成三部分： （（i））mP（（有效负载，如卫星）有效负载，如卫星）（（ii））mF（（燃料质量）燃料质量）（（iii））mS（（结构质量结构质量——如外壳、燃料容器及推进器）如外壳、燃料容器及推进器） 最终质量为最终质量为mP + mS ，，初始速度为初始速度为0，，所以末速度：所以末速度：根据目前的技术条件和燃料性根据目前的技术条件和燃料性能，能，u只能达到只能达到3公里公里/秒，即使秒，即使发射空壳火箭，其末速度也不发射空壳火箭，其末速度也不超过超过6.6公里公里/秒 目前根本不目前根本不可能用一级火箭发射人造卫星可能用一级火箭发射人造卫星火箭推进力在加速整个火箭时，其火箭推进力在加速整个火箭时，其实际效益越来越低如果将结构质实际效益越来越低如果将结构质量在燃料燃烧过程中量在燃料燃烧过程中不断减少，那不断减少，那么末速度能达到要求吗？么末速度能达到要求吗？ 2 2、理想火箭模型、理想火箭模型 假设：假设： 记结构质量记结构质量mS在在mS + mF中占的比例为中占的比例为λ，，假设火假设火箭能随时抛弃无用的结构，结构质量与燃料质量以箭能随时抛弃无用的结构，结构质量与燃料质量以λ与与（（1-λ））的比例同时减少。 的比例同时减少 建模建模: : 由由 得到：得到：解得：解得： 理想火箭与一级火箭最大的区别在于，当火箭燃料理想火箭与一级火箭最大的区别在于，当火箭燃料耗尽时，结构质量也逐渐抛尽，它的最终质量为耗尽时，结构质量也逐渐抛尽，它的最终质量为mP，， 所以最终速度为：所以最终速度为： 只要只要m0足够大，我们可以足够大，我们可以使卫星达到我们希望它具使卫星达到我们希望它具有的任意速度有的任意速度考虑到空气阻力和重力等因素，估考虑到空气阻力和重力等因素，估计（按比例的粗略估计）发射卫星计（按比例的粗略估计）发射卫星要使要使υ=10.5公里公里/秒才行，则可推秒才行，则可推算出算出m0/ mp约为约为51,即发射一吨重即发射一吨重的卫星大约需要的卫星大约需要50吨重的理想火箭吨重的理想火箭 3 3、理想过程的实际逼近、理想过程的实际逼近——多级火箭卫星系统多级火箭卫星系统 记火箭级数为记火箭级数为n，，当第当第i级火箭的燃料烧尽时，第级火箭的燃料烧尽时，第i+1级火级火箭立即自动点火，并抛弃已经无用的第箭立即自动点火，并抛弃已经无用的第i级火箭用级火箭用mi表示第表示第i级火箭的质量，级火箭的质量，mP表示有效负载。 表示有效负载 先作如下假设：先作如下假设： （（i））设各级火箭具有相同的设各级火箭具有相同的λ ,即即i级火箭中级火箭中λmi为结构为结构质量，（质量，（1-λ））mi为燃料质量为燃料质量 （（ii））设燃烧级初始质量与其负载质量之比保持不变，设燃烧级初始质量与其负载质量之比保持不变，并记比值为并记比值为k k 考虑二级火箭：考虑二级火箭： 由由3.11式，当第一级火箭燃烧完时，其末速度为：式，当第一级火箭燃烧完时，其末速度为： 当第二级火箭燃尽时，末速度为：当第二级火箭燃尽时，末速度为： 该该假设有点强加假设有点强加的味道，先权作的味道，先权作讨论的方便吧讨论的方便吧 又由假设（又由假设（ii），），m2=kmP，，m1=k(m2+mP)，，代入上式，代入上式，仍设仍设u=3公里公里/秒，且为了计算方便，近似取秒，且为了计算方便，近似取λ=0.1，，则可则可得：得： 要使要使υ2=10.5公里公里/秒，则应使秒，则应使: 即即k≈11.2，，而而: 类似地，可以推算出三级火箭：类似地，可以推算出三级火箭： 在同样假设下在同样假设下: : 要使要使υ3=10.5公里公里/秒，则秒，则(k+1)/(0.1k+1)≈3.21，，k≈3.25，，而而（m1+ m2+ m3+ mP）/ mP≈77。 三级火箭比二级火箭三级火箭比二级火箭几乎节省了一半几乎节省了一半 是否三级火箭就是最省呢是否三级火箭就是最省呢？最简单的方法就是对四？最简单的方法就是对四级、五级等火箭进行讨论级、五级等火箭进行讨论 考虑考虑N N级火箭：级火箭： 记记n级火箭的总质量（包含有效负载级火箭的总质量（包含有效负载mP）为）为m0 ，，在在相同的假设下可以计算出相应的相同的假设下可以计算出相应的m0/ mP的的值，见表值，见表3-2n（（级数）级数）1 2 3 4 5 …… ∞∞（理想）（理想） 火箭质量（吨）火箭质量（吨）/ 149 77 65 60 …… 50表3-2由于工艺的复杂性及每节火箭由于工艺的复杂性及每节火箭都需配备一个推进器，所以使都需配备一个推进器，所以使用四级或四级以上火箭是不合用四级或四级以上火箭是不合算的，三级火箭提供了一个最算的，三级火箭提供了一个最好的方案好的方案当然若燃料的价钱很便宜当然若燃料的价钱很便宜而推进器的价钱很贵而推进器的价钱很贵,且且制作工艺非常复杂的话，制作工艺非常复杂的话，也可选择二级火箭。 也可选择二级火箭 4 4、火箭结构的优化设计、火箭结构的优化设计 3 3中中已经能说过假设已经能说过假设(ii)(ii)有点强加的味道；现去掉该有点强加的味道；现去掉该假设假设，，在各级火箭具有相同在各级火箭具有相同λλ的粗糙假设下，来讨论火箭的粗糙假设下，来讨论火箭结构的最优设计结构的最优设计 W1=m1+……+ mn+ mP W2=m2+……+ mn+ mP……Wn= mn+ mPWn+1= mP记记应用（应用（3.113.11）可求得末速度：）可求得末速度： 记记则则又又问题化为，在问题化为，在υn一定的条件下，求使一定的条件下，求使k1 k2……kn最小最小 解条件极值问题：解条件极值问题： 或等价地求解无约束极值问题：或等价地求解无约束极值问题： 可以解出最优结构设计应满足：可以解出最优结构设计应满足： 火箭结构优化设计讨论中火箭结构优化设计讨论中我们得到与假设（我们得到与假设（ii））相相符的结果，这说明前面的符的结果，这说明前面的讨论都是有效的！讨论都是有效的！ §3.23.2 药物在体内的分布药物在体内的分布 何为房室系统？何为房室系统？ 在用微分方程研究实际问题时，人们常常采用一种叫在用微分方程研究实际问题时，人们常常采用一种叫“房室房室系统系统”的观点来考察问题。 根据研究对象的特征或研究的不同精的观点来考察问题根据研究对象的特征或研究的不同精度要求，我们把研究对象看成一个整体（单房室系统）或将其剖度要求，我们把研究对象看成一个整体（单房室系统）或将其剖分成若干个相互存在着某种联系的部分（多房室系统）分成若干个相互存在着某种联系的部分（多房室系统） 房室具有以下特征：它由考察对象均匀分布而房室具有以下特征：它由考察对象均匀分布而成，房室中考察对象的数量或浓度（密度）的变成，房室中考察对象的数量或浓度（密度）的变化率与外部环境有关，这种关系被称为化率与外部环境有关，这种关系被称为“交换交换”且交换满足着总量守衡在本节中，我们将用房且交换满足着总量守衡在本节中，我们将用房室系统的方法来研究药物在体内的分布室系统的方法来研究药物在体内的分布交换环境内部单房室系统均匀分布 药物的分解与排泄（输出）速率通常被认为是与药物当前的药物的分解与排泄（输出）速率通常被认为是与药物当前的浓度成正比的，即：浓度成正比的，即： 药物分布的单房室模型药物分布的单房室模型 单房室模型是最简单的模型，它假设：体内药物在任一时刻单房室模型是最简单的模型，它假设：体内药物在任一时刻都是均匀分布的，设都是均匀分布的，设t时刻体内药物的总量为时刻体内药物的总量为x(t)；；系统处于一种系统处于一种动态平衡中，即成立着关系式：动态平衡中，即成立着关系式： 药物的输入规律与给药的方式有药物的输入规律与给药的方式有关。 下面，我们来研究一下在几种常关下面，我们来研究一下在几种常见的给药方式下体内药体的变化规律见的给药方式下体内药体的变化规律 机体环境药物总量图3-8 假设药物均匀分布 情况情况1 1 快速静脉注射快速静脉注射机体环境只输出不输入房室其解其解为：为：药物的浓度：药物的浓度： 与放射性物质类似，医学上将血浆药物浓度衰减一半所需与放射性物质类似，医学上将血浆药物浓度衰减一半所需的时间称为药物的血浆半衰期：的时间称为药物的血浆半衰期： 负增长率的Malthus模型 在快速静脉注射时，总量为在快速静脉注射时，总量为D的药物在瞬间被注入体内的药物在瞬间被注入体内设机体的体积为设机体的体积为V，，则我们可以近似地将系统看成初始总量为则我们可以近似地将系统看成初始总量为D，，浓度为浓度为D/V，，只输出不输入的房室，即系统可看成近似地只输出不输入的房室，即系统可看成近似地满足微分方程：满足微分方程： （3.12） 情况情况2 2 恒速静脉点滴恒速静脉点滴 机体环境恒定速率输入房室药物似恒速点滴方式进入体内，即药物似恒速点滴方式进入体内，即: : 则体内药物总量满足：则体内药物总量满足： （x(0)=0） （3.13） 这是一个一阶常系数线性方程，其解为：这是一个一阶常系数线性方程，其解为： 或或易见易见：：称为稳态血药浓度 对于多次点滴，设点滴时间为对于多次点滴，设点滴时间为T1，，两次点滴之间的间隔时两次点滴之间的间隔时间设为间设为T2，，则在第一次点滴结束时病人体内的药物浓度可由上则在第一次点滴结束时病人体内的药物浓度可由上式得出。 其后式得出其后T2时间内为情况时间内为情况1 1故：（第一次） 0≤t≤T1 T1≤t≤T1 +T2 类似可讨论以后各次点滴时类似可讨论以后各次点滴时的情况，区别只在初值上的的情况，区别只在初值上的不同第二次点滴起，患者不同第二次点滴起，患者 体内的初始药物浓度不为零体内的初始药物浓度不为零 情况情况3 3 口服药或肌注口服药或肌注 y(t)x(t)K1yK1x环境机体外部药物 口服药或肌肉注射时，药物的吸收方式与点滴时不同，药口服药或肌肉注射时，药物的吸收方式与点滴时不同，药物虽然瞬间进入了体内，但它一般都集中与身体的某一部位，物虽然瞬间进入了体内，但它一般都集中与身体的某一部位，靠其表面与肌体接触而逐步被吸收设药物被吸收的速率与存靠其表面与肌体接触而逐步被吸收设药物被吸收的速率与存量药物的数量成正比，记比例系数为量药物的数量成正比，记比例系数为K1，，即若记即若记t时刻残留药物时刻残留药物量为量为y(t)，则，则y满足：满足： D为口服或肌注药物总量 因而因而：：所以所以：：解得解得：：从而药物浓度从而药物浓度：： 图图3-93-9给出了上述三种情况下体内血药浓度的变化曲线。 给出了上述三种情况下体内血药浓度的变化曲线容易看出，快速静脉注射能使血药浓度立即达到峰值，常用于容易看出，快速静脉注射能使血药浓度立即达到峰值，常用于急救等紧急情况；口服、肌注与点滴也有一定的差异，主要表急救等紧急情况；口服、肌注与点滴也有一定的差异，主要表现在血药浓度的峰值出现在不同的时刻，血药的有效浓度保持现在血药浓度的峰值出现在不同的时刻，血药的有效浓度保持时间也不尽相同时间也不尽相同图3-9 我们已求得三种常见给药方式下的血药浓度我们已求得三种常见给药方式下的血药浓度C C( (t t) )，，当然也当然也容易求得血药浓度的峰值及出现峰值的时间，因而，也不难根容易求得血药浓度的峰值及出现峰值的时间，因而，也不难根据不同疾病的治疗要求找出最佳治疗方案据不同疾病的治疗要求找出最佳治疗方案 上述研究是将机体看成一个均匀分布的同质单元，故上述研究是将机体看成一个均匀分布的同质单元，故被称单房室模型，但机体事实上并不是这样药物进入血被称单房室模型，但机体事实上并不是这样药物进入血液，通过血液循环药物被带到身体的各个部位，又通过交液，通过血液循环药物被带到身体的各个部位，又通过交换进入各个器官。 因此，要建立更接近实际情况的数学模换进入各个器官因此，要建立更接近实际情况的数学模型就必须正视机体部位之间的差异及相互之间的关联关系，型就必须正视机体部位之间的差异及相互之间的关联关系，这就需要多房室系统模型这就需要多房室系统模型IIIk12k21两房室系统图3-10 图图3-103-10表示的是一种常见的两房室模型，表示的是一种常见的两房室模型，其间的其间的k k1212表示由室表示由室I I渗透到室渗透到室IIII的变化率前的的变化率前的系数，而系数，而k k2121则表示由室则表示由室IIII返回室返回室I I的变化率前的变化率前的系数，它们刻划了两室间的内在联系，其的系数，它们刻划了两室间的内在联系，其值应当用实验测定，使之尽可能地接近实际值应当用实验测定，使之尽可能地接近实际情况 当差异较大的部分较多当差异较大的部分较多时，可以类似建立多房时，可以类似建立多房室系统室系统，即，即N N房室系统房室系统 离散优化模型 比较算法的好坏，从不同的角度出发，有各种不同的标准在这里，比较算法的好坏，从不同的角度出发，有各种不同的标准在这里，我们仅就算法的计算速度作一个十分粗略的比较。 我们仅就算法的计算速度作一个十分粗略的比较 例例1 1 （整理问题）给定（整理问题）给定n n个实数个实数a a1 1, , a a2 2,…,,…, a an n，要求将它整理成由小到大，要求将它整理成由小到大排列（或由大到小排列）的顺序：排列（或由大到小排列）的顺序：b b1 1, , b b2 2,…,,…, b bn n，，b b1 1≤≤ b b2 2≤≤……≤≤ b bn n算法（算法1 1）） 取出取出a a1 1, , a a2 2,…,,…, a an n中的最小者，令其为中的最小者，令其为b b1 1从a a1 1, , a a2 2,…,,…, a an n中中去除去除b b1 1，在余下的，在余下的n n——1 1个数中选出最小者，令其为个数中选出最小者，令其为b b2 2，，……，直至得到，直至得到b b1 1, , b b2 2,…,,…, b bn n容易看出，为了排出容易看出，为了排出b b1 1, , b b2 2,…,,…, b bn n，算法工作了，算法工作了 次比较算法（算法2 2））步步0 0 b b1 1←←a a1 1 步步1 1 设已有设已有b b1 1, ,……, ,b bk k (1 (1≤≤k 若）若k k+1<+1 可见在较大规模的秒可见在较大规模的整理问题时，算法整理问题时，算法2 2明显优于算法明显优于算法1 1 现设一台每小时能作现设一台每小时能作M M次运算的计算机，并设求解某一问题已有了两个不次运算的计算机，并设求解某一问题已有了两个不同的算法算法同的算法算法A A对规模为对规模为n n的实例约需作的实例约需作n n2 2次运算而算法次运算而算法B B则约需作则约需作2 2n n次运次运算于是，运用算法算于是，运用算法A A在一小时内可解一个规模为在一小时内可解一个规模为 的实例，而算法的实例，而算法B B则则只能解一个规模为只能解一个规模为loglog2 2M M的实例两者的差别究竟有多大呢？让我们来对的实例两者的差别究竟有多大呢？让我们来对比一下假如计算机每秒可作比一下假如计算机每秒可作100100万次运算，则算法万次运算，则算法A A一小时大约可解一个一小时大约可解一个规模为规模为n n=60000=60000的实例，而算法的实例，而算法B B在一小时内只能解一个规模在一小时内只能解一个规模 的实的实例且例且n n每增加每增加1 1，算法，算法B B求解时所化的时间就将增加求解时所化的时间就将增加1 1倍。 例如算法倍例如算法B B求解一求解一个个n n=50=50的实例将连续运算的实例将连续运算357357年多现设计算速度提高了现设计算速度提高了 100100倍，这对计算机来讲已是一个相当大的改进此倍，这对计算机来讲已是一个相当大的改进此时算法时算法A A可解问题的规模增大了可解问题的规模增大了1010倍，而算法倍，而算法B B可解问题的规模只增加了可解问题的规模只增加了loglog2 2100<7100<7前者可解问题的规模成倍成倍地增加而后者则几乎没有什么前者可解问题的规模成倍成倍地增加而后者则几乎没有什么改变，今天无法求解的问题将来也很少有希望解决改变，今天无法求解的问题将来也很少有希望解决 由于这一原因，人们对由于这一原因，人们对算法作了如下的分类：算法作了如下的分类： （多项式算法）设（多项式算法）设A A是求解某一问题是求解某一问题D D的一个算法，对的一个算法，对D D的一个规模为的一个规模为n n的实的实例，用例，用f fA A( (D D, ,n n) )表示用算法表示用算法A A求解这一实例所作的初等运算的次数若存在求解这一实例所作的初等运算的次数若存在一个多项式一个多项式P P( (n n) )和一个正整数和一个正整数N N，当，当n n≥≥N N时总有时总有 f fA A( (D D, ,n n) )≤≤P P( (n n) )（不论求解（不论求解D D的怎样的实例，只要其规模为的怎样的实例，只要其规模为n n），则称算法），则称算法A A为求解问题为求解问题D D的一个多项式的一个多项式时间算法，简称多项式算法。 时间算法，简称多项式算法定义定义4.14.1定义定义4.24.2（指数算法）设（指数算法）设B B是求解某一问题是求解某一问题D D的一个算法，的一个算法，f f B B( (D D, ,n n) )为用算法为用算法B B求解求解D D的一个规模为的一个规模为n n的实例时所用的初等运算次数，若存在一个常数的实例时所用的初等运算次数，若存在一个常数k k>0>0，对任，对任意正整数意正整数N N，总可找到一个大于，总可找到一个大于N N的正整数的正整数n n及及D D的一个规模为的一个规模为n n的实例，对的实例，对这一实例有这一实例有f fB B( (D D, ,n n)=)=O O (2 (2knkn) )，则称，则称B B为求解问题为求解问题D D的一个指数算法的一个指数算法 多项式算法被称为是多项式算法被称为是“好好”的算法即所谓有效算法，而指数算法则一般被的算法即所谓有效算法，而指数算法则一般被认为是认为是“坏坏”算法，因为它只能求解规模极小的实例算法，因为它只能求解规模极小的实例 下面的表下面的表1 1列出了在规模大约为列出了在规模大约为n n时各类算法的计算量，而表时各类算法的计算量，而表2 2则反映了当计算机速度提高则反映了当计算机速度提高1010倍、倍、100100倍时，各类算法在倍时，各类算法在1 1小时小时内可求解的问题的模型的增长情况，（前三个是多项式时间的，内可求解的问题的模型的增长情况，（前三个是多项式时间的，后两个是指数时间的）后两个是指数时间的） 表表1 1 几个多项式和指数时间复杂性函数增长情况几个多项式和指数时间复杂性函数增长情况算法要求的算法要求的计计算量算量规规模模n的近似的近似值值101001000n101001000nlogn336649966n31031061092n10241.27×10301.05×10301n!3628800101584×102567 表表2 12 1小时内可解的问题实例的规模小时内可解的问题实例的规模算法要求的算法要求的计计算量算量用用现现在在计计算机算机用快用快10倍倍计计算机算机用快用快100倍倍计计算机算机nN110N1100N1nlognN28.2N267N2n3N32.15N34.64N32nN4N4+3.32N4+6.64n!N5≤N5+2≤N5+3由定义由定义1 1知知4 4n n2 2与与2 2n n2 2都是都是O O( (n n2 2) )，，n nloglogn n+3+3n n是是O O( (n nloglogn n) )，我们在以后分析时，我们在以后分析时间复杂性函数时也往往忽略常系数和增长速度较慢的项，因为前者可通间复杂性函数时也往往忽略常系数和增长速度较慢的项，因为前者可通过提高计算机速度来提高效率，而后者增长速度最快的项才是决定效率过提高计算机速度来提高效率，而后者增长速度最快的项才是决定效率的关键因素。 的关键因素 2004年浙江大学数学建模竞赛年浙江大学数学建模竞赛（（B题）通讯卫星上的开关设置题）通讯卫星上的开关设置 地面上存在着n个接收站与n个发送站，而在通讯卫星上则设置了若干种开关模式开关模式可用矩阵P=(pij)来表示，若卫星可接收发送站i发射的信息并将信息传送回地面的接收站j时，矩阵中的元素pij =1，否则pij =0通讯卫星上的接收发送任务也可以用一个矩阵T=（tij）来表示，其元素tij为需经通讯卫星传递的由i发点发送到j接收点的信息量的传送时间长度由于技术上的原因，当发送站i在发送给接收站j信息时，它不能同时发送给别的接收站信息；同样，当接收站j在接收发送站i的信息时，也不能同时接收其他发送站发送的信息你的任务是： n设计一组开关模式，k=1, …,r（注：r应当尽可能小），使得对任意给定的任务矩阵T，卫星开关设置均能完成要求的发接收任务n设计一个算法，在发接收任务T给出后，可根据你设计的开关模式（k=1,…,r）求出各开关的使用时间λk，使得在完成预定传送任务的前提下使用各开关模式的总时间最短n同样由于技术上的原因，开关模式的总数r有一个上限当需要传送的任务数数量较大时，仍无法分派任务。 请你想一些办法来解决这一困难，（当然，这时你可能要作出一些牺牲，即传送时间可能会增加一些） 问题及模型问题及模型问题的标准形式为：在地面上存在着问题的标准形式为：在地面上存在着n个收站与个收站与n个发战，而在通讯卫星上个发战，而在通讯卫星上则设置了若干种开关模式开关模式可用矩阵则设置了若干种开关模式开关模式可用矩阵P=(pij)来表示，若卫星可接来表示，若卫星可接收发射站收发射站i发射的信息并将信息传送回地面的接收站发射的信息并将信息传送回地面的接收站j时，矩阵元素时，矩阵元素pij =1，否，否则则pij =0通讯卫星的接发任务也可用一矩阵通讯卫星的接发任务也可用一矩阵T=（（tij）来表示，其元素）来表示，其元素tij为需为需经通讯卫星传递的由经通讯卫星传递的由i发点发送到发点发送到j接受点的信息量的传送时间长度问题要接受点的信息量的传送时间长度问题要求求求求r并设计一组开关模式并设计一组开关模式Pk，，k=1, …,r及模式及模式Pk的使用时间的使用时间λk，使得在完，使得在完成预定传送任务的前提下各开关模式使用的总时间最短，即要求求解下面成预定传送任务的前提下各开关模式使用的总时间最短，即要求求解下面的问题：的问题： 例例1 1 设设 这是一个有这是一个有3个发送站与个发送站与3个接收站的实例，个接收站的实例，tij在矩阵中已给出，例如在矩阵中已给出，例如由发站由发站1传送到收站传送到收站1的通讯量为的通讯量为3单位时间等。 单位时间等分析分析 容易看出，三个发站需传送的时间分别为容易看出，三个发站需传送的时间分别为6、、5、、5；而三个收站需接；而三个收站需接收的时间分别为收的时间分别为6、、3、、7为完成全部传送任务，通讯卫星总传送时间至少为完成全部传送任务，通讯卫星总传送时间至少应为应为7单位时间，即的下界为单位时间，即的下界为7由于技术上的原因，当发站由于技术上的原因，当发站i在发送给收站在发送给收站j信息时，它不能同时发送给别的信息时，它不能同时发送给别的收站信息；同样，当收站收站信息；同样，当收站j在接收发站在接收发站i的信息时，也不能同时接收其他发站的信息时，也不能同时接收其他发站发送的信息这一要求说明，任一开关模式发送的信息这一要求说明，任一开关模式Pk应具有以下性质：（应具有以下性质：（1））Pk的的每一行中有且只有一个每一行中有且只有一个1，每一列中也有且只有一个，每一列中也有且只有一个1；（；（2）所有的）所有的1均位均位于不同的行列中于不同的行列中满足（满足（1）、（）、（2）的矩阵）的矩阵 被称为置换矩阵，被称为置换矩阵，n阶置换矩阵阶置换矩阵Pk共有共有n!个，当个，当n较大时，我们不可能在通讯卫星上设置这么多种不同的开关模式。 因而，为较大时，我们不可能在通讯卫星上设置这么多种不同的开关模式因而，为了设计出切实可行的开关模式，我们还得另想办法了设计出切实可行的开关模式，我们还得另想办法 （设计方法（设计方法1））注意到注意到Pk每行（或列）元素之和均为每行（或列）元素之和均为1，故不管如何指派开关的使用时间，故不管如何指派开关的使用时间（即不论如何取（即不论如何取λk），矩阵），矩阵均具有某些特殊的性质，例如其行和（及列和）均为同一常数这样的矩均具有某些特殊的性质，例如其行和（及列和）均为同一常数这样的矩阵构成一个线性空间（参见逻辑模型第一节阵构成一个线性空间（参见逻辑模型第一节 Dürer魔方），为减少开关模魔方），为减少开关模式的种类，可取此空间的一组基底作为开关模式在使用这种开关模式时，式的种类，可取此空间的一组基底作为开关模式在使用这种开关模式时，无论无论T的元素的元素tij怎么取，通讯卫星对每一发（收）点的开通时间总和是恒定怎么取，通讯卫星对每一发（收）点的开通时间总和是恒定的在这种开关模式下，可按如下方式指派各开关模式的使用时间：的在这种开关模式下，可按如下方式指派各开关模式的使用时间： 步步1 先将先将T改变为改变为 ，， 满足：满足：（（1）） ≥T（（2）记）记 ，，步步2 用用Pk表示表示 ，即将，即将 分解为分解为（（r为空间的维数）为空间的维数） 将将T化为化为 的方法一般有无穷多种，如可如下化法：的方法一般有无穷多种，如可如下化法：令令 事实上，事实上， ，（即通讯卫星传送总时间的下界）。 即通讯卫星传送总时间的下界）令令 其中其中 用这种方法化例中的用这种方法化例中的T，得到，得到的任一行（或列）中元素之和均为的任一行（或列）中元素之和均为7 定义定义1 1称行和、列称行和、列%%和均相等的矩阵为双随机矩阵（和均相等的矩阵为双随机矩阵（Doubly stochastic matrix））定理定理1 1 （（Birkhoff定理，定理，1944）任一）任一n阶双随机矩阵均可写成至多阶双随机矩阵均可写成至多 （（n－－1））2+1个置换矩阵的非线性组合个置换矩阵的非线性组合的分解可如下进行：的分解可如下进行：步步1 选取由选取由Pij>0可推出可推出 >0的置换矩阵的置换矩阵P步步2 确定确定 步步3 取取 ，用，用 －－ 代替代替步步4 若若 =0，停；否则，返回步，停；否则，返回步1例例2. 2. 为方便起见，我们来分解一个元素均为非负整数的为方便起见，我们来分解一个元素均为非负整数的3阶双随机矩阵，阶双随机矩阵，（由（由Birkhoff定理，定理，r≤5）） 解：取解：取 ，，λ=min {1, 3, 3 } =1分解成分解成，再取，再取 因因min{ 5, 5, 3} = 3，又有，又有，取，取 于是又有于是又有 易得分解结果为：易得分解结果为： 尚需解决的问题是如何求尚需解决的问题是如何求P，使得，使得Pij>0必有必有 。 读者不难发现，此问读者不难发现，此问题可以通过求解一个两分图上的最大流（或最大匹配）来实现，计算量题可以通过求解一个两分图上的最大流（或最大匹配）来实现，计算量为为O(n4)，是多项式时间可解的具体方法为：作一两分图，若，是多项式时间可解的具体方法为：作一两分图，若 ，，则作边（则作边（i, j），令边容量为），令边容量为1，这样，可作出，这样，可作出P的充要条件是该最大流问的充要条件是该最大流问题的最大流量为题的最大流量为n对例9.33，，n=3由于所有由于所有 ，先取，先取 ，， －－P1为为 相应的两分图为：相应的两分图为：于是又可求得于是又可求得 ，相应的两分图为：，相应的两分图为： 又可得又可得 ，，…，如此下去，直到作不出，如此下去，直到作不出P为至，为至，由于由于 的特殊性质及的特殊性质及Birkhoff定理，上述分解必能在不超过定理，上述分解必能在不超过r= (n－－1)2 + 1步内终止步内终止上述开关设计方法要求在通讯卫星上设置上述开关设计方法要求在通讯卫星上设置(n－－1)2 + 1种不同的开关模式种不同的开关模式（即（即Pk），当），当n稍大时，稍大时，(n－－1)2 + 1仍显得太大而使得使用时不便。 例如，仍显得太大而使得使用时不便例如，当当n=41时，时，(n－－1)2 + 1=| 60 |为实用方便，人们研究了限止开关模式个为实用方便，人们研究了限止开关模式个数的相应问题数的相应问题 若要求若要求r≤n，即要求通讯卫星上至多设置，即要求通讯卫星上至多设置n种开关模式，则问题化为令种开关模式，则问题化为令r≤n，求不超过，求不超过n个置换矩阵个置换矩阵Pk及及λk，使之满足：，使之满足： min S.t 为了使任意一对发射法与接收站之间的传送均为可能实现的，自然应要求为了使任意一对发射法与接收站之间的传送均为可能实现的，自然应要求Pk满足满足（5.1） （5.2） （右面的矩阵有（右面的矩阵有n2个值为个值为1的分量，每一的分量，每一Pk 恰有恰有n个个1分量）故分量）故r=n容易看出，（容易看出，（5.1）隐含着）隐含着T的每一元素只能被唯一的的每一元素只能被唯一的P复盖，即复盖，即T的元素在的元素在分解中是不可分割的，这当然是一个好性质，使实际操作时较为方便，但分解中是不可分割的，这当然是一个好性质，使实际操作时较为方便，但可惜的是对一般的双随机矩阵，分解很可能无解可惜的是对一般的双随机矩阵，分解很可能无解。 例例3 3 若取若取， ， （注意：（注意：T已是双随机矩阵，行和列和均为已是双随机矩阵，行和列和均为10）） 则min S.t 的解为的解为λ1=3，，λ2=4，，λ3=5 （大于（大于10）而）而 但等号经常并不成立但等号经常并不成立1985年，年，F·Rendel证明，在给定满足（证明，在给定满足（5.2）的置）的置换矩阵换矩阵P1,…,Pn后，求解问题（后，求解问题（5.1）是）是NP难的，从而不可能存在多项式难的，从而不可能存在多项式时间算法，除非时间算法，除非P=NP现要求现要求r≤2n一种自然而方便的开关设置为引入两组各有一种自然而方便的开关设置为引入两组各有n个开关模式的置换矩阵个开关模式的置换矩阵P1,…,Pn，，Q1,…,Qn，满足下面的（，满足下面的（5.3)式：式：例如，当例如，当n=3时，可令：时，可令： （注：这种设置方法保持了其内在的对称性，不失为一种明智的做法注：这种设置方法保持了其内在的对称性，不失为一种明智的做法现在，我们来分解例现在，我们来分解例9.33中的双随机矩阵中的双随机矩阵 ，令，令 =，得方程组，得方程组 求出各对角线与反对角线上的三个元素之和，并作一些简单的消去运算；求出各对角线与反对角线上的三个元素之和，并作一些简单的消去运算；将矩阵的所有元素相加，可得下面的方程组：将矩阵的所有元素相加，可得下面的方程组：注意到（注意到（5.3），易证空间），易证空间 的维数为的维数为 5，，故故 之一可任取，（稍加注意即可保持非负性），之一可任取，（稍加注意即可保持非负性），例如，令例如，令μ3=0，求得，求得 ，故有，故有 读者不难验证，上述方法可推广到读者不难验证，上述方法可推广到n是奇数的一般情况。 事实，由各对角线是奇数的一般情况事实，由各对角线元素之和可导出元素之和可导出n－－1个方程，由各反对角线元素之和又可导出个方程，由各反对角线元素之和又可导出n－－1个方程，个方程，加上矩阵所有元素之和导出的等式，共计可导出加上矩阵所有元素之和导出的等式，共计可导出2 n－－1个方程，并易知它个方程，并易知它们是独立的另一方面空间们是独立的另一方面空间 的维数恰为的维数恰为2 n－－1，故，故 之一可任取，而通过方程组解得所有的之一可任取，而通过方程组解得所有的 ，，（只须注意保持其非负性即可）（只须注意保持其非负性即可）但当但当n为偶数时，情况就不大相同了让我们先来观察一下为偶数时，情况就不大相同了让我们先来观察一下n=4的情况当当n=4时，时， 易见，易见， 具有非常特殊的结构，一般的偶数阶双随机矩具有非常特殊的结构，一般的偶数阶双随机矩阵，即使其元素是非负整数，也无法用阵，即使其元素是非负整数，也无法用Pk、、Qk来分解当当 具有上述结构时，能否用具有上述结构时，能否用Pk和和Qk来分解呢？易见，由各对角线元素之来分解呢？易见，由各对角线元素之和可导出：和可导出： 另外，由反对角线元素之和又可导出另外，由反对角线元素之和又可导出上述方程中只有上述方程中只有6个是独立的，且已不可能再得出新的独立方程，（读者可个是独立的，且已不可能再得出新的独立方程，（读者可自行分析之）故可选取其中自行分析之）故可选取其中2个的值，进而可解出其余。 例如，若令个的值，进而可解出其余例如，若令 4= 3=0，可得，可得 2=1，， 1=0，进而可求得，进而可求得 1=2，， 4=3，， 3=3及及 2=4，， 已达到下界已达到下界易见，易见，P1 + P3 = Q2 + Q4，，P2 + P4 = Q1 + Q3，故空间，故空间 的维的维数为数为6，与上面的分析是一致的与上面的分析是一致的读者可将上述讨论推广到读者可将上述讨论推广到n为一般偶数的情况，分析方法是完全类似的为一般偶数的情况，分析方法是完全类似的 当当n是偶数时，我们虽无法将一般的双随机矩阵分解为是偶数时，我们虽无法将一般的双随机矩阵分解为Pk 、、Qk的非负组合，的非负组合，但上述讨论仍然是十分有意义的首先，要求完成的任务矩阵是但上述讨论仍然是十分有意义的首先，要求完成的任务矩阵是T，在将，在将T转换成时我们可尽量使具有上述的特殊结构（有兴趣的读者可自行研究这转换成时我们可尽量使具有上述的特殊结构（有兴趣的读者可自行研究这一问题），只要能做到这一点，即可给出一个达到下界的开关模式的指派一问题），只要能做到这一点，即可给出一个达到下界的开关模式的指派方式。 其次，即使这样的努力没有成功，也容易给出一个具有上述特殊结方式其次，即使这样的努力没有成功，也容易给出一个具有上述特殊结构构 矩阵，矩阵，并使并使 尽可能地小，即给出一种开关指派的近似最佳方尽可能地小，即给出一种开关指派的近似最佳方法，由此可设计出效果较好的近似算法法，由此可设计出效果较好的近似算法由于技术水平的提高，目前通讯卫星传送信息已允许一个发射站同时向多由于技术水平的提高，目前通讯卫星传送信息已允许一个发射站同时向多个接收站发送信息，当然，同时发送的信息条数具有某一上限，例如上限个接收站发送信息，当然，同时发送的信息条数具有某一上限，例如上限为为v1987年，年，J.L.Lewandowski和和C.L.Liu研究了如下更一般的问题：研究了如下更一般的问题： 给定一正整数给定一正整数v，（，（v为通讯卫星传送容量的总限止），求开关模式为通讯卫星传送容量的总限止），求开关模式M：：= ：：={ ; (0, 1) | , i= 1, …, m; ,i= 1, …, n, }的设计，要求所用的开关模式总数的设计，要求所用的开关模式总数量量r尽可能小，并使尽可能小，并使有解，其中有解，其中T为信息传送量矩阵（需满足一定要求），为信息传送量矩阵（需满足一定要求），ak为开关模式为开关模式Mk的使的使用时间。 他们设计了一个求解此问题的用时间他们设计了一个求解此问题的O(n5)算法，有兴趣的读者可直接阅算法，有兴趣的读者可直接阅读他们的论文读他们的论文 The End 。
