大学数学实验清华.doc

Exp01 数学建模初步
人口指数增长模型p1
人口阻滞增长模型p2
Exp02 差分方程和数值微分
濒危物种(Florida沙丘鹤)自然演变和人工孵化p5
一年生植物的繁殖p6
汽车租赁公司的运营p7
按年龄分组的种群增长p8
离散形式的阻滞增长模型（Logistic模型）p10
寄主—寄生模型p12
Exp03 插值与数值积分

Exp04 常微分方程数值解

Exp05 线性代数方程组的数值解

Exp06 非线性方程求解

Exp07 无约束优化

Exp08 约束优化

Exp09 整数规划

Exp10 数据的统计与分析

Exp11 统计推断

Exp12 回归分析

人口指数增长模型
%Exp01m 文件名：popEXPdemo.m
% 1790-1900美国人口数据
x=[3.9 5.3 7.2 9.6 12.9 17.1 23.2 ...
31.4 38.6 50.2 62.9 76.0];
t=0:length(x)-1;
r=0.2743;x0=4.1884;
plot(t,x,'b+'); hold on;
fplot('popexpfun',[0,length(x)],'r',[],[],r,x0)
xlabel('t');
ylabel('x');
hold off;
pause
% 加上1910-2000美国人口数据
x1=[92.0 106.5 123.2 131.7 150.7 179.3 204.0 226.5 251.4 281.4];
x=[x,x1];
t=0:length(x)-1;
r=0.2022;x0=6.0450;
plot(t,x,'b+'); hold on;
fplot('popexpfun',[0,length(x)+3],'r',[],[],r,x0)
xlabel('t','fontsize',20,'color','red');
ylabel('x','fontsize',20,'color','b');
hold off;
命令式m文件：popEXPdemo.m
图1：

fplot, 函数popexpfun.m

图2：

fplot, 函数popexpfun.m
%Exp01m 文件名：PopExpFun.m
function x = PopExpFun(t,r,x0)
x=x0*exp(r*t);
函数式m文件：PopExpFun.m

输出图1：指数增长模型拟合图形（1790年至1900年）

输出图2：指数增长模型模型拟合图形（1790年至2000年）

人口阻滞增长模型（Logistic模型）
%Exp01m 文件名：popdemo.m
% Logistics Model

fplot('popinc',[0,1]);
xlabel('x','fontsize',20);
ylabel('dx/dt','fontsize',20);
pause;

r=0.29;
xm=1;
x0=0.05;
fplot('popcal',[0,20],[],[],[],r,x0,xm);
xlabel('t','fontsize',20);
ylabel('x','fontsize',20);
pause

% 1790-1900美国人口数据
x=[3.9 5.3 7.2 9.6 12.9 17.1 23.2 31.4 38.6 50.2 ...
62.9 76.0 92.0 106.5 123.2 131.7 150.7 179.3 ...
204.0 226.5 251.4 281.4];
t=0:length(x)-1;
r=0.2557;x0=x(1);xm=392.0886;
plot(t,x,'b+'); hold on;
fplot('popcal',[0,length(x)+3],'r',[],[],r,x0,xm)
xlabel('t','fontsize',20,'color','red');
ylabel('x','fontsize',20,'color','b');
hold off;
命令式m文件：popdemo.m
阻滞增长模型（Logistic模型）
图1：输出Logistic模型dx/dt~x曲线
（r=0.29, xm=1）
fplot, 函数popinc.m

图2：Logistic模型x~t曲线（r=0.29，xm=1，x0=0.05）

fplot, 函数popcal.m

图3：输出阻滞增长模型拟合图形
（r=0.2557，x0=x(1)，xm=392.0886）

fplot, 函数popcal.m
%Exp01m 文件名：popinc.m
function y = popinc(x)
r=0.29;
xm=1;
y=r*x-r*x.^2/xm;
函数式m文件：popinc.m
阻滞增长模型

%Exp01m 文件名：popcal.m
function x = popcal(t,r,x0,xm)
x=xm/(1+(xm/x0-1)*exp(-r*t));
函数式m文件：popcal.m
阻滞增长模型

输出图1：Logistic模型dx/dt~x曲线（r=0.29, xm=1）

输出图2：Logistic模型x~t曲线（r=0.29，xm=1，x0=0.05）

输出图3：输出阻滞增长模型拟合图形（r=0.2557，x0=x(1)，xm=392.0886）

