数学建模第三次作业.docx

数学建模第三次作业.docx

《数学建模第三次作业.docx》由会员分享，可在线阅读，更多相关《数学建模第三次作业.docx（13页珍藏版）》请在冰豆网上搜索。

数学建模第三次作业

院系：

数学学院

专业：

信息与计算科学

年级：

2014级

学生姓名：

王继禹

学号：

201401050335

教师姓名：

徐霞

6.6习题

3.一个慢跑者在平面上沿着他喜欢的路径跑步，突然一只狗攻击他，这只狗以恒定速率跑向慢跑者，狗的运动方向始终指向慢跑者，计算并画出狗的轨迹。

解：

（1）模型分析建立:

狗的轨迹：

在任意时刻，狗的速度向量都指向它的目标慢跑者。

假设1：

慢跑者在某路径上跑步，他的运动由两个函数X（t）和Y（t）描述。

假设2：

当t=0时，狗是在点（x0,y0）处，在时刻t时，它的位置是（x（t）,y（t））

那么下列方程成立：

（1）狗以恒定速率跑:

X’2+y’2=w2

（2）狗的速度向量平行于慢跑者与狗的位置的差向量：

将上述方程带入等式：

可得：

再将λ代入第二个方程，可得狗的轨迹的微分方程：

（2）程序及结果

dog函数

[dog.m]

function[zs,isterminal,direction]=dog（t,z,flag）

globalw;%w=speedofthedog

X=jogger（t）;

h=X-z;

nh=norm（h）;

ifnargin<3||isempty（flag）

zs=（w/nh）*h;

else

switch（flag）

case'events'

zs=nh-1e-3;

isterminal=1;

direction=0;

otherwise

