数模实验报告

1、袀重庆 交 通 大 学艿学生 实 验 报 告芆肁虿实验课程名称数学建模莈蚇开课实验室数学实验室螃蚂学院信息年级 09级专业班计科 3 班蒈螄学生 姓 名郭凯强蒅蒁 学号09060304薈膅开 课 时 间 2011至 2012学年第 1期羃芀蚈薆膂罿蚄成绩肁一 .小行螈实验指导教师羇杜绍洪星的轨道问题袈一天文学家要确定一颗小行星的轨道，他在轨道平面内建立以太阳为原点的直角坐标系，在两坐标轴上取天文测量单位 （一天文单位为地球到太阳的平均距离： 1.49597871011 m ）。在 5 个不同的时间对小行星作了5 次观察 ,测得轨道上5 个点的坐标数据如表5。膃表5坐标数据袄螀x1袈x2x3节x

2、4蕿x5薄羈 x 坐 标 羅 5.764羄 6.286薂 6.759肇 7.168莆 7.408蒂莁y1膇 y2螇 y3膄 y4膀 y5芇 y 坐 标 袄 0.648蚂 1.202衿 1.823莇 2.526芅 3.360莃由开普勒第一定律，小行星轨道为一椭圆。（注：椭圆的一般方程可表示为：羂 a1 x22a2 xya3 y 22a4 x2a5 y10蒇 1、请建立求椭圆轨道方程的数学模型。蚅 2、查找资料并求小行星的近日点距离和远日点距离以及椭圆轨道周长近似值。螁 3、画出小行星的轨迹图。蚀 解： (1).需要确定系数 ai ,i1,2,3, 4,5;蒇 利用已知的数据，不妨设xiyii1,

3、2,3,4,5;肆 欲确定系数 ai 等价于求解一个线性方程组：a1 x122a2 x1 y1a3 y122a4 x12a5 y11 0a1 x222a2 x2 y2a3 y222a4 x22a5 y210薃 a1x322a2 x3 y3a3 y322a4 x32a5 y310 可写成矩阵的形式： ax ba1 x422a2 x4 y4a3 y422a4 x42a5 y410a1 x522a2 x5 y5a3 y522a4 x52a5 y510x122x1 y1y122x12 y1a11x222x2 y2y222x22 y2a21葿 其中， a x22x y3y22x2 y3， xa3 ， b

4、13333x22x4yy22x2 y4a414444x522x5 y5y522x52 y5a5133.22377.47010.419911.52801.2960a1139.513815.11151.444812.57202.4040a21薇 45.684124.64333.323313.51803.6460a3151.380236.21276.380714.33605.0520a4154.878549.781811.289614.81606.7200a51蒇可见，解答上述问题就是对线性方程进行求解。羁 模型求解：利用 matlab 求解 ,程序代码如下：薂 %该函数用于求解椭圆方程的各系数蚆

