数学模型课程设计.doc

数学建模实践数学建模课程设计（程序设计和论文）题目 1.求微分方程的数值解 2.傅立叶级数 3.确定罪犯藏身地点问题4.曲线拟合与回归分析 5.麦克劳林多项式 6.酒杯的生产 班级 / 学号 14140101/2011041401011 学 生 姓 名 黄中武 指 导 教 师 单锋 朱丽梅 沈阳航空航天大学课 程 设 计 任 务 书课 程 名 称 数学建模实践 院（系） 理学院 专业 信息与计算科学 班级 14140101 学号2011041401011 姓名 黄中武 课程设计题目1.求微分方程的数值解 2.傅立叶级数 3.确定罪犯藏身地点问题 4.曲线拟合与回归分析 5.麦克劳林多项式 6.酒杯的生产 课程设计时间: 2013年 6月 17日至 2013年 7月 4日 要求1、学习态度要认真，要积极参与课程设计，锻炼独立思考能力；2、严格遵守上机时间安排；3、按照MATLAB编程训练的任务要求来编写程序；4、根据任务来完成数学建模论文；5、报告书写格式要求按照沈阳航空航天大学“课程设计报告撰写规范”；7、报告上交时间：课程设计结时上交报告。8、严谨抄袭行为。课程设计的内容及要求：内容1. 求微分方程的数值解其中初始条件为。其他参数为计算在时刻处是函数值。2. 设函数是以为周期的函数，（1）编写表示函数的函数M文件y=fd(x)，并绘图；（2）求函数的函数傅立叶级数：（3）对任意的x（x可以是数组）和n编写三角多项式的函数M文件y=fly(x,n)（4）对任意的数组x和n，在同一平面内画出函数的图形，并进行比较。3. 确定罪犯藏身地点问题已知罪犯系列犯罪的地点坐标如下：PlaceDimensionLongitude153497.73N13230.83W2534828.74N13158.44W3535028.40N12932.95W4534838.71N14549.71W553494.08N13155.91W6532557.38N2153.14W7534759.94N14619.78W8533915.67N14646.64W9532740.55N21326.73W10534241.97N15226.05W11534721.45N14547.33W12534830.40N14016.68W13534920.17N13440.65W（1）将角度制的坐标转换为弧度制坐标；（2）将系列犯罪的地点（13个）的球面坐标转换为平面坐标；并画出这些点的图形（要求带适当宽度的网格，行宽1.25*103，列宽h=2.5*103）。（3）求这13个点之间距离的最大值和最小值；分别求这13个点相邻点横坐标和纵坐标之差绝对值的最大值和最小值； （4）对所得13个点的每个点，在其周围找16个点（可疑点），其中格宽为d=500。将这16个点与其余的12个犯罪的地点进行比较，若与某犯罪的地点距离小于等于d，则去掉该可疑点。求出所有剩余的点（构成可疑点集合）。（5）对于d=500:500:60000，重复上述过程，对于不同的d，求可疑点的个数，画出d与可疑点个数的图形，并确定d为何值使，可疑点个数最小。（6）对上述所得的d，类似问题4，求可疑点YD。（7）对于YD中每个点（如第j个点），计算该点到问题2中的13个点的距离，建立该点的衰减函数：，求的最大值及YD中对应的点罪犯所在的可能性最大的点，并将该点加在问题2的图形上。（8）对于d= (21:30)*500，重复上述过程，对于不同的d，求罪犯所在的可能性最大的点，并将这些点加在问题2的图形上。4. 曲线拟合与回归分析(1)以函数M文件形式编写用1次多项式进行曲线拟合的最小二乘法通用程序，并根据如下数据计算系数.(2) 再用回归分析法进行线性回归分析 33组的值序号123 4 5 6780120731808012512581.1901133.02731808012512581.1902129.63731808012512581.1903158.77731808012512581.1904145.32731808012512581.190512078.5961808012512581.190612075.451808012512581.190712090.4871808012512581.190812083.8481808012512581.190912073231.398012512581.1901012073198.488012512581.1901112073212.648012512581.1901212073190.558012512581.190131207318075.85712512581.190141207318065.95812512581.190151207318087.25812512581.190161207318097.82412512581.190171207318080150.7112581.190181207318080141.5812581.190191207318080132.3712581.190201207318080156.9312581.190211207318080125138.8881.190221207318080125131.2181.190231207318080125141.7181.190241207318080125149.2981.19025120731808012512560.5829026120731808012512570.9629027120731808012512564.8549028120731808012512575.5299029120731808012512581.1104.8430120731808012512581.1111.2231120731808012512581.198.09232120731808012512581.1120.44 对应的33组的值 序号1234560164.78140.87-144.25119.09135.44157.691165.81140.13-145.14118.63135.37160.762165.51140.25-144.92118.7135.33159.983167.93138.71-146.91117.72135.41166.814166.79139.45-145.92118.13135.41163.645164.94141.5-143.84118.43136.72157.226164.8141.13-144.07118.82136.02157.57165.59143.03-143.16117.24139.66156.598165.21142.28-143.49117.96137.98156.969167.43140.82-152.26129.58132.04153.610165.71140.82-147.08122.85134.21156.2311166.45140.82-149.33125.75133.28155.0912165.23140.85-145.82121.16134.75156.7713164.23140.73-144.18119.12135.57157.214163.04140.34-144.03119.31135.97156.3115165.54141.1-144.32118.84135.06158.2616166.88141.4-144.34118.67134.67159.2817164.07143.03-140.97118.75133.75158.8318164.27142.29-142.15118.85134.27158.3719164.57141.44-143.3119134.88158.0120163.89143.61-140.25118.64133.28159.1221166.35139.29-144.2119.1136.33157.5922165.54140.14-144.19119.09135.81157.6723166.75138.95-144.17119.15136.55157.5924167.69138.07-144.14119.19137.11157.6525162.21141.21-144.13116.03135.5154.2626163.54141-144.16117.56135.44155.9327162.7141.14-144.21116.74135.4154.8828164.06140.94-144.18118.24135.4156.6829164.66142.27-147.2120.21135.28157.6530164.7142.94-148.45120.68135.16157.6331164.67141.56-145.88119.68135.29157.6132164.69143.84-150.34121.34135.12157.645.（1）编写表示函数的函数M文件y=fd(x)，并绘图；（2）求函数（3）对任意的x（x可以是数组）和n编写多项式的函数M文件y=fly(x,n)编写对任意固定的n计算多项式函数值的函数M文件（4）对任意的数组x和n，在同一平面内画出函数的图形，并进行比较。（5）用的近似值，要求误差不超过6. Production of drinking glassesThe main activity of a company in northern France is the production of drinking glasses. It currently sells six different types (V1 to V6), that are produced in batches of 1000 glasses, and wishes to plan its production for the next 12 weeks. The batches may be incomplete (fewer than 1000 glasses). The demand in thousands for the 12 coming weeks and for every glass type is given in the following table.Table : Demands for the planning period (batches of 1000 glasses)For every glass type the initial stock is known, as well as the required final stock level (in thousands). Per batch of every glass type, the production and storage costs in BC are given, together with the required working time for workers and machines (in hours), and the required storage space (measured in numbers of trays).The number of working hours of the personnel is limited to 390 hours per week, and the machines have a weekly capacity of 850 hours. Storage space for up to 1000 trays is available. Which quantities of the different glass types need to be produced in every period to minimize the total cost of production and storage?Table 2: Data for the six glass types指导教师 年 月 日负责教师 年 月 日学生签字 年 月 日沈阳航空航天大学课 程 设 计 成 绩 评 定 单课 程 名 称 数学建模实践 院（系） 理学院 专业 信息与计算科学 课程设计题目1.求微分方程的数值解 2.傅立叶级数 3.确定罪犯藏身地点问题 4.曲线拟合与回归分析 5.麦克劳林多项式 6.酒杯的生产 学号 2011041401011 姓名 黄中武 指导教师评语：课程设计成绩 指导教师签字 年 月 日目 录第一章 编程任务10摘要10正文12问题一：求微分方程的数值解121.1问题的重述121.2问题的分析121.3程序框图131.4问题的结果14问题二：傅立叶级数142.1问题重述142.2问题分析142.3程序框图152.4问题的结果172.4.1问题一的图象172.4.2问题二的傅立叶级数证明172.4.3问题三的图像182.4.4问题四的图像19问题三：确定罪犯藏身地点193.1问题重述193.2问题分析203.3问题的程序框图213.4问题求解253.4.1问题一的图象253.4.2问题二的结果263.4.3问题三的结果263.4.4问题四的结果273.4.5问题五的图像273.4.6问题六的图像(部分)283.4.7问题七的图像293.4.8问题八的图像29问题四：曲线拟合与回归分析304.1问题重述304.2问题分析324.3程序框图324.4问题求解334.4.1问题一的求解334.4.2问题二的求解33问题五：麦克劳林多项式展开345.1问题重述345.2问题分析345.3程序框图355.4问题的求解375.41问题一的图像375.42问题二的麦克劳林多项式375.43问题三的结果385.44问题四的结果385.45问题五的结果39第二章 建模任务40问题六：Production of drinking glasses406.1问题重述406.2问题分析406.2问题假设和符号说明426.3模型建立436.4模型求解446.5结果分析与检验456.6模型的的优缺点46参考文献48第一章 编程任务摘要对于问题一，求解微分方程组，建立M文件，将题中给出的微分方程组及其所含参数值输入，通过调用函数ode45函数（4-5阶runge-kutta算法，中精度）来实现问题的求解，在命令窗口输入所求范围及初值，从而求出结果。对于问题二，求函数的傅立叶级数时，首先判断其是否为周期函数，若不是周期函数，必须先进行周期延拓，但在具体计算过程中，并非用到延拓后的函数。在问题（1）中，通过判断x所在的周期来作出图形；问题（2）在求函数傅立叶变换时，根据在一个函数周期中的函数表达式求出傅立叶级数的系数；问题（3）将问题（2）所求出的傅立叶级数的系数画出图形；问题（4）将问题（1）和问题（2）所作图形结合在一起，进行比较。对于问题三，对于问题(1)，我们把度，分，秒独立开来，通过角度之间的换算，最终全部以度作为最后的单位，再将结果除以180/pi，从而将角度值转换为弧度制。对于问题(2)，试题中给我们十三组坐标，我们通过建立M文件，将这十三组坐标列出，通过plot函数输出。对于问题(3)，求最小距离和最大距离。我们首先建立M文件，将十三组数据输出，通过for循环，将所有循环出来，再根据两点之间的距离公式算出距离，最后输出最小距离和最大距离。对于问题(4)，首先建立矩阵存储13给点的坐标，再以中央点为中心构建4*4的小区域，利用循环语句扣除不必要的点，将剩下的点进行距离上的比较，求得记结果。对于问题(5)，是在上一小题的基础上求解，在此，我们通过距离d作为循环量，通过不同的d来画出与可疑点的图形，并在此确定一个d，是可以点数目最少。对于问题(6)，对于问题(5)中的d，按照问题(4)的方法，得到结果。对于问题(7)，我们按照两点距离公式，计算问题(6)中YD和问题二中点的距离，并通过建立衰减函数，求出衰减函数的最大值和罪犯可能的所在点。对于问题(8)，在问题(7)的基础上，对不同的d，以距离d作为循环量，按照(7)的方法，算出罪犯可能的所在点。对于问题四，本题考查的是根据多组数据进行最小二乘法曲线拟合，在求取系数时，我们采用高等代数上关于最小二乘法的求解方法，通过矩阵之间的运算求出系数。对于问题五，根据麦克劳林的推导公式，我们可以用matlab求出它的次麦克劳林展式并将其以表达式的形式输出出来。然后，编写M文件求对于给定的次数和求出的值。最后对于第三小问，用plot命令将函数的值与其的麦克劳林展式在区间上进行比较。关键词：微分方程数值解；傅立叶级数；角度与弧度；最小二乘法；正文问题一：求微分方程的数值解1.1问题的重述其中初始条件为。其他参数为计算在时刻处是函数值。1.2问题的分析求解微分方程组，建立M文件，将题中给出的微分方程组及其所含参数值输入，通过调用函数ode45函数（4-5阶runge-kutta算法，中精度）来实现问题的求解，在命令窗口输入所求范围及初值，从而求出结果。1.3程序框图建立M文件输入函数表达式输入变量值在命令窗口输入初值输出结果开始结束图1.3 问题一的程序框图1.4问题的结果图1.4 问题一的结果问题二：傅立叶级数2.1问题重述设函数是以为周期的函数，（1）编写表示函数的函数M文件y=fd(x)，其中可以是任意实数，也可以是任意数组，并在任意区间上绘出函数图形；（2）推倒函数的函数傅立叶级数：(3)对任意的x（x可以是数组）和n编写三角多项式的函数M文件y=fly(x,n)（4）对任意的数组x和n，在任意区间上在同一平面内画出函数的图形，并进行比较。2.2问题分析对于问题一，通过建立M文件，输入自变量的范围和函数的表达式，通过plot函数画图。对于问题二，对于傅立叶基数的证明，可参照数学分析上的证明。对于问题三，该函数为周期函数，其图像为周期图像，首先建立M文件，输入自变量的范围，通过for语句的循环，来实现傅立叶级数程序的运行，最后通过plot函数，来实现图形的绘出。对于问题四，建立M文件，在M文件中分别写上两个函数图像的运行程序，通过hold on来实现统一窗口绘出两个函数的图像。 2.3程序框图开始输入x，i=1m=lenth(x)x(i)-(2*m+1)*pini=i+1否是输入plot(x,y)结束图2.3-1 问题一y=fd(x)程序框图图2.3-2 问题三y=fly(x,n)程序框图开始调用函数fd（x，n）输入x，n调用函数调用plot函数结束图2.3-3 同时绘制与其傅里叶级数图像的程序框图2.4问题的结果2.4.1问题一的图象图2.3.1 fd(x)的函数图像2.4.2问题二的傅立叶级数证明求函数的函数傅立叶级数：证明：函数展开成傅里叶级数，设是定义在上的有界周期函数，则能展开成三角级数 综上求得结果为：2.4.3问题三的图像图2.3.3 y=fly(x,n)函数图像 2.4.4问题四的图像图2.3.4 20阶傅里叶图像与原图像比较问题三：确定罪犯藏身地点3.1问题重述已知罪犯系列犯罪的地点坐标如下：PlaceDimensionLongitude153497.73N13230.83W2534828.74N13158.44W3535028.40N12932.95W4534838.71N14549.71W553494.08N13155.91W6532557.38N2153.14W7534759.94N14619.78W8533915.67N14646.64W9532740.55N21326.73W10534241.97N15226.05W11534721.45N14547.33W12534830.40N14016.68W13534920.17N13440.65W（1）将角度制的坐标转换为弧度制坐标；（2）将系列犯罪的地点（13个）的球面坐标转换为平面坐标；并画出这些点的图形（要求带适当宽度的网格，行宽1.25*103，列宽h=2.5*103）。（3）求这13个点之间距离的最大值和最小值；分别求这13个点相邻点横坐标和纵坐标之差绝对值的最大值和最小值；（4）对所得13个点的每个点，在其周围找16个点（可疑点），其中格宽为d=500。将这16个点与其余的12个犯罪的地点进行比较，若与某犯罪的地点距离小于等于d，则去掉该可疑点。求出所有剩余的点（构成可疑点集合）。（5）对于d=500:500:60000，重复上述过程，对于不同的d，求可疑点的个数，画出d与可疑点个数的图形，并确定d为何值使，可疑点个数最小。（6）对上述所得的d，类似问题4，求可疑点YD。（7）对于YD中每个点（如第j个点），计算该点到问题2中的13个点的距离，建立该点的衰减函数：，求的最大值及YD中对应的点罪犯所在的可能性最大的点，并将该点加在问题2的图形上。（8）对于d= (21:30)*500，重复上述过程，对于不同的d，求罪犯所在的可能性最大的点，并将这些点加在问题2的图形上。3.2问题分析对于问题(1)，我们把度，分，秒独立开来，通过角度之间的换算，最终全部以度作为最后的单位，再将结果除以180/pi，从而将角度值转换为弧度制。对于问题(2)，试题中给我们十三组坐标，我们通过建立M文件，将这十三组坐标列出，通过plot函数输出。对于问题(3)，求最小距离和最大距离。我们首先建立M文件，将十三组数据输出，通过for循环，将所有循环出来，再根据两点之间的距离公式算出距离，最后输出最小距离和最大距离。对于问题(4)，首先建立矩阵存储13给点的坐标，再以中央点为中心构建4*4的小区域，利用循环语句扣除不必要的点，将剩下的点进行距离上的比较，求得记结果。对于问题(5)，是在上一小题的基础上求解，在此，我们通过距离d作为循环量，通过不同的d来画出与可疑点的图形，并在此确定一个d，是可以点数目最少。对于问题(6)，对于问题(5)中的d，按照问题(4)的方法，得到结果。对于问题(7)，我们按照两点距离公式，计算问题(6)中YD和问题二中点的距离，并通过建立衰减函数，求出衰减函数的最大值和罪犯可能的所在点。对于问题(8)，在问题(7)的基础上，对不同的d，以距离d作为循环量，按照(7)的方法，算出罪犯可能的所在点。3.3问题的程序框图开始输入经纬度结束经纬度数*pi/180输出弧度制坐标in求坐标的距离输出最小距离否 是i=i+1图3.3-3 问题三的最值流开始输入坐标x,y结束输出符合条件的剩余点在图中做出犯罪点相邻的16个点做出的点与其余12个犯罪点比较除去距离小于d的点图3.3-4 问题四的距离上的比较流程图结束开始求出所对应的满足题意的可疑点计算这些可疑点到原先13个点的距离将距离代入衰减函数中求出衰减函数最大值所对应的可疑点，即为罪犯所在的可能性最大的点d=17d= (21:30)*500图3.3-5 求解总个数最小时所对应d的流程图开始ts=0.01:0.01:10导入myfun3_4结束直接将d的值输入即可图3.3-6 求可疑点流程图开始ts=0.01:0.01:10找出所求点并标注在图上结束把第二问的图画上图3.3-7 画图求点流程附加说明：由于问题(5)和问题(6)算法一样，所以流程图可以共用一个3.4问题求解3.4.1问题一的图象图3.4.1 问题一的结果3.4.2问题二的结果问题二的图象：图3.4.2 问题二的图像3.4.3问题三的结果图3.4.3 问题三的图像3.4.4问题四的结果图3.4.4 问题四的结果3.4.5问题五的图像图3.4.5 问题五的图像3.4.6问题六的图像(部分)图3.4.6 问题六的结果3.4.7问题七的图像图3.4.7 平面坐标上的可疑点3.4.8问题八的图像图3.4.8 罪犯所在的可能性最大的点问题四：曲线拟合与回归分析4.1问题重述(1)以函数M文件形式编写用1次多项式进行曲线拟合的最小二乘法通用程序，并根据如下数据计算系数.(2)再用回归分析法进行线性回归分析。 33组的值序号123 4 5 6780120731808012512581.1901133.02731808012512581.1902129.63731808012512581.1903158.77731808012512581.1904145.32731808012512581.190512078.5961808012512581.190612075.451808012512581.190712090.4871808012512581.190812083.8481808012512581.190912073231.398012512581.1901012073198.488012512581.1901112073212.648012512581.1901212073190.558012512581.190131207318075.85712512581.190141207318065.95812512581.190151207318087.25812512581.190161207318097.82412512581.190171207318080150.7112581.190181207318080141.5812581.190191207318080132.3712581.190201207318080156.9312581.190211207318080125138.8881.190221207318080125131.2181.190231207318080125141.7181.190241207318080125149.2981.19025120731808012512560.5829026120731808012512570.9629027120731808012512564.8549028120731808012512575.5299029120731808012512581.1104.8430120731808012512581.1111.2231120731808012512581.198.09232120731808012512581.1120.44 对应的33组的值 序号1234560164.78140.87-144.25119.09135.44157.691165.81140.13-145.14118.63135.37160.762165.51140.25-144.92118.7135.33159.983167.93138.71-146.91117.72135.41166.814166.79139.45-145.92118.13135.41163.645164.94141.5-143.84118.43136.72157.226164.8141.13-144.07118.82136.02157.57165.59143.03-143.16117.24139.66156.598165.21142.28-143.49117.96137.98156.969167.43140.82-152.26129.58132.04153.610165.71140.82-147.08122.85134.21156.2311166.45140.82-149.33125.75133.28155.0912165.23140.85-145.82121.16134.75156.7713164.23140.73-144.18119.12135.57157.214163.04140.34-144.03119.31135.97156.3115165.54141.1-144.32118.84135.06158.2616166.88141.4-144.34118.67134.67159.2817164.07143.03-140.97118.75133.75158.8318164.27142.29-142.15118.85134.27158.3719164.57141.44-143.3119134.88158.0120163.89143.61-140.25118.64133.28159.1221166.35139.29-144.2119.1136.33157.5922165.54140.14-144.19119.09135.81157.6723166.75138.95-144.17119.15136.55157.5924167.69138.07-144.14119.19137.11157.6525162.21141.21-144.13116.03135.5154.2626163.54141-144.16117.56135.44155.9327162.7141.14-144.21116.74135.4154.8828164.06140.94-144.18118.24135.4156.6829164.66142.27-147.2120.21135.28157.6530164.7142.94-148.45120.68135.16157.6331164.67141.56-145.88119.68135.29157.6132164.69143.84-150.34121.34135.12157.644.2问题分析本题考查的是根据多组数据进行最小二乘法曲线拟合，在求取系数时，我们采用高等代数上关于最小二乘法的求解方法，通过矩阵之间的运算求出系数。4.3程序框图开始输入矩阵A，B，i=1i6b=B(:,i); z=A*b; y=inv(A*A); a(:,i)=y*z;输出系数aij结束图4.3 最小二乘法曲线拟合的流程图4.4问题求解4.4.1问题一的求解 图4.4.1 问题一的结果4.4.2问题二的求解图4.4.2 问题二的结果问题五：麦克劳林多项式展开5.1问题重述（1）编写表示函数的函数M文件y=fd(x)，并绘图；（2）求函数（3）对任意的x（x可以是数组）和n编写多项式的函数M文件y=fly(x,n)编写对任意固定的n计算多项式函数值的函数M文件（4）对任意的数组x和n，在同一平面内画出函数的图形，并进行比较。（5）用的近似值，要求误差不超过。5.2问题分析针对本题，我们需要将转化为，在此，我们对和分别进行麦克劳林展开为和，则展开后可表示为，编写M文件求出其表达式形式即可。对于得到的表达式结果用利用字符表达式带入值很容易求出结果。在第三问中，输入的值将其带入和中，并求出它们的差值。为了方便画图，将程序做一修改，使其可以输入一个向量，可以得到一个结果向量。然后对其在区间画图进行比较。5.3程序框图开始输入x，y输出plot(x,y)结束图5.3.1 问题一的函数流程图开始输入n，x，s1=0，s2=0i=1ins1=s1+(-1)(i-1)*(xi)/is2=s2+(-1)*(xi)/i输出s=s2-s1结束否是图5.3.2 麦克劳林函数流程图开始输入x，mi=