error（['Unknowflag:

'flag]）;

end

end

慢跑者的运动轨迹方程，水平向右

[jogger.m]

functions=jogger（t）;

s=[8*t;0];

标记的函数

[cross.m]

functioncross（Cx,Cy,v）

Kx=[CxCxCxCx-vCx+v];

Ky=[CyCy+2.5*vCy+1.5*vCy+1.5*vCy+1.5*v]

plot（Kx,Ky）;

plot（Cx,Cy,'o'）;

主程序：

静态显示

[main1.m]

globalw

y0=[60;70];

w=10;

options=odeset（'RelTol',1e-5,'Events','on'）;

[t,Y]=ode23（'dog',[0,20],y0,options）;

clf;holdon;

axis（[-10,100,-10,70]）;

plot（Y（:

1）,Y（:

2））;

J=[];

forh=1:

length（t）,

w=jogger（t（h））;

J=[J;w'];

end

plot（J（:

1）,J（:

2）,':

'）;

p=max（size（Y））;

cross（Y（p,1）,Y（p,2）,2）

holdoff;

动态显示

[main2.m]

globalw;

y0=[60;70];

w=10;

options=odeset（'RelTol',1e-5,'Events','on'）;

[t,Y]=ode23（'dog',[0,20],y0,options）;J=[];

forh=1:

length（t）;

w=jogger（t（h））;

J=[J;w'];

end

xmin=min（min（Y（:

1））,min（J（:

1）））;

xmax=max（max（Y（:

1））,max（J（:

1）））;

ymin=min（min（Y（:

2））,min（J（:

2）））;

ymax=max（max（Y（:

2））,max（J（:

2）））;

clf;holdon;

axis（[xmin-10xmaxymin-10ymax]）;

title（'ThejoggerandtheDog'）;

forh=1:

length（t）-1,

plot（[Y（h,1）,Y（h+1,1）],[Y（h,2）,Y（h+1,2）],'-','Color','red','EraseMode','none'）;plot（[J（h,1）,J（h+1,1）],[J（h,2）,J（h+1,2）],'-','Color','green','EraseMode','none'）;

drawnow;

pause（0.1）;

end

plot（J（:

1）,J（:

2）,':

'）;

p=max（size（Y））;

cross（Y（p,1）,Y（p,2）,2）

holdoff;

结果

t=12.2761812635281，在12.27秒后狗追上慢跑者。

慢跑者轨迹是椭圆轨迹

[jogger2.m]

functions=jogger2（t）

s=[10+20*cos（t）20+15*sin（t）];

狗的微分方程

[dog.m]

function[zs,isterminal,direction]=dog（t,z,flag）

globalw;%w=speedofthedog

X=jogger2（t）;

h=X-z;

nh=norm（h）;

ifnargin<3||isempty（flag）

zs=（w/nh）*h;

else

switch（flag）

case'events'

zs=nh-1e-3;

isterminal=1;

direction=0;

otherwise

error（['Unknowflag:

'flag]）;

end

end

主程序

[main3.m]

globalw;

y0=[60;70];

w=10;

options=odeset（'RelTol',1e-5,'Events','on'）;

[t,Y]=ode23（'dog',[0,20],y0,options）;J=[];

forh=1:

length（t）;

w=jogger2（t（h））;

J=[J;w'];

end

xmin=min（min（Y（:

1））,min（J（:

1）））;

xmax=max（max（Y（:

1））,max（J（:

1）））;

ymin=min（min（Y（:

2））,min（J（:

2）））;

ymax=max（max（Y（:

2））,max（J（:

2）））;

clf;holdon;

axis（[xmin-10xmaxymin-10ymax]）;

title（'ThejoggerandtheDog'）;

forh=1:

length（t）-1,

plot（[Y（h,1）,Y（h+1,1）],[Y（h,2）,Y（h+1,2）],'-','Color','red','EraseMode','none'）;plot（[J（h,1）,J（h+1,1）],[J（h,2）,J（h+1,2）],'-','Color','green','EraseMode','none'）;

drawnow;

pause（0.1）;

end

plot（J（:

1）,J（:

2）,':

'）;

p=max（size（Y））;

cross（Y（p,1）,Y（p,2）,2）

holdoff;

结果取w=25有

t=4.017776368842910，经过4秒左右狗追上慢跑者。

8.平面上有n（n>=2）个圆，任何两个圆都相交但无3个圆共点。

试问n个圆把平面划分成多少个不连通的区域？

解：

∵一个圆将平面分为2份

两个圆相交将平面分为4=2+2份，

三个圆相交将平面分为8=2+2+4份，

四个圆相交将平面分为14=2+2+4+6份，

…

平面内n个圆，其中每两个圆都相交于两点，且任意三个圆不相交于同一点，

则该n个圆分平面区域数f（n）=2+（n-1）n=n2-n+2

证明：

（1）当n=1时，一个圆把平面分成两个区域，而12-1+2=2，命题成立．

（2）假设n=k（k≥1）时，命题成立，即k个圆把平面分成k2-k+2个区域．

当n=k+1时，第k+1个圆与原有的k个圆有2k个交点，这些交点把第k+1个圆分成了2k段弧，

而其中的每一段弧都把它所在的区域分成了两部分，因此增加了2k个区域，

共有k2-k+2+2k=（k+1）2-（k+1）+2个区域．

∴n=k+1时，命题也成立．

由

（1）、

（2）知，对任意的n∈N*，命题都成立．

9.某人有

元钱，他每天买一次物品，每次买物品的品种很单调，或者买一元钱的甲物品，或者买二元钱的乙物品，问他花完这

元钱有多少不同的方式？

解：

设an表示花完这n元钱的方案种数，

若n=1，则只能买甲，有一种方法，故a1=1，

若n=2，则可以买2个甲，或者1个乙或1个丙，即a2=3，

当n≥3时，花钱的方式由购买甲和购买乙购买丙的种数之和构成，

即an=an-1+an-2+an-2=an-1+2an-2

则当n≥3时，an+an-1=2（an-1+an-2），

即{an+1+an}是公比q=2的等比数列，首项为a2+a1=1+3=4，

则an+1+an=4•2n-1=2n+1，

∴an+an-1=2n，

两式相减得an+1-an-1=2n+1-2=2，（n≥2），

若n是奇数，an=2n-1+2n-3+…+22+a1=（2n+1-1）/3

若n是偶数，an=2n-1+2n-3+…+23+a2=（2n+1+1）/3．

7.6习题

1.在化工生产中常常需要知道丙烷在各种温度

和压力

下的导热系数

。

下面是实验得到的一组数据：

/

68

68

87

87

106

106

140

140

/

KPa

9.7981

13.324

9.0078

13.355

9.7918

14.277

9.6563

12.463

K

0.0848

0.0897

0.0762

0.0807

0.0696

0.0753

0.0611

0.0651

试求

=99

和

=10.3x

KPa下的K。

解：

找出温度T相等时，导热系数K与压力P的关系。

由于在不同温度时，仅给出两个K、P的值，因此采用线性近似，把K、P看作是线性关系。

建立M文件：

functiony=y_lagr1（x0,y0,x）

n=length（x0）;m=length（x）;

fori=1:

m

z=x（i）;

s=0.0;

fork=1:

n

p=1.0;

forj=1:

n

ifj~=k

p=p*（z-x0（j））/（x0（k）-x0（j））;

end

end

s=p*y0（k）+s;

end

y（i）=s;

end

主程序：

p1=[9.7981,13.324];k1=[0.0848,0.0897];%T=68℃

p2=[9.0078,13.355];k2=[0.0762,0.0807];%T=87℃

p3=[9.7918,14.277];k3=[0.0696,0.0753];%T=106℃

p4=[9.6563,12.463];k4=[0.0611,0.0651];%T=140℃

a2=polyfit（p2,k2,1）;a3=polyfit（p3,k3,1）;

x1=polyval（a2,10.3）;x2=polyval（a3,10.3）;

%x1，x2分别是P=10.3*10^3kPa下87℃和106℃的K值

plot（10.3,x1,'k+',10.3,x2,'k+',p1,k1,p2,k2,p3,k3,p4,k4）

xlabel（'丙烷压力P'）

ylabel（'丙烷导热系数K'）

title（'在不同温度下丙烷导热系数与压力的关系图'）

gtext（'T=68℃'）,gtext（'T=87℃'）,gtext（'T=106℃'）,gtext（'T=140℃'）

运行后图中所标点为P=10.3*10^3kPa时，T=87℃和T=106℃对应的导热系数K值。

在T=87℃和T=106℃之间仍采用线性近似来求T=99℃时的导热系数K。

程序如下：

x=[87,106];

y=[x1,x2];

a=polyfit（x,y,1）;

z=polyval（a,99）

z=0.0729

plot（99,z,'k+',x,y）

grid

xlabel（'丙烷温度T'）

ylabel（'丙烷导热系数K'）

title（'压力P=10.3*10^3kPa时丙烷导热系数与温度的关系'）

运行结果:

T=99℃、P=10.3*10^3KPa时K=0.0729。

4.用电压V=10伏的电池给电容器充电，电容器上t时刻的电压为v（t）=V-（V-V0）

其中v0是电容器的初始电压，

是充电常数。

试由下面一组t，v数据确定V0和

。

t/s

0.5

1

2

3

4

5

7

9

V/伏

6.36

6.48

7.26

8.22

8.66

8.99

9.34

9.63

解：

建立M文件

functionf=dianya（x,t）

f=10-（10-x

（1））*exp（-t/x

（2））

%x

（1）=V0；x

（2）=τ

主程序：

t=[0.51234579];

V=[6.366.487.268.228.668.999.439.63];

x0=[0.2,0.05];

x=lsqcurvefit（'dianya',x0,t,V）%x

（1）代表V0，x

（2）代表τ

f=dianya（x,t）

结果

x=5.55773.5002,即初始电压v0=5.5577，充电常数为τ=3.5002．