最新数学建模第三次作业.docx

精品文档
院系：数学学院
专业：信息与计算科学
年级：2014 级
学生姓名：王继禹
学号：201401050335
教师姓名：徐霞
6.6 习题
3.一个慢跑者在平面上沿着他喜欢的路径跑步，突然一只狗攻击他，这只狗以恒定速率跑向慢跑者，狗的运动方向始终指向慢跑者，计算并画出狗的轨迹。

解：
（1）模型分析建立:
狗的轨迹：在任意时刻，狗的速度向量都指向它的目标慢跑者。

假设 1：慢跑者在某路径上跑步，他的运动由两个函数 X(t)和 Y(t)描述。

假设 2：
当 t=0 时，狗是在点 (x0,y0)处，在时刻 t 时，它的位置是 (x(t),y(t)) 那么下列方程
成立：
222
(1)狗以恒定速率跑:X’+y’=w
(2)狗的速度向量平行于慢跑者与狗的位置的差向量：
将上述方程带入等式：,可得：
再将λ代入第二个方程，可得狗的轨迹的微分方程：
（2）程序及结果
dog 函数
[dog.m]
function[zs,isterminal,direction] = dog(t,z,flag)
global w; % w=speed of the dog
X=jogger(t);
h = X-z;
nh=norm(h);
if nargin<3 || isempty(flag)
zs=(w/nh)*h;
else
switch (flag)
case 'events'
zs = nh-1e-3;
isterminal = 1;
direction = 0;
otherwise
error(['Unknow flag:'flag]);
end
end
慢跑者的运动轨迹方程，水平向右
[jogger.m]
function s = jogger(t);
s = [8*t;0];
标记的函数
[cross.m]
function cross(Cx,Cy,v)
Kx = [Cx Cx Cx Cx-v Cx+v];
Ky = [Cy Cy+2.5*v Cy+1.5*v Cy+1.5*v Cy+1.5*v]
plot(Kx,Ky);
plot(Cx,Cy,'o' );
主程序：
静态显示
[main1.m]
global w
y0 = [60;70];
w=10;
options = odeset('RelTol',1e-5,'Events', 'on' ); [t,Y] = ode23('dog' ,[0,20],y0,options);
clf;hold on ;
axis([-10,100,-10,70]);
plot(Y(:,1),Y(:,2));
J=[];
for h=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)
hold off ;
动态显示
[main2.m]
global w;
y0=[60;70];
w=10;
options = odeset('RelTol',1e-5,'Events', 'on' ); [t,Y]=ode23('dog' ,[0,20],y0,options); J=[];
for h=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;hold on ;
axis([xmin-10 xmax ymin-10 ymax]);
title('The jogger and the Dog');
for h = 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', 'EraseMo de', 'none');
drawnow;
pause(0.1);
end
plot(J(:,1),J(:,2),':' );
p = max(size(Y));
cross(Y(p,1),Y(p,2),2)
hold off;
结果
t=12.2761812635281，在 12.27 秒后狗追上慢跑者。

慢跑者轨迹是椭圆轨迹
[jogger2.m]
function s=jogger2(t)
s=[10+20*cos(t) 20+15*sin(t)];
狗的微分方程
[dog.m]
function[zs,isterminal,direction] = dog(t,z,flag)
global w; % w=speed of the dog
X=jogger2(t);
h = X-z;
nh=norm(h);
if nargin<3 || isempty(flag)
zs=(w/nh)*h;
else
switch (flag)
case 'events'
zs = nh-1e-3;
isterminal = 1;
direction = 0;
otherwise
error(['Unknow flag:'flag]);
end
end
主程序
[main3.m]
global w;
y0=[60;70];
w=10;
options = odeset('RelTol',1e-5,'Events', 'on' ); [t,Y]=ode23('dog' ,[0,20],y0,options); J=[];
for h=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;hold on ;
axis([xmin-10 xmax ymin-10 ymax]);
title('The jogger and the Dog');
for h = 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', 'EraseMo de', 'none');
drawnow;
pause(0.1);
end
plot(J(:,1),J(:,2),':' );
p = max(size(Y));
cross(Y(p,1),Y(p,2),2)
hold off;
结果取 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.某人有n元，他每天一次物品，每次物品的品种很，或者一元的甲物品，
或者二元的乙物品，他花完n 元有多少不同的方式？
解： an 表示花完n 元的方案种数，
若n=1，只能甲，有一种方法，故a1=1，
若n=2，可以 2 个甲，或者 1 个乙或 1 个丙，即 a2 =3，当 n≥3 ，花
的方式由甲和乙丙的种数之和构成，
即 a n=a n-1+a n-2+a n-2=a n-1+2a n-2
当 n≥ 3 ， a n +a n-1=2（ a n-1+a n-2），
即{a n+1 +a n}是公比 q=2 的等比数列，首a2 +a1=1+3=4，
a n+1+a n=4?2n-1=2n+1，
n
∴a n+a n-1=2 ，
两式相减得a n+1-a n-1=2n+1-2-=2- ，（n≥ 2），
若n 是奇数， a n=2n-1+2n-3+⋯ +22+a1=(2n+1-1)/3
n-1 n-33n+1
若 n 是偶数， a n=2 +2+⋯ +2 +a2=(2+1)/3
．
7.6 习题
1.在化工生产中常常需要知道丙烷在各种温度T 和压力 P 下的导热系数K 。

下面是实验得
到的一组数据：
T /C68688787106106140140
9.798113.3249.007813.3559.791814.2779.656312.463
P /103
KPa
K0.08480.08970.07620.08070.06960.07530.06110.0651试求 T =99 C 和P=10.3x103KPa下的 K。

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

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

建立 M 文件：
function y=y_lagr1(x0,y0,x)
n=length(x0);m=length(x);
for i=1:m
z=x(i);
s=0.0;
for k=1:n
p=1.0;
for j=1:n
if j~=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。

t
4.用电压 V=10 伏的电池给电容器充电，电容器上t 时刻的电压为v(t)=V-(V-V0) e T,其中 v0是电容器的初始电压，是充电常数。

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

t/s0.51234579 V/伏 6.36 6.487.268.228.668.999.349.63解：
建立 M 文件
function f=dianya(x,t)
精品文档
f=10-(10-x(1))*exp(-t/x(2))
%x(1)=V0； x(2)=τ
主程序：
t=[0.5 1 2 3 4 5 7 9];
V=[6.36 6.48 7.26 8.22 8.66 8.99 9.43 9.63];
x0=[0.2,0.05];
x=lsqcurvefit('dianya',x0,t,V)%x(1)代表 V0， x(2)代表τ
f=dianya(x,t)
结果
x=5.5577 3.5002,即初始电压v0=5.5577 ，充电常数为τ =3.5002．精品文档。