本题可用
x=fminbnd（'-0.1*x*（1-x）',0,1）
y=0.1*x*（1-x）
3）指示最大值坐标
用线性绘图函数plot，调用格式如下：
plot（x1,y1,'颜色线型数据点图标',x2,y2,'颜色线型数据点图标',…）
本题可用
holdon;%在上面的同一张图上画线（同坐标系）
plot（[0,x],[y,y],':
',[x,x],[0,y],':
'）;
4）图形的标注
使用文本标注函数text，调用格式如下：
格式1
text（x,y,文本标识内容,'HorizontalAlignment','字符串1'）
x,y给定标注文本在图中添加的位置。
'HorizontalAlignment'为水平控制属性，控制文本标识起点位于点（x,y）同一水平线上。
'字符串1'为水平控制属性值，取三个值之一：
'left'，点（x,y）位于文本标识的左边。
'center'，点（x,y）位于文本标识的中心点。
'right'，点（x,y）位于文本标识的右边。
格式2
text（x,y,文本标识内容,'VerticalAlignment','字符串2'）
x,y给定标注文本在图中添加的位置。
'VerticalAlignment'为垂直控制属性，控制文本标识起点位于点（x,y）同一垂直线上。
'字符串1'为垂直控制属性值，取四个值之一：
'middle'，'top'，'cap'，'baseline'，'bottom'。
（对应位置可在命令窗口应用确定）
本题可用
text（0,y,'（di/dt）m','VerticalAlignment','bottom'）;
text（x,-0.001,num2str（x）,'HorizontalAlignment','center'）;
5）坐标轴标注
调用函数xlabel，ylabel和title
本题可用
title（'SI模型di/dt~i曲线'）;
xlabel（'i'）;
ylabel（'di/dt'）;
☆程序运行结果（比较[138]图2）：
（在图形窗口菜单选择Edit/CopyFigure，复制图形）
1.2画i~t曲线图p136~138
求出微分方程的解析解i（t），画出i~t曲线（i（0）=0.15,k=0.2,t=0~30）（见[138]图1比较）。
参考程序：
%5.1传染病模型——模型2
%文件名：
p136fig1.m
%di/dt=ki（1-i）,i（0）=i0
clear;clc;
x=dsolve（'Dx=k*x*（1-x）','x（0）=x0'）%求微分方程的解析解，为符号表达式
x0=0.15;k=0.2;%xi对应i，xi0对应i0，k对应λ
tt=0:
0.1:
30;%时间单位为天
fors=1:
length（tt）%x的表达式中没有点运算，按标量运算取值xx
t=tt（s）;
xx（s）=eval（x）;%给出xi0=0.2,k=0.2,t,求符号表达式xi的对应值
end%xx为复数表示
plot（tt,xx）;
axis（[03101.1]）;
title（'图1SI模型的i~t曲线'）;
xlabel（'t（天）'）;ylabel（'i（病人所占比例）'）;
提示：
1）求解微分方程dsolve]，见提示；
2）画出i~t曲线（i（0）=0.15,λ=0.2,t=0~30）
用for循环，函数length,eval,plot,axis,title,xlabel,ylabel。
☆程序运行结果（见[138]图1）：
命令窗口中的结果：
图形窗口中的结果（比较[138]图1）：
2.（编程）传染病模型3（SIS模型）
已知传染病模型3（SIS模型）：
其中，
i（t）是第t天病人在总人数中所占的比例。
λ是每个病人每天有效接触的平均人数（日接触率）。
i0是初始时刻（t=0）病人的比例。
σ是整个传染期内每个病人有效接触的平均人数（接触数）。
2.1画
曲线图p138~139
取λ=0.1，σ=1.5，画出如下所示的
曲线图。
试编写一个m文件来实现。
（在图形窗口菜单选择Edit/CopyFigure，复制图形）
（注：
p139图3）
提示：
用fplot函数画出
的曲线图；
在上图上用plot函数画一条过原点的水平线；
用title,xlabel,ylabel标注。
★编写的M文件和运行结果（见[139]图3）：
%传染病模型3（SIS模型）的di/dt~i曲线图
%文件名：
p138fig3.m
%λ=0.1,σ=1.5
clear;clc;
fplot（'-0.1*x*（x-（1-1/1.5））',[00.4-0.00050.003]）;
holdon;
plot（[0,0.4],[0,0],':
'）;
title（'SIS模型的di/dt~i曲线'）;
xlabel（'i'）;
ylabel（'di/dt'）;
2.2画i~t曲线图p138~139
要求：
求出微分方程的解析解i（t）。
取λ=0.2,σ=3,t=0~40，画出如下所示的图形。
试编写一个m文件来实现。
（注：
p139图4）
其中
蓝色实线为i（0）=0.2时的i~t曲线（第1条）；
黑色虚点线为过点（0,1-1/σ）的水平线（第2条）；
红色虚线为i（0）=0.9时的i~t曲线（第3条）。
提示
图例标注可用
legend（'i（0）=0.2','1-1/¦σ','i（0）=0.9'）;
★编写的M文件和运行结果（比较[139]图4）：
解法一：
程序：
%传染病模型3（SIS模型）的i~t曲线图
%文件名：
p138fig4.m
clear;clc;
%λ=0.2,σ=3，i（0）=0.2,x代表i
x=dsolve（'Dx=-0.2*x*（x-（1-1/3））','x（0）=0.2'）%求微分方程的解析解，为符号表达式
tt=0:
0.1:
40;%时间单位为天
fori=1:
length（tt）
t=tt（i）;
xx（i）=eval（x）;
end
plot（tt,xx）;
holdon;
plot（[0,41],[1-1/3,1-1/3],'-.k'）;
%λ=0.2,σ=3，i（0）=0.9
x=dsolve（'Dx=-0.2*x*（x-（1-1/3））','x（0）=0.9'）
tt=0:
0.1:
40;%时间单位为天
fori=1:
length（tt）
t=tt（i）;
xx（i）=eval（x）;
end
plot（tt,xx,':
r'）;
legend（'i（0）=0.2','1-1/σ','i（0）=0.9'）;
axis（[04001]）;
title（'图1SI模型的i~t曲线（λ=0.2,σ=3）'）;
xlabel（'t（天）'）;ylabel（'i（病人所占比例）'）;
命令窗口的结果：
图形窗口的结果：
解法二：
程序：
%传染病模型3（SIS模型）的i~t曲线图
%文件名：
p138fig4.m
clear;clc;
%λ=0.2,σ=3，x代表i
x=dsolve（'Dx=-0.2*x*（x-（1-1/3））','x（0）=x0'）%求微分方程的解析解，为符号表达式
tt=0:
0.1:
40;%时间单位为天
forx0=[0.2,0.9]%i（0）=0.2,0.9
fort=tt
xx（2-（x0==0.2）,round（t/0.1）+1）=eval（x）;
end
end
plot（tt,xx（1,:
）,'-b',[0,41],[1-1/3,1-1/3],'-.k',tt,xx（2,:
）,':
r'）;
legend（'i（0）=0.2','1-1/σ','i（0）=0.9'）;
axis（[04001]）;
title（'图1SI模型的i~t曲线（λ=0.2,σ=3）'）;
xlabel（'t（天）'）;ylabel（'i（病人所占比例）'）;
命令窗口的结果：
图形窗口的结果：
与解法一相同
解法三：
程序
%传染病模型3（SIS模型）的i~t曲线图
%文件名：
p138fig4.m
clear;clc;
x=dsolve（'Dx=-lam*x*（x-（1-1/si））','x（0）=x0'）%求微分方程的解析解，为符号表达式
tt=0:
0.1:
40;%时间单位为天
lam=0.2;si=3;%λ=0.2,σ=3，x代表i
forx0=[0.2,0.9]%i（0）=0.2,0.9
fort=tt
xx（2-（x0==0.2）,round（t/0.1）+1）=eval（x）;
end
end
plot（tt,xx（1,:
）,'-b',[0,41],[1-1/3,1-1/3],'-.k',tt,xx（2,:
）,':
r'）;
legend（'i（0）=0.2','1-1/σ','i（0）=0.9'）;
axis（[04001]）;
title（'图1SI模型的i~t曲线（λ=0.2,σ=3）'）;
xlabel（'t（天）'）;ylabel（'i（病人所占比例）'）;
命令窗口的结果：
图形窗口的结果：
与解法一相同
3.（验证）传染病模型4（SIR模型）p140~141
SIR模型的方程：
设λ=1，μ=0.3，i（0）=0.02，s（0）=0.98。
输入p140的程序并运行，结果与教材p141的图7和图8比较。
ode45,pause的用法见提示。
☆2个M文件（见[140]）和运行结果（比较[141]图7、图8）：
函数M文件：
%5.1传染病模型——模型4（SIR模型）
%文件名：
ill.m
functiony=ill（t,x）
a=1;b=0.3;%λ用a表示，μ用b表示
y=[a*x
（1）*x
（2）-b*x
（1）,-a*x
（1）*x
（2）]';%i用x
（1）表示，s用x
（2）表示
命令M文件：
%5.1传染病模型——模型4（SIR模型）
%文件名：
p140.m
clear;clc;
ts=0:
50;
x0=[0.02,0.98];
[t,x]=ode45（'ill',ts,x0）;[t,x]
plot（t,x（:
1）,t,x（:
2））,grid,pause%图7i（t）,s（t）图形
plot（x（:
2）,x（:
1））;grid,%图8i~s图形（相轨线）
i（t）,s（t）图形（比较[141]图7）：
i~s图形（相轨线）（比较[141]图8）：
4.（验证）人口指数增长模型参数估计及结果分析（美国1790-2000年人口）p163~164
美国1790-2000年人口统计数据（以百万为单位）
年
1790
1800
1810
1820
1830
1840
1850
1860
1870
1880
1890
人口
3.9
5.3
7.2
9.6
12.9
17.1
23.2
31.4
38.6
50.2
62.9
年
1900
1910
1920
1930
1940
1950
1960
1970
1980
1990
2000
人口
76.0
92.0
106.5
123.2
131.7
150.7
179.3
204.0
226.5
251.4
281.4
人口指数增长模型：
x（t）=x0ert
（1）用表中数据进行数据拟合求参数r，x0。
将x（t）=x0ert两边取对数，可得
y=rt+a
其中，y=lnx,a=lnx0，即x0=exp（a）。
采用线性最小二乘法进行数据拟合，用MATLAB中的函数polyfit计算。
（r为每10年估计的人口增长率，x0为1790年估计的初始人口数）
以下是M文件：
clc;formatcompact;
x=[3.95.37.29.612.917.123.231.4...
38.650.262.976.092.0106.5123.2131.7...
150.7179.3204.0226.5251.4281.4];
t=0:
length（x）-1;%t=0为1790年，t=1为1800年，…
y=log（x）;%取1790-2000年的数据
ra=polyfit（t,y,1）;
disp（'用1790-2000年数据估计的参数为：
'）
r=ra
（1）
x0=exp（ra
（2））
☆
（1）运行程序并给出结果（见[167]）：
（2）人口指数增长模型计算结果与实际数据比较（数据表）
以下是M文件：
clc;formatcompact
x=[3.95.37.29.612.917.123.231.4...
38.650.262.976.092.0106.5123.2131.7...
150.7179.3204.0226.5251.4281.4];
t=0:
length（x）-1;
y=log（x）;%取1790-2000年的数据
ra=polyfit（t,y,1）;
r=ra
（1）;x0=exp（ra
（2））;
x2=x0*exp（r*t）;
disp（'指数增长模型拟合美国人口数据的结果（x2）'）
formatshortg;
disp（[1790+10*t;x;round（10*x2）/10]'）;
☆
（2）运行程序并给出结果（见[167]表4中x2列）：
（3）人口指数增长模型计算结果与实际数据比较（拟合图形）
以下是M文件：
x=[3.95.37.29.612.917.123.231.4...
38.650.262.976.092.0106.5123.2131.7...
150.7179.3204.0226.5251.4281.4];
t=0:
length（x）-1;
y=log（x）;%取1790-2000年的数据
ra=polyfit（t,y,1）;
r=ra
（1）;x0=exp（ra
（2））;
t2=linspace（0,length（x）-1,30）;
x2=x0*exp（r*t2）;
plot（t,x,'r+',t2,x2,'b'）;
title（'指数增长模型拟合图形'）;
☆（3）运行程序并给出结果（见[168]图3（b））：
5.（验证，编程）估计阻滞增长模型的参数和绘制图形p165~168
美国1790-2000年人口统计数据（以百万为单位）
年
1790
1800
1810
1820
1830
1840
1850
1860
1870
1880
1890
人口
3.9
5.3
7.2
9.6
12.9
17.1
23.2
31.4
38.6
50.2
62.9
年
1900
1910
1920
1930
1940
1950
1960
1970
1980
1990
2000
人口
76.0
92.0
106.5
123.2
131.7
150.7
179.3
204.0
226.5
251.4
281.4
人口阻滞增长模型：
或
（1）（验证）用1860-1990年的数据拟合估计参数r,xm
用下面方程
估计参数r,xm。
程序如下：
%用数值微分的三点公式计算美国人口增长率（/10年）
clear;clc;
x=[3.95.37.29.612.917.123.231.4...
38.650.262.976.092.0106.5123.2131.7...
150.7179.3204.0226.5251.4];%1790-1990年的人口
xx=x（1,8:
end）;%取1860-1990年的数据
len=length（xx）;h=1;
dx=[-3*xx
（1）+4*xx
（2）-xx（3）,xx（3:
len）-xx（1:
len-2）,...%求数值微分dx
3*xx（len）-4*xx（len-1）+xx（len-2）]/（2*h）;
y=dx./xx;
sr=polyfit（xx,y,1）;
r=sr
（2）
xm=-r/sr
（1）
☆
（1）运行程序并给出结果（见[168]差异大，书上结果数据处理过）：
（2）（编程）计算阻滞增长模型拟合美国人口数据（1790-1990）
取参数r=0.2557,xm=392.0886,x0=3.9，编写程序计算1790-1990年的人口，参考的输出如下（第1列为年，第2列为实际人口，第3列为计算人口）：
★
（2）给出程序及其运行结果（见[167]表4中x列）：
clear;clc;formatcompact;
x=[3.95.37.29.612.917.123.231.4...
38.650.262.976.092.0106.5123.2131.7...
150.7179.3204.0226.5251.4];%1790-1990年的人口
t=0:
length（x）-1;
x0=3.9;r=0.2557;xm=392.0886;%代入参数值
xx=xm./（1+（xm/x0-1）*exp（-r*t））;%计算人口
formatshortg;
[1790+10*t;x;round（10*xx）/10]'
（3）（编程）绘制阻滞增长模型拟合图形（1790-1990）
取参数r=0.2557,xm=392.0886,x0=3.9，编写程序绘制1790-1990年的人口数据拟合图形。
参考图形如下，实线为模型曲线，+为实际数据值。
★（3）给出程序及其运行结果（见[168]图4）：
clear;clc;formatcompact;
x=[3.95.37.29.612.917.123.231.4...
38.650.262.976.092.0106.5123.2131.7...
150.7179.3204.0226.5251.4];%1790-1990年的人口
t=0:
length（x）-1;
x0=3.9;r=0.2557;xm=392.0886;
t3=linspace（0,length（x）-1,30）;
x3=xm./（1+（xm/x0-1）*exp（-r*t3））;
plot（t,x,'r+',t3,x3,'b-'）;
附1：
实验提示
第1.2题
提示：
求解微分方程dsolve
1）求解微分方程
常微分方程符号解用函数dsolve,调用格式如下：
dsolve（'equ1','equ2',…,'变量名'）
以代表微分方程及初始条件的符号方程为输入参数，多个方程或初始条件可在一个输入变量内联立输入，且以逗号分隔。
默认的独立变量为t，也可把t变为其他的符号变量。
字符D代表对独立变量的微分，通常指d/dt。
本题可用
x=dsolve（'Dx=k*x*（1-x）','x（0）=x0'）
2）画出i~t曲线（i（0）=0.15,λ=0.2,t=0~30）
用for循环，函数length,eval,plot,axis,title,xlabel,ylabel。
第3题
ode45,pause的用法
1）求解微分方程的数值解函数ode45，格式如下：
[t,x]=ode45（'fun',ts,x0）
fun是由一个或多个待解方程写成的函数式m文件；
ts=[t0,tf]表示此微分方程的积分限是从t0到tf，也可以是一些离散的点，形式为ts=[t0,t1,…,tf]；
x0为初值条件。
2）等待用户反应命令pause：
程序执行到该命令时暂停，直到用户按任意键后继续（处在命令窗口有效）。
附2：
第5章微分方程模型
[136]5.1传染病模型
[138]题1.1、1.2答案
[139]题2.1
升级会员 