濒危物种(Florida沙丘鹤)自然演变和人工孵化
%Exp02m 文件名：exam0201.m
clear all
% x0:初始值；b:人工孵化数
x0=100;n=20;b=0;
k=(0:n)';
y1=exf0201(x0,n,0.0194,b);% 给定x0,n,r，b,调用exf0201计算
y2=exf0201(x0,n,-0.0324,b);
y3=exf0201(x0,n,-0.0382b);
round([k,y1',y2',y3']), % 对结果四舍五入取整
plot(k,y1,k,y2,':',k,y3,'--'), % 将3条线画在一个图上
gtext('r=0.0194'), gtext('r=-0.0324'), gtext('r=-0.0382'),
% 在图上做标记
命令式文件：exam0201.m

调用函数exf0201.m
% Exp02m 文件名：exf0201.m
function x=exf11(x0,n,r,b)
% 建立名为exf11的函数M文件， x0,n,r,b可调节
a=1+r;
x=x0; % 赋初值
for k=1:n
x(k+1)=a*x(k)+b; % 按照（3）迭代计算
end
函数式文件：exf0201.m
输出结果：
ans =
0 100 100 100
1 102 97 96
2 104 94 93
3 106 91 89
4 108 88 86
5 110 85 82
6 112 82 79
7 114 79 76
8 117 77 73
9 119 74 70
10 121 72 68
11 124 70 65
12 126 67 63
13 128 65 60
14 131 63 58
15 133 61 56
16 136 59 54
17 139 57 52
18 141 55 50
19 144 53 48
20 147 52 46

一年生植物的繁殖
%Exp02m 文件名：exam0202.m
clear all
n=20;
k=(0:n)';
y1=exf0202(100, n+1,0.18); % 给定x0, n, b，调用exf0202计算
y2=exf0202(100, n+1,0.19);
y3=exf0202(100, n+1,0.20);
round([k,y1',y2',y3']),
plot(k,y1,k,y2,':',k,y3,'--'),
gtext('b=0.18'),gtext('b=0.19'),gtext('b=0.20')
命令式文件：exam0202.m

调用函数exf0202.m
%Exp02m 文件名：exf0202.m
function x=exf12(x0,n,b)
% 建立名为exf12的函数M文件， x0, n, b可调节
c=10;a1=0.5;a2=0.25;
p=-a1*b*c;q=-a2*(1-a1)*c*b^2;
x(1)=x0; % 赋初值
x(2)=-p*x(1); % 按照（10）迭代计算
for k=3:n
x(k)=-p*x(k-1)-q*x(k-2); % 按照（10）迭代计算
end
函数式文件：exf0202.m
输出结果：
ans =
0 100 100 100
1 90 95 100
2 85 95 105
3 80 94 110
4 76 94 115
5 71 93 121
6 67 93 127
7 63 93 133
8 60 92 139
9 56 92 146
10 53 91 152
11 50 91 160
12 47 90 167
13 45 90 175
14 42 90 184
15 40 89 192
16 37 89 202
17 35 88 211
18 33 88 221
19 31 88 232
20 30 87 243

汽车租赁公司的运营
%Exp02m 文件名：exam0204.m
clear all
A=[0.6,0.2,0.1;0.3,0.7,0.3;0.1,0.1,0.6];
x(:,1)=[200,200,200]’; % 赋初值
n=10;
for k=1:n
x(:,k+1)=A*x(:,k); % 按照（31）迭代计算
end
round(x),
k=0:10;
plot(k,x),grid,
gtext(‘x1(k)’),gtext(‘x2(k)’),gtext(‘x3(k)’),
命令式文件：exam0204.m
输出结果：
ans =
200 180 176 176 178 179 179 180 180 180 180
200 260 284 294 297 299 300 300 300 300 300
200 160 140 130 125 122 121 121 120 120 120

按年龄分组的种群增长
%Exp02m 文件名：exam0205m
% 按年龄分组的种群增长
clear all
b=[0,0.2,1.8,0.8,0.2];
s=diag([0.5,0.8,0.8,0.1]); % 对角阵，对角元素为0.5,0.8,0.8,0.1
L=[b;s,zeros(4,1)]; % 按照（36）构造矩阵L
x(:,1)=100*ones(5,1); % 赋初值
K=30;
for k=1:K
x(:,k+1)=L*x(:,k); % 按照（37）迭代计算
end
round (x),
y=diag([1./sum(x)]); % 为向量x归一化做的计算
z=x*y, % z是向量x的归一化
k=0:K;
subplot(1,2,1), plot(k,x),grid % 在一个图形窗内画两张图
subplot(1,2,2), plot(k,z),grid
命令式文件：exam0205.m
输出结果：
120 117 120 127 128 130 136 138 141 146 149 152 157 161 165
89 96 94 96 102 102 104 109 111 113 117 120 122 126 129
8 9 10 9 10 10 10 10 11 11 11 12 12 12 13
z =
Columns 1 through 10
0.2000 0.5769 0.4564 0.3658 0.5003 0.4781 0.4050 0.4712 0.4766 0.4305
02000 0.0962 0.3112 0.2599 0.1462 0.2519 0.2594 0.1799 0.2273 0.2488
0.2000 0.1538 0.0830 0.2836 0.1662 0.1178 0.2187 0.1844 0.1388 0.1898
0.2000 0.1538 0.1328 0.0756 0.1813 0.1339 0.1023 0.1555 0.1423 0.1160
0.2000 0.0192 0.0166 0.0151 0.0060 0.0183 0.0145 0.0091 0.0150 0.0149
Columns 11 through 20
0.4575 0.4708 0.4445 0.4529 0.4651 0.4517 0.4522 0.4610 0.4549 0.4530
0.2021 0.2184 0.2387 0.2143 0.2167 0.2313 0.2202 0.2177 0.2267 0.2225
0.1869 0.1544 0.1771 0.1841 0.1641 0.1724 0.1804 0.1695 0.1713 0.1774
0.1426 0.1428 0.1252 0.1366 0.1410 0.1306 0.1345 0.1389 0.1334 0.1340
0.0109 0.0136 0.0145 0.0121 0.0131 0.0140 0.0127 0.0129 0.0137 00131
Columns 21 through 30
0.4583 0.4561 0.4540 0.4569 0.4564 0.4548 0.4562 0.4564 0.4553 0.4559
0.2192 0.2242 0.2231 0.2204 0.2229 0.2231 0.2213 0.2224 0.2229 0.2218
0.1722 0.1715 0.1755 0.1733 0.1721 0.1743 0.1737 01726 0.1737 0.1737
0.1373 0.1348 0.1342 0.1363 0.1353 0.1346 0.1357 0.1355 0.1349 0.1354
0.0130 0.0134 0.0132 0.0130 0.0133 0.0132 0.0131 0.0132 0.0132 0.0131
Column 31
0.4562
0.2222
0.1730
0.1355
0.0132

离散形式的阻滞增长模型（Logistic模型）
%Exp02m 文件名：exam0206.m
r=[0.3,1.8, 2.5];
x=0.1; % 赋初值
n=40;
for j=1:3
R=r(j); % 取r值
for i=1:n
x(i+1)=x(i)+R*x(i)*(1-x(i)); % 按照（48）迭代计算
end
xx(:,j)=x';
end
k=(0:40)';
[k,xx] % 输出结果
subplot(1,3,1),plotk,xx(:,1)), % 在一个图形窗内画3张图
subplot(1,3,2),plot(k,xx(:,2)),
subplot(1,3,3),plot(k,xx(:,3)),
命令式文件：exam0206.m
输出结果：
ans =
0 0.1000 0.1000 0.1000
1.0000 0.1270 0.2620 0.3250
2.0000 0.1603 0.6100 0.8734
3.0000 0.2006 1.0382 1.1498
4.0000 0.2487 0.9668 0.7192
5.0000 0.3048 1.0246 1.2241
6.0000 0.3684 0.9792 0.5384
7.0000 0.4382 1.0158 1.1597
8.0000 0.5120 0.9869 0.6967
9.0000 0.5870 1.0102 1.2250
10.0000 0.6597 0.9917 0.5360
11.0000 0.7271 1.0065 1.1578
12.0000 0.7866 0.9947 0.7011
13.0000 0.8370 1.0042 1.2250
14.0000 0.8779 0.9966 0.5359
15.0000 0.9101 1.0027 1.1577
16.0000 0.9346 0.9978 0.7012
17.0000 0.9529 1.0017 1.2250
18.0000 0.9664 0.9986 0.5359
19.0000 0.9761 1.0011 1.1577
20.0000 0.9831 0.9991 0.7012
21.0000 0.9881 1.0007 1.2250
22.0000 0.9916 0.9994 0.5359
23.0000 0.9941 1.0005 1.1577
24.0000 0.9959 0.9996 0.7012
25.0000 0.9971 1.0003 1.2250
26.0000 0.9980 0.9998 0.5359
27.0000 0.9986 1.0002 1.1577
28.0000 0.9990 0.9999 0.7012
29.0000 0.9993 1.0001 1.2250
30.0000 0.9995 0.9999 0.5359
31.0000 0.9997 1.0001 1.1577
32.0000 0.9998 0.9999 0.7012
33.0000 0.9998 1.0000 1.2250
34.0000 0.9999 1.0000 0.5359
35.0000 0.9999 1.0000 1.1577
36.0000 0.9999 1.0000 0.7012
37.0000 1.0000 1.0000 1.2250
38.0000 1.0000 1.0000 0.5359
39.0000 1.0000 1.0000 1.1577
40.0000 1.0000 1.0000 0.7012

寄主—寄生模型
%Exp02m 文件名：exam0207.m
% 寄主—寄生模型
clear all
n=100;
% 参数：x0,y0,n, r, N, a, b
z=exf0207(50,10,n,1.5,100,0.025,0.02);
% 给定参数值，调用exf0207计算
k=0:n;
plot(k,z(:,1),k,z(:,2)),grid % 作xk，yk的图形
命令式文件：exam0207.m

调用函数exf0207.m
%Exp02m 文件名：exf0207.m
function z=exf0207(x0,y0,n,r,N,a,b)
% 建立名为exf0207的函数M文件，各参数可调节
x=x0;y=y0; % 赋初值
for k=1:n
x(k+1)=x(k)+r*x(k)*(1-x(k)/N)-a*x(k)*y(k);% 按照（61）迭代计算
y(k+1)= b*x(k)*y(k);
end
z=[x',y'];
函数式文件：exf0207.m
输出结果：

13