5、function xishu=claua()薄 a= 5.764 0.648;蚃 6.286 1.202;芁 6.759 1.823;螆 7.168 2.526;羅 7.408 3.360莅 ;肀 x2=a(:,1).2;袆 m,n=size(a);莆 xy2=;袃 for i=1:m衿 xy2=xy2;2.*(a(i,1)*a(i,2);羆 end螇 y2=a(:,2).2;薅 doublex=2.*a(:,2);袂 doubley=2.*a(:,1);羆 b=-ones(m,1);羄 aa=x2,xy2,y2,doublex,doubley;肂 xishu=inv(aa)*b;蚁程序运行结

6、果：肆 xishu =莄 0.0507螄 -0.0351荿 0.0381蒀 -0.2265螅 0.1321膂所以轨道方程是：蒂 0.0507 x22*0.0351 xy0.0381y 22*0.2265 x2*0.1321 y102、3、 薀 查找资料并求小行星的近日点距离和远日点距离以及椭圆轨道周长近似值。a1a2a1a2a4d a2a3a5膆令 ca2a3a4a5a1羄 利用 matlab指令 eig(c) 指令求出 c 的特征值 10.080120. 0088 。利用 det(c) 和 det(d) 计算二阶和三阶行列式 |c|=7.0286e-004和|d|=-0.0011 。利用 公

7、 式 长 半 轴 add=20.1289 ， 焦 半 径=183.8359 短 半 轴 b2 c1 cca2b 2 =182.7306 , 得近日点 h=a-c=1.1053，远日点 h=2c+h=366.5665 。椭圆周长的近似值为l 3= 770.0552膁a bab2蚀 用 matlab 软件，求解上述数值。程序代码如下：薇 function a,b,c,root,det_c,det_d,near_sun,far_sun,l=get_abc()莂 xishu=claua();羀 c=xishu(1),xishu(2);虿 xishu(2),xishu(3)蚄 ;肄 d=xishu(1)

8、,xishu(2),xishu(4);蝿 xishu(2),xishu(3),xishu(5);蝿 xishu(4),xishu(5),xishu(1)肅 ;薂 root=eig(c);螂 det_c=det(c);衿 det_d=det(d);蒆 a=abs(det_d/(det_c*root(1);芄 b=abs(det_d/(det_c*root(2);薁 c=sqrt(a.2-b.2);罿 near_sun=a-c;袇 far_sun=a+c;螁 l=pi*(1.5*(a+b)-sqrt(a*b);4、5、 芀用运 matlab 软件画出小行星的轨迹图如下：聿肃相关程序代码如下：蒃肈

9、function draw_cir腿 a,b,c=get_abc;蒄 t=0:0.001:2*3.1415926;袁 x=a*cos(t);肁 y=b*sin(t);艿 plot(x,y);袅 title(- 行星轨迹图 -);薃 hold on;袀 grid on;艿二、陈酒出售的最佳时机问题芆某酒厂有批新酿的好酒， 如果现在就出售， 可得总收入 r050 万元（rmb ），如果窖藏起来待来日（第 n 年）按陈酒价格出售，第 n 年末可得总收入为：1n（ 1）肁 r r0 e6（万元）虿 当银行利率为 r 0.05 时，试分析这批好酒窖藏多少年后出售可使总收入的现值最大？莈 第一种方案：如果

10、现在出售这批好酒， 可得本金 50 万元。 由于银行利率为r 0.05 ，按照复利计算公式，第 n 年连本带利资金积累为：蚇 b( n)50(1 0.05) n（ 2）螃第二中方案：如果窖藏起来，待第n 年出售。原来的 50 万元到第 n 年时增值为：1蚂 r(n)n（ 3）50e6蒈 1、利用公式（ 2）和（ 3）分别计算出 16 年内采用两种方案， 50 万元增值的数目，并填写表 1 和表 2。螄蒅蒁表 1第一种方案薈第 1 年膅第 2 年羃第 3 年芀第 4 年蚈第 5 年薆第 6 年蚄第 7 年罿第 8 年螈 52.5羇 55.125膂 57.88肁 60.77袈 63.81膃 67.

11、00袄 70.36螀 73.81第 9 年薄第 10节第 11蕿第 12羈第 13羅第 14羄第 15薂第 16袈年年年年年年年肇 77.57莆 81.44蒂 85.52莁 89.79膇 94.28螇 99.00膄 103.95膀 109.14芇袄表 2第二种方案蚂第 1 年衿第 2 年莇第 3 年芅第 4 年莃第 5 年羂第 6 年蒇第 7 年蚅第 8 年螁 59.07蚀 63.29蒇 66.73肆 69.78薃 72.58葿 75.21薇 77.71蒇 80.11羁 第 9 年薂 10 年蚆第 11薄第 12蚃第 13芁第 14螆第 15羅第 16年年年年年年莅 82.44肀 84.70袆

12、 86.90莆 89.07袃 91.19衿 93.28羆 95.35螇 97.39薅比较表 1 和表 2 的数据。考虑如下问题（应简述原因） ：袂（1）如果酒厂希望在 2 年后投资扩建酒厂， 应选择哪一种方案使这批好酒所具有的价值发挥最大作用？羆（ 2）如果酒厂希望在6 年后将资金用作其它投资，应该选择哪一种方案？羄 解：（1）第二种，由表可知，在第二年的时候，第一种方案只有 55.125 而第二种有 63.29，故采用第二种方案。肂 （2）第二种，由表可知，六年之前都是第二种的获利多，故选用第二种。蚁运用 matlab 求解表一和表二的数据：肆 function bn,rn=plan12莄

13、r0=50;%本金螄 r=0.05;%银行利率荿 for i=1:16蒀 bn(i)=r0*(1+r)i;%保存方案一的数据螅 rn(i)=r0*exp(1/6*sqrt(i);%保存方案二的数据膂 end蒂 2、 现在酒厂有一笔资金，数额为 x （万元）将其存入银行，到第 n 年时增值为 r(n) （万元）。已知复利公式为：薀 r(n) x (n)(1 0.05) n（ 4）膆 其中 x (n) 为 r(n) 的现值。故 r(n) 的现值计算公式为：r(n)（ 5）羄 x (n)0.05) n(1膁将（ 3）代入（ 5），可得酒厂将这批好酒窖藏起来作为陈酒在第n 年后出售所得总收入的现值为：

14、1 n50e6（ 6）蚀 x (n)0.05) n(1薇利用公式（ 6），计算出 16 年内陈酒出售后总收入r( n) 的现值并填写表3：莂羀虿表3 陈酒出售后的现值蚄第 1 年肄第 2 年蝿第 3 年蝿第 4 年肅第 5 年薂第 6 年螂第 7 年衿第 8 年蒆 56.26芄 57.41薁 57.65罿 57.41袇 56.87螁 56.12芀 55.23聿 54.22第 9 年蒃第 10肈第 11腿第 12蒄第 13袁第 14肁第 15艿第 16肃年年年年年年年袅 53.14薃 52.00袀 50.81艿 49.60芆 48.36肁 47.11虿 45.86莈 44.61蚇根据表 3 中的

15、数据考虑如下问题（应简述原因） ：螃（1）如果酒厂打算将这批好酒出售所得收入用于 8 年后的另外投资， 应选择哪一年作为出售陈酒的最佳时期？蚂（ 2）如果综合考虑银行利率，将出售陈酒后所得总收入再存入银行，使得 8 年后资金增值最大，又应该作何选择？蒈 解：（1）第八年，由上表可知，到第八年出售可以增值。螄 （2）由于第 i 年出售酒后，要存入银行生利，故得n 年后的收入表达式：蒅 yx(1o.o5)8 n1n8蒁 由 matlab 求解可得薈for n=1:8膅 b(n)=50*exp(sqrt(n)/6)*(1+0.05)(8-n)羃 end芀运行结果如下：蚈 b =薆 56.2553 5

16、7.4058 57.6464 57.4087 56.8690 56.1221 55.2269 54.2231蚄 由上可知，第三年时，总收入最大为入银行，此时可使得 8 年后资金增值最大。57.6464，故应该三年后出售再存罿 3考虑银行利率按连续复利公式计算：螈 r(t)x (t )e0.05t ，或 x (t )r(t)e 0.05t（ 7）而酒厂将这批好酒窖藏到第n 年，作为陈酒出售总收入为：r(t)1 t（ 8）50e6结合公式（ 7）和（ 8），将 t 年后陈酒出售总收入的现值x 视为时间 t 的函数。试写出函数 x (t ) 的表达式，并利用求一元函数最大值的方法求出酒厂将这批好酒作

17、为陈酒出售的最佳时机。解：由 (7)、(8)公式得，x(t )50e1/6 t e 0.05 t ，应使得其取得最大值。求出 x(t) 的唯一驻点 t0=1/（ 144*0.052）=3 年所以陈酒出售的最佳时期为第三年。matlab 程序代码：%用于求解陈酒出售最佳时机和对应的最大收益function y,p=xt_max()f=inline(50.*exp(1/6.*sqrt(t)-0.05*t);x=f(1:16);t=1:16;plot(t,x);%画出收益随年份的变化图像hold on;f=inline(-50.*exp(1/6.*sqrt(t)-0.05*t);y,p=fminbn

18、d(f,1,16);p=-p;运行结果，图像如下示。并求得， y=2.7778，p=57.4498.由上述的求解结果以及函数图像，可知在第三年出售时为最佳时机。三、电力公司生产经营计划某电力公司经营两座发电站，发电站分别位于两个水库之上，位置如下图所示：水源 b水源 a水库 a发电站 a水库 b发电站 b已知发电站 a 可以将水库 a 的 1 万立方米的水转换为400 千度电能，发电站b 只能将 1 万立方米的转为 200 千度电能。发电站 a ，b 每个月的最大发电能力分别是 60000 千度， 35000 千度。每个月最多有50000 千度的电能够以200 元/千度的价格售出，多余的电能只

19、能够以140 元/千度的价格售出。水库a ，b 的其他有关数据如下（单位：万立方米） ：水库 a水库 b水库最大储水量20001500水源流入水量本月20040下月13015水库最小储水量1200800水库目前储水量1900850根据上述资料制定该电力公司本月和下月的生产经营计划?解： (1) 问题分析本题要解决的是水库所在何种情况下怎样调度水量，能使所放出的水来产生的效益。其实，问题本身就是一个线性规划问题，它是要求我们利用线性代数的有关知识来解决放多少水， 蓄多少水才能使水电公司收益最大。这样我们可以根据题中的数据，利用题中的一系列约束条件建立目标函数，即线性规划方程。最后利用 matla

20、b 软件对其求解，可得最优解，即为公司最大的收益时，各未知量的解。(2) 模型假设1、假如水资源流入水量是在每个月开始时发生的。2、假如水库中的水允许不发电而放走。3、假如不发生洪涝、干旱等灾害。4、假如电力公司发电系统均能正常运作。5、假如所有发生电量都能售出。6、假如售出不受价值规律的约束。7、假如水库水流一直处在动态过程。(3)符号说明1、 ma1：a 水库本月发电用水量。（单位：万 m3）2、 ma2：a 水库下月发电用水量。3、 mb1 ：b 水库本月发电用水量。4、 mb2 ：b 水库下月发电用水量。5、 na1： a 水库本月直接放走的水量。6、 na2： a 水库下月直接放走的

21、水量。7、 nb1： b 水库本月直接放走的水量。8、 nb2： b 水库下月直接放走的水量。9、 wa1：本月结束时 a 水库的水量。10、wa2：下月结束时 a 水库的水量。11、wb1：本月结束时 b 水库的水量。12、wb2：下月结束时 b 水库的水量。13、x1 ：本月以高价出售的电量。14、x2 ：下月以高价出售的电量。15、y1 ：本月以低价出售的电量。16、y2 ：下月以低价出售的电量。17、s：水库发电所获得的最大效益。 （单位：万元）(4) 模型建立根据线性规划模型，结合数学建模知识，我们可以得出如下模型，即线性优化方程组：maxs=200（x1+x2 ）+140（y1+y

22、2 ）约束条件有：1）每个月的发电量高于当月卖出的电量：400 ma1+200mb1 = x1+y1400 ma2+ 200mb2=x2+y22）水量永恒约束：ma1+ na1+wa1=1900+200mb1+ nb1+wb1=850+40+ ma1+ na1ma2+ na2 +wa2= wa1+130mb1+ nb2+wb2= wb1+15+ ma2+ na23）发电能力限制：400ma160000；400ma260000200 mb135000；200 mb2 350004）水库蓄水量限制：1200 wa12000；1200wa2 2000800wb11500；800 wb2 15005）

23、高价电量的限制：x1 50000；x2 50000(5)模型求解运用 matlab 软件求解，相关代码如下：f =0 0 0 0 0 0 0 0 0 0 0 0 -200 -200 -140 -140;a=;b=;aeq=400 0 200 0 0 0 0 0 0 0 0 0 -1 0 -1 00 400 0 200 0 0 0 0 0 0 0 0 0 -1 0 -11 0 0 0 1 0 0 0 1 0 0 0 0 0 0 0-1 0 1 0 -1 0 1 0 0 0 1 0 0 0 0 00 1 0 0 0 1 0 0 -1 1 0 0 0 0 0 00 -1 0 1 0 -1 0 1 0 0 -1 1 0 0 0 0;beq=0 0 2100 890 130 15;vlb=0;0;0;0;0;0;0;0;1200;1200;800;800;0;0;0;0;vub=150 150 175 175 150 150 175 175 2000 2000 1500 1500 50000 50000 45000 45000;x,fval=linprog(f,a,b,aeq,beq,vlb,vub);运行 matlab ，得出来的结果如下：x =1.0e+004