《MATLAB基础应用》配套教学课件

单击此处添加副标题内容MATLAB基础应用第一章Matlab中的数组操作matlab中的运算和操作是以数组为对象的，数组又包括：数值数组、字符数组、元胞数组等。数值数组：（1）n元数值向量（行向量与列向量）

（2）数值矩阵

（3）由数值矩阵构成的元胞数组几个标点符号的作用：逗号：用来将数组中的元素分开。（可用空格代替）分号：用来将矩阵中的行分开。（可用回车键代替）冒号：相当于文字中的省略号。中括号：界定数组的首与尾。一、数组的建立1.直接输入法

matlab在创立数组时以逗号或空格表示分列，分号或回车表示分行。数组开头“[”、结尾“]”

行数组：如a=[1,2,3,8,-1]

列数组：b=[1;2;3;8;-1]或a’

矩阵：A=[2,4,1;8:-2:4;2,4,6]2.通过数组编辑器生成矩阵步骤：先建立空矩阵a=[],

然后在工作空间(workspace)中点开a进入数组编辑器，输入元素。3.用函数创建数组定步长生成法：

x=a:t:b(t步长，省略是为1)；定数线性采样法：

x=linspace(a,b,n)，

a与b是数组的第一个和最后一个元素，

n是采样的总点数。x=linspace(2,5,6)x=2.00002.60003.20003.80004.40005.0000

zeros(m):m阶全零方阵

zeros(m,n):m×n阶全零方阵

eye(m):m阶单位阵

ones(m):m阶全1方阵

ones(m,n):m×n阶全1方阵

rand(m):m阶均匀分布随机方阵

randn(m):m阶正态分布随机方阵4.元胞数组的创建

元胞数组是MATLAB的一种特殊数据类型，可以将元胞数组看做一种无所不包的通用矩阵，或者叫做广义矩阵。组成元胞数组的元素可以是任何一种数据类型的量，每一个元素也可以具有不同的尺寸，每一个元素的内容也可以完全不同，元胞数组的元素叫做元胞。建立元胞数组：{}a={'matlab',20;ones(2,3),1:10}a='matlab'[20][2x3double][1x10double]二、数组的操作数组的编址：数组a建立后，a中各元素的编址方法如下：

单下标编址：a(1)表示a的第1个元素，a(n)表示a的第n个元素，对于二元数组按列优先原则进行单下标编址。双下标编址：a(2,3)表示矩阵a的第2行第3列元素。1.数组元素与子数组的提取

提取数组a的第3个元素：y=a(3)

提取a的第3到7个元素：y=a(3:7),a=linspace(1,20,6)a=1.004.808.6012.4016.2020.00

提取a的第1,3,5个元素构成数组b:b=a(1:2:5)

提取a的第2到5个元素，并反转次序构成数组b1:b1=a(5:-1:2)

按条件提取子数组:

提取a的元素值大于10的元素构成数组b2b2=a(find(a>10))b1=16.200012.40008.60004.8000

b2=12.4016.2020.00二维数组A的元素的提取：由于数组A有两种编址方法，matlab会根据接受的指令，先判断是哪一种编址方法，然后再进行元素的提取。如：

A=[1,2,3,4,5;0,1,2,3,4;-1,0,1,2,3;-2,-1,0,1,2]1234501234-10123-2-1012b1=A(5)b2=A(2,3)b3=A(2,:)b4=A([2,3])b5=[A(2,:);A(1,:)]b6=A([2,4],[3,5])

b1=2b2=2b3=[01234]b4=[0,-1]b5=0123412345b6=2402b7=A([2,4],[3,5,1])

b7=24002-2元胞数组元素的提取：（）和

{}有着本质的区别，

{}用于表示元胞的内容，（）小括号表示指定的元胞。a={'matlab',20;ones(2,3),1:10;ones(4,5),eye(4)}a='matlab'[20][2x3double][1x10double][4x5double][4x4double]a='matlab'[20][2x3double][1x10double][4x5double][4x4double]b=[4x4double]b=a(3,2)b=a{3,2}b=1000010000100001a{2,3}=cell(2)改变元胞数组元素的元胞：赋值

a='matlab'[20][][2x3double][1x10double]{2x2cell}[4x5double][4x4double][]2.数组拼接与数组中的元素值的改变在命令窗口实验:a=1:2:11a(1)=0a(1:4)=[2,-1,-2,-3]a([2,5])=[1.5,0.5]x=0:2:10,y=[-2,-5,-8],xx=[x,y]

yy=[xx([2,5]);y(2:3)]a=1357911a=0357911a=2-1-2-3911a=21.5-2-30.511xx=0246810-2-5-8

yy=28-5-8x=0246810y=-2-5-8空数组的使用：建立空数组A：A=[]，空数组大小任意。可用空数组删除已有数组中的元素B=1:8B(1:2:5)=[]B=24678A=[2,3,4,5,6;1,2,3,4,5;0,1,2,3,4;-1,0,1,2,3]删除矩阵A第3行：

A(3,:)=[]删除矩阵A第2列：

A(:,2)=[]3.常用的数组操作命令(1)确定数组大小命令A=[2,5,1,-2,7;0,3,-1,-2,4;0,0,5,8,8]求数组A行数与列数的最大值：n=length(A)提取数组A的行数与列数：[m,n]=size(A)(2)排序命令将一维数组x的元素排序：x=[3,-1,2,5,7,4,6,11,13,9,2,0,7,8]b=sort(x),[b,k]=sort(x)b=-1022345677891113k=2123111647513141089(3)改变数组形状的命令

x=[3,-1,2,5,7,4,6,11,13,9,2,8]

将一维数组x按条件转化为矩阵：

B=reshape(x,3,4)(4)数组的复制：

c=[1,2,5]c1=repmat(c,4,1)c2=repmat(c,1,4)c3=repmat(c,3,2)c4=repmat(c,3)B=3569-1711224138c3=125125125125125125(5)稀疏矩阵与满矩阵的转化：稀疏矩阵生成命令：sparse(a,b,c)

数组a,b,c的大小必须相同

数组a与b分别指定元素的行标与列标，数组c指定元素的值

A=sparse([2,4,18],[3,12,20],[-5,-3,-8])

创建稀疏矩阵A，A的(2,3),(4,12),(18,20)元素分别为-5,-3,-8，其余元素为零，A为18×20阶矩阵。

将稀疏矩阵x变回满矩阵:A=full(A)A=(2,3)-5(4,12)-3(18,20)-8(6)sum(A):矩阵A按列求和，返回一个行向量；

sum(A,2):矩阵A按行求和，返回一个列向量。

max(A):返回由矩阵各列的最大值构成的向量。

max(A,B):返回A与B对应元素最大值构成的矩阵

min(A)，min(A,B)类似创建18×20稀疏矩阵A，使A的(2,3),(4,12),(15,16)元素分别为-5,-3,-8，其余元素为零？A=sparse([2,4,15,18],[3,12,16,20],[-5,-3,-8,0])A=[1,2,3,4,5;0,6,2,3,4;-1,0,7,2,3;2,-1,0,1,2]1234506234-107232-1012B=sum(A,2)B=1515114b1=max(A)b1=26745返回由矩阵A各行的最大值构成的列向量？b2=max(A')'b3=max(A,2)?b3=22345262342272322222b2=5672

b=diag(A):提取方阵A的对角线元素构成列向量b

A=diag(b):用一维数组b的元素生成对角方阵A

(7)diag命令：

A=diag(b,k):b为一维数组，k为整数将b元素作为偏离主对角线的第k条对角生成方阵Ab=[2,3,-1,5,6],A=diag(b,1)B=diag(b,-2)A=020000003000000-100000050000006000000B=000000000000002000000030000000-1000000050000000600例1.1输入n阶矩阵n=input('输入方阵阶数n=')a1=4*ones(n,1);a2=2*ones(n-1,1);a3=ones(n-2,1);A1=diag(a1)+diag(a2,1)+diag(a3,2);A2=diag(a2,-1)+diag(a3,-2);A=A1+A2

(8)find命令：

find(A)找出A的不为0的元素的下标

find(A,k)找出A的前k个不为0的元素的下标

find(A,k,’last’)找出A的后k个不为0的元素的下标

find(g(A))，其中g(A)是数组A的逻辑表达式，返回数组A中满足条件g(A)的元素下标。

A=[0,0,2,-1,3,0,0,5,0,6,-7,0,0,9];b1=find(A,3)b2=find(A,2,'last')

b1=345b2=1114B=[0,1,0;2,3,0;4,0,0]c1=find(B)[m,n]=find(B)B=010230400c1=2345m=2312n=1122B=[0,1,0;2,3,0;4,0,0]

t=find(B>2)[x,y]=find(B>2)[m,n]=find(B>=1&B<=3)B=010230400t=35x=32y=12m=212n=122求连乘积

c1=prod(4:6)c2=cumprod(4:6)

求组合nchoosek(a,k)D=nchoosek(5,2)

D=10

M=323137212717N=

wxy

wxz

wyzxyzb='wxyz'N=nchoosek(b,3)(9)排列组合c1=120c2=420120a=[3,2,1,7];M=nchoosek(a,2)求向量的全排列

perms([2,1,8])三、数组的运算1.数值运算数组的运算，也称点运算，是同阶数组对应分量的运算。包括点乘、点除和点乘方，对应的运算符号为.*,./.^这些运算符的公共特点是在算符前加上一个英文句号，以便与对应的矩阵运算相区别。设A与B为同阶数的数组，k为常数，

A+B，A-B，k*AA.*B,A./B,A.^n

2.关系运算与逻辑运算关系运算：等于：==，不等于：~=，小于：<，大于：>，小于等于：<=，大于等于：>=。逻辑运算：与：&；或：|；非：~。关系运算返回值为0或1

如a=3>5，得a=0例1.2已知数组A=[2,5,1,-2,7;0,3,-1,-2,4;0,0,5,8,8]，要求将A中大于4的元素减去2，小于0的元素加上1，其余元素不变构成矩阵B。A=[2,5,1,-2,7;0,3,-1,-2,4;0,0,5,8,8]A=251-2703-1-2400588B=231-15030-1400366B=A+(-2)*(A>4)+(A<0)例1.3矩阵A=[2,6,1,-12,7;0,3,4,-5,4;1,0,5,8,9]，要求将A中能被3整除的元素保留其余元素不变构成矩阵B。A=[2,6,1,-12,7;0,3,4,-5,4;1,0,5,8,9]A=261-127034-5410589B=060-1200300000009B=A.*(mod(A,3)==0)3.数组的集合运算

a=1:12;b=2:10;

setdiff(a,b)（a与b的差集）得：1，11，12intersect(a,b)（a与b的交集）得：2，3，4，5，6，7，8，9，10union(a,b)（a与b的并集）四、矩阵的运算1.矩阵的基本运算数k与矩阵A的运算：k+Ak*A

加法运算：A+B，A-B，（A与B为同阶矩阵）乘法运算：A*B（A为m×k矩阵，B为k×n矩阵）

右除运算：A/B（A乘B的逆，B为可逆矩阵）左除运算：A\B（

A的逆乘B，A为可逆矩阵）幂运算：若A为方阵，A^n为矩阵A的n次幂2.基本矩阵函数

转置A’特征值eig(A),行列式det(A),秩rank(A),逆inv(A),迹trace(A)，条件数cond(A)

练习题输入一个矩阵A，提取A的第2行第1列的元素；提取A的第1,3,4列所有元素；让A的第1列和第3列互换；删除A的第1列。将矩阵A的第2行元素扩大2倍再增加3后作为A的第3行。求出数组x中所有奇数的和与所有的偶数和

x=[8,3,5,-1,6,10,7,2,9,11,-4,-5,-1,12,11];4.生成5阶随机数矩阵，其中的随机数服从均值为20标准差为2的正态分布

5.生成一个15×20阶的均匀分布随机矩阵A，并将A中大于0.5的元素换成0.56.生成一个20阶的方阵A，A中(2,3),(4,5),(8,6),(5,1),(12,2),(14,7),(15,12),(18,7)位置的均为3，其余元素为零.7.计算和式：

8.查寻fliplr命令的用法，并建立20阶的方阵9.已知矩阵

求（1）A的逆矩阵（2）解方程组Ax=b

（3）求A的特征值与特征向量34第二章Matlab程序设计一、自定义函数1.匿名函数格式：函数句柄

=@(自变量列表)函数表达式

例2.1定义函数，并计算f(x)在点x=-2,1,2.5,3,5.2的值。f=@(x)x.^2+3*x+5x=[-2,1,2.5,3,5.2]y1=f(x)y1=3.009.0018.7523.0047.6435

例2.2定义函数，并计算g(3,2)g=@(x,y)sin(x.^2+3*y)y=g(3,2)y=0.6503

例2.3定义复合函数并计算函数z在(x,y)=(2,1),(3,2),(4,1),(4,3)处的值u=@(x,y)x./yv=@(x,y)3*x-2*yf=@(x,y)u(x,y).^2.*log(v(x,y))x=[2,3,4,4]y=[1,2,1,3]z=f(x,y)z=5.54523.621236.84143.1854362.m-文件函数在matlab界面上先点击新建文件按钮，进入窗口Editor-Untitled，按如下格式建立函数并保存。格式：function[y1,y2]=ff(x1,x2)……y1=…y2=…输入变量：x1,x2，输出变量：y1,y2注：输入、输出变量可以为1个或多个。37例2.4建立函数文件，实现对矩阵x的非零元素取倒数，零元素不变，然后调用该函数。（在窗口Editor-Untitled写入）functiony=fun204(a)y=spfun(@(x)1./x,a);（保存）（在窗口Command-Window写入）a=[0,0,2;3,0,1;4,0,0]b=fun204(a)c=full(b)c=000.50000.333301.00000.25000038二、m-文件的建立（程序设计）编辑m-文件：在窗口Editor-Untitled按要求编辑程序,输完程序后保存,在对话框中输入文件名.运行m-文件:在窗口Editor-Untitled单击Debug→run,或按F5键。或在窗口Command-Window输入文件名并回车。39三、流程控制语句（一）分支结构

1.if条件语句

if表达式执行语句

end

if表达式执行语句1else

执行语句2end

if表达式1

执行语句1

elseif

表达式2

执行语句2else

执行语句3end40例2.5定义函数并计算f(x)在x=-5,-3,-1.5,-1,0.5,1,2处的值functiony=ff(x)ifx<-2y=x.^2-2;elseifx<0y=x.^2.*exp(x);elsey=log(x+10);endx=[-5,-3,-1.5,-1,0.5,1,2]y=ff(x)y=1.60941.94592.14012.19722.35142.39792.484941

2.switch分值选择语句

switch表达式

case常量表达式1

语句块1case常量表达式2

语句块2……case常量表达式n

语句块notherwise

语句块n+1end42function[m,y]=ff(x)s=sum(x)m=mod(s,5)switchmcase1t=mod(x,2);y=sum(t.*x);case2t=1-mod(x,2);y=sum(t.*x);otherwise

y=max(x)-min(x);end例2.6定义函数，对数组x做如下处理：记数组x所有元素和为s，（1）若s除以5余数为1，求出x中所有奇数的和。（2）若s除以5余数为2，求出x中所有偶数的和。（3）其余情况求出x中最大值与最小值的差。x=[2,-4,3,6,7,11][m,y]=ff(x)m=0y=1543（二）循环结构

1.for循环语句

for变量=数组执行语句1……

执行语句nend若数组为行向量，依次取数组元素执行循环体。若数组为矩阵，依次取矩阵的列执行循环体。44

例2.7随机生成15个0~99的整数构成数组x，并构造下列数组：(1)将数组x的第偶数个元素平方，其余元素不变。(2)将数组x的第偶数个元素平方，其余元素置0。

x=floor(rand(1,15)*100);y1=x;y2=zeros(1,15);

fori=2:2:14y1(i)=x(i)^2;y2(i)=x(i)^2;

end

[x;y1;y2]45

例2.8数组x=(8,3,5,-1,6,10,7,2,9,11),求出其中所有除以5余1的元素之和以及所有除以5余4的元素之和.x=[8,3,5,-1,6,10,7,2,9,11];s1=0;s2=0;

form=x

switchmod(m,5)

case1s1=s1+m;

case4s2=s2+m;

endend

s1,s2

不用循环完成

x=[8,3,5,-1,6,10,7,2,9,11];i1=find(mod(x,5)==1);s1=sum(x(i1))i2=find(mod(x,5)==4);s2=sum(x(i2))46例2.9对矩阵A的各列按例2.6的要求处理。-1-3-14-52-3-4-218411-537-2113121413792115-5174877823A=[…];i=1;a=[],b=[];

fork=A[m,y]=ff(k);

a(i)=m;

b(i)=y;i=i+1;

endD=[a;b]D=2140324921142647

2.while-end循环语句

while表达式执行语句1……

执行语句nend注：表达式一般是由逻辑运算和关系运算组成的表达式，表达式的值非零继续循环，表达式值为0终止循环。48例2.10记三位正整数x的个位为a，十位为b，百位为c求出满足ab-2c=4a+b(c+2)的第一个三位数。N=100;m=1;whilema=mod(N,10);b=mod((N-a)/10,10);c=(N-a-10*b)/100;ifa*b-2*c==4*a+b*(c+2)m=0;N1=N;

endN=N+1;endN1运行结果：N1=28949

（三）

程序的流程控制命令：

continue

终止一次循环

break

结束循环

(1)a=3;b=2;

fori=1:4b=b+1;

ifi<2

continueenda=a+2;

end

a,b(2)a=3;b=2;

fori=1:4b=b+1;

ifi>2

breakenda=a+2;

end

a,b50例2.11在某地区开掘水渠，已知该水渠经过的若干点坐标(X,Y)，X=[7650,7950,8250,8550,8850,9150,9450,9750,10050,10350,10650,10950,11250,11550,11850,12150,12450,12750,13050,13350,13650,13950,14250,14550];Y=[4650,4890,5100,5190,5220,5190,5130,4890,4650,4470,4320,4200,4080,3930,3840,3720,3660,3630,3570,3540,3600,3690,3780,3900];假设该水渠是以这些点为节点的折线，（1）作出水渠的折线图，求出水渠的总长度。（2）A(7650,4650)为水渠起点，B(xn,yn)为水渠上任一点，编写matlab函数要求输入xn，返回水渠上A到B段长度及yn。（3）求x1=8100与x2=13500所对应的水渠上的点的坐标，并求这两点间水渠的长度。（4）若将水渠分为等长的6段，求出分点坐标51（1）求水渠的总长度记已知节点的坐标：Ai(xi,yi),(i=1,2,…,n)水渠总长度：X=[7650,7950,8250,8550,8850,9150,9450,9750,10050,10350,10650,10950,11250,11550,11850,12150,12450,12750,13050,13350,13650,13950,14250,14550];Y=[4650,4890,5100,5190,5220,5190,5130,4890,4650,4470,4320,4200,4080,3930,3840,3720,3660,3630,3570,3540,3600,3690,3780,3900];

plot(X,Y,'*-');n=length(X);x1=X(2:n)-X(1:n-1);y1=Y(2:n)-Y(1:n-1);L=sum((x1.^2+y1.^2).^0.5)%水渠长度52L=7.5145e+00353（2）编写函数ABAmAm+1先找到B(xn,yn)所在的线段AmAm+1搜索节点Ai(xi,yi)的横坐标Ai(xi,yi),(i=1,2,…,n)找到第一个xi≥xn的xi，得下标m得到线段AmAm+1的方程，可求出yn水渠上A到B段长度：54function[p,xn,yn]=fun0000(X,Y,xn)n=length(X);X0=X(1:n-1);Y0=Y(1:n-1);X1=X(2:n);Y1=Y(2:n);D=sqrt((X1-X0).^2+(Y1-Y0).^2);%各线段的长度，fori=1:n

if

X(i)>=xnm=i-1;

breakendendk=(Y(m+1)-Y(m))/(X(m+1)-X(m));yn=Y(m)+k*(xn-X(m));s=sum(D(1:m-1));p=s+sqrt((xn-X(m)).^2+(yn-Y(m)).^2);55（3）求指定段的水渠长度

xc=8100;xd=13500;[z1,xc,yc]=fun0000(X,Y,xc)[z2,xd,yd]=fun0000(X,Y,xd)

zd=z2-z156

练习题1.定义函数并计算f(x)在x=1,1.5,2,2.5,3,…,9.5,10处的值2.定义函数在单位圆周上等距取6个点，并计算f(x,y)在这些点的值。3.编程产生裴波那契数列:1,1,2,3,5,8,13,21,…要求数列最后一个数不超过1000。（提示：Fn+2=Fn+1+Fn，F1=F2=1)4.计算级数的近似值，要求误差小于10-5576.某钢管厂的原料钢管长为19m，现要求将原料钢管切割成4m,6m和8m钢管，编程求出有多少种切割模式可使余料小于4m，并给出这些切割模式的方案。5.随机生成5行15列的矩阵A=(aij)，其元素为0~9的整数，要求用双重循环和单重循环两种方法构造矩阵B=(bij)。当i+j为奇数时当i+j为偶数时其中8.完成例2.11第(4)问

7.随机生成20行15列的矩阵A=(aij)，其元素为0~9的整数，求出其中所有除以5余1的元素之和以及所有除以5余4的元素之和.第三章matlab图形功能

绘图原理是描点法；一、二维图形1．plot函数给定图形上的n个数据点的x坐标与y坐标，格式：plot(x,y)或plot(x,y,S)用字符串S设置曲线的颜色、线形和点的形状。将这n个点依次连接起来构成折线。颜

色线

型点

型r红色—实线p五角星k黑色—

—虚线h六边形b蓝色—.点划线xX形g绿色：点线o圆圈y黄色**m洋红.点w白色d菱形字符串S所设定的颜色、线形状和点的形状：控制方法：plot(x,y,‘颜色线型点型’)例3.1作出的图形clear,clcx=0:pi/20:2*pi;y=sin(x).*cos(x).^2;plot(x,y,'r--h')2.subplot函数

格式：subplot(m,n,p)将图形窗口分成mn个坐标窗口，并打开第p个窗口供图形操作。例3.2将图形窗口分成3个坐标窗口，分别绘制曲线subplot(3,1,1)x=0:pi/20:2*pi;y=sin(x).*cos(x).^2;plot(x,y,'k-*')subplot(3,1,2)y=exp(x).*sin(x);plot(x,y,'b--o')subplot(3,1,3)x1=cos(x);y1=1/sqrt(2)*sin(x);plot(x1,y1,‘m:.')axisequal3．多重线（在同一个画面上画多条曲线）例3.3x=0:pi/15:2*pi;y1=sin(x);y2=cos(x);plot(x,y1,‘b:*,x,y2,‘r-.p’)多重线还可以用holdon,holdoff来控制完成x=0:pi/15:2*pi;y1=sin(x);plot(x,y1,‘b:*’)holdony2=cos(x);plot(x,y2,‘r:p’)holdoff图形窗口的控制

Tools—editplot—showpropertyeditor例3.4在同一坐标窗口绘制曲线族(n=1,2,3,4,5,6)x=0:pi/40:2*pi;y=sin(x).*cos(x);plot(x,y,'linewidth',2.5)colors=['rkgmy'];holdon

fork=1:5y=y.*cos(x);plot(x,y,colors(k),'linewidth',2.5)

endtitle('y=sin(x)*(cos(x))^n')%加标题

legend('n=1','n=2','n=3','n=4','n=5','n=6')%加图例4.图形的控制grid命令％在图形上加网格；xlabel(‘x轴’)%在x轴加标记；ylabel(‘y轴’)%在y轴加标记；title(‘正弦、余弦曲线’)%给图形加标题；text(x0,y0,‘字符串’)%在图上(x0,y0)处加上字符串；gtext%用鼠标控制在图形上加字符串。axis([xminxmaxyminymax])%确定x，y轴的范围；axisequal%设定x轴y轴单位相同；axissquare%设定图框成方形；axisoff%清除坐标刻度。4.数学函数的简易作图ezplot函数:适用于绘制参数曲线和隐函数确定的曲线。调用格式为ezplot(@(x)fun,[a,b])%绘制显函数曲线

ezplot(@(x,y)fun,[xmin,xmax,ymin,ymax])%绘制隐函数曲线(3)ezplot(@(x)fun1,@fun2)%绘制参数方程曲线

例3.5作出下列函数的图形ezplot(@(x)cos(x).^3)figure(2)ezplot(@(x,y)(x.^4+y.^4-1),[-1.2,1.2,-1.2,1.2])gridonfigure(3)ezplot(@(t)(t-sin(t)),@(t)(1-cos(t)),[-2*pi,4*pi])gridonfplot函数：fplot(fun,[a,b],S)%数值a,b界定绘图区间，字符串S设定曲线颜色和线形fplot(@sin,[-pi/2,3*pi],'r--')5.极坐标作图作图函数polar(t,r)t为极角，r为极半径。

例3.6下列极坐标曲线图t=0:pi/50:2*pi;r1=2*(1+cos(t));polar(t,r1)r2=cos(2*t);figure(2)polar(t,r2)r3=sin(3*t);figure(3)polar(t,r3)二、三维图形（一）空间曲线作图

空间曲线作图命令：plot3(x,y,z,S)

其中x,y,z是同阶向量（曲线上的点列坐标），

S是字符串，用于设置曲线的颜色、线形和点。

例3.7作曲线图:t=0:pi/50:2*pi;x=sin(t).*cos(t).^2;y=cos(t).*sin(t).^2;z=t.*(2*pi-t);plot3(x,y,z,'r--*’)gridon（二）曲面作图原理：matlab的曲面作图是对参数方程设计的。在数值处理中，通过对参数的离散点得到自变量的网格剖分，形成栅格点。栅格数据点的产生用meshgrid命令格式：meshgrid(x,y)t=-5:0.5:5;s=t;[S,T]=meshgrid(s,t);plot(S,T,'.')求出栅格数据点(S,T)处(X,Y,Z)的值然后执行命令：mesh(X,Y,Z)或surf(X,Y,Z)meshc(X,Y,Z)与surfc(X,Y,Z)分别为带等高线的网格曲面和阴影曲面图。例3.8作曲面z=f(x,y)的图形x=-7.5:0.5:7.5;y=x;[X,Y]=meshgrid(x,y);R=sqrt(X.^2+Y.^2+1);Z=sin(R)./R;mesh(X,Y,Z)surf(X,Y,Z)例3.9画抛物柱面y=2x2的图形clc,clearx=-1:0.1:1;

z=0:0.1:2;[X,Z]=meshgrid(x,z);Y=2*X.^2;surf(X,Y,Z)axisequalfigure(2)surf(X,Y,Z)view(-37.5+180,30)axisequal例3.10画圆柱面x2+y2=9的图形clear,clct=0:pi/50:2*pi;s=0:0.1:4;[T,S]=meshgrid(t,s);x=3*cos(T);y=3*sin(T);z=S;mesh(x,y,z)例3.11(动画设计)comet3命令a=10;b=8;t=[0:0.001:10*pi];f=sqrt(a^2-b^2);%地球与另一焦点的距离th=12.5*pi/180;%未经轨道与x-y平面的倾角E=exp(-0.05*t);%轨道收缩率x=E.*(a*cos(t)-f);y=E.*(b*cos(th)*sin(t));z=E.*(b*sin(th)*cos(t));plot3(x,y,z,'b');%画全程轨线pause(2),holdonsphere(15);%画地球axisoff,axisequalcomet3(x,y,z,0.01);%画运动轨线holdoff（三）等高线图1.contour命令：绘制曲面的等高线图用法：contour(z,v)或contour(x,y,z,v)z是一个矩阵，由z=f(x,y)在一组点处的函数值给出。若v是正整数，给出需要画的等高线的条数。若v是向量，给出需要画等高线的指定高度。2.clabel命令：标明等高线的高度值3.contourf命令：绘制填充等高线图4.contour3命令：绘制立体等高线图例3.12画出曲面z=sin(xy)的图形和等高线图x=-2:0.05:2;y=x;[X,Y]=meshgrid(x,y);Z=sin(X.*Y);mesh(X,Y,Z)figure(2)v=[-0.8,-0.5,-0.3,0.2,0.5,0.8];[C,h]=contour(Z,v);clabel(C,h)figure(3)[C,h]=contourf(Z,v);clabel(C,h)

colormapcool%定义等高线的色图

colorbar%定义等高线的色标

练习题1.画出曲线及其包络线2.画出分段函数的图形

的图形。t的取值范围是3.画出极坐标系下函数的图形4.对在同一坐标系下画出一组5.分别画出旋转抛物面的图形（1）在区域上（2）在区域上椭圆的图形，标明不同颜色，x,y轴的坐标刻度相同。6.画出单叶双曲面的图形7.画出椭圆抛物面的等高线图第四章方程与微分方程一、代数方程的求根多项式求根（1）matlab中的多项式是用行向量表示，如：p=[2,3,-1,4,5]

表示多项式：（2）求多项式在某点的值的函数：polyval(p,x)p：多项式对应的行向量,x：待求值的点如：求在x=1,2,-3,5点的值p=[2,3,-1,4,5]x=[1,2,-3,5]y=polyval(p,x)y=1365651625

（3）polyder(p)

求多项式p的导数。（4）poly(x)

构造以x为根的多项式。（5）m=roots(p)

求多项式p的的根如：x=[-1,2,3]q=poly(x)如：p=[-1,2,3,1]x=roots(p)q=1-416x=3.0796-0.5398+0.1826i-0.5398-0.1826i2.求方程根符号解sovle函数格式：solve(f,x)功能：对变量x解方程f=0

注：使用时先用symsx声明x为符号变量symsxy='2*x^2-3*x-8'm=solve(y,x)m1=vpa(m,7)如：求方程2x2-3x-8=0的根m=3/4-73^(1/2)/473^(1/2)/4+3/4m1=-1.3860012.886001如：解方程组

symsxyy1='2*x^2-3*y-8'y2='x-2*y+1'[m,n]=solve(y1,y2)m1=vpa(m,7)n1=vpa(n,7)m=313^(1/2)/8+3/83/8-313^(1/2)/8n=313^(1/2)/16+11/1611/16-313^(1/2)/16m1=2.586476-1.836476n1=1.793238-0.4182379

3.求方程根的数值解大多数非线性方程的求根问题非常困难，通常需要求方程的数值根。常用求数值解的方法有：二分法、不动点迭代法、牛顿迭代法、牛顿下山法、割线法等。（1）求方程数值解的函数fzero格式：fzero(f,x0)功能：求方程f=0在x0附近的一个根例4.1求方程函数sin(x+1)-0.2x=0的根。

x=-5:0.05:5;y=sin(x+1)-0.2*x;

plot(x,y)gridonpausen=input('输入根的个数n=')

fork=1:nb0=input('输入初值b0=')a=fzero('sin(x+1)-0.2*x',b0)

end例4.2求解非线性方程组：（2）求非线性方程（组）数值解的函数fsolve格式：fsolve(f,x0)（x0为初值）h1=ezplot(@(x,y)5*x.^2+y.^2.*sin(y)-3)set(h1,'color',[1,0,0])holdonh2=ezplot(@(x,y)y.^2+10*cos(x.*y)+4*x+2*y-7);set(h2,'color',[0,0,1])holdoffgridon

1,-10.5，-3

-0.5，-3-1，0-0.8，3-1.2，4

functionf=fun402(x)f1=5*x.^2+y.^2.*sin(y)-3;f2=1*y.^2+10*cos(x.*y)+4*x+2*y-7;f=[f1;f2]x=[1,0.5,-0.5,-1,-0.8,-1.2];y=[-1,-3,-3,0,3,4];m=length(x);r=[];

fork=1:m

r(:,k)=fsolve(@fun402,[x(k);y(k)]);

endrr=0.93600.4483-0.3502-0.7746-0.7121-1.1380-1.2139-3.3232-3.35520.05183.09303.4397

1,-10.5，-3

-0.5，-3-1，0-0.8，3-1.2，4将求得的根标注在图形上，观察是否正确h1=ezplot(@(x,y)5*x.^2+y.^2.*sin(y)-3)set(h1,'color',[1,0,0])holdonh2=ezplot(@(x,y)y.^2+10*cos(x.*y)+4*x+2*y-7);set(h2,'color',[0,0,1])r=[0.9360,0.4483,-0.3502,-0.7746,-0.7121,-1.1380;-1.2139,-3.3232,-3.3552,0.0518,3.0930,3.4397];u=r(1,:);v=r(2,:);plot(u,v,'mh')holdoffgridon二、微分方程求解1.求微分方程解析解

dsolve(‘方程1’,‘方程2’,…,‘方程n’,‘初始条件,‘自变量’）

初始条件省缺时，是求微分方程的通解。表示微分方程时，Dy代表y的导数，D2y代表y的二阶导数，D3y代表y的三阶导数……例4.3求解微分方程：

y=dsolve('D2y-2*Dy+y-x^2=0','x')y=4*x+C48*exp(x)+x^2+C49*x*exp(x)+6y=dsolve('D2y+4*Dy+29*y','y(0)=0,Dy(0)=15','x')

结果：y=(3*sin(5*x))/exp(2*x)y=dsolve('D2y-1000*(1-y^2)*Dy-y=0','y(0)=2,Dy(0)=0','x')y=[emptysym]例4.4求微分方程组的通解[x,y,z]=dsolve('Dx=2*x-3*y+3*z','Dy=4*x-5*y+3*z','Dz=4*x-4*y+2*z')2.求微分方程数值解

多数微分方程是没有解析解的，即使一些看上去形式非常简单的微分方程，如这个微分方程不能用初等函数及其积分来表达它的解!

从实际问题中抽象出来的微分方程，通常主要依靠数值解法来解决！微分方程有解的一个充分条件：微分方程初值问题：（1）若f(x,y)在区域上连续，且满足：则初值问题（1）在[a,b]上有唯一解y=y(x)简单的说微分方程数值解法，就是将连续问题离散化，用数值方法给出它的解在某些离散节点上的近似值。求数值解的方法有欧拉法、龙格-库塔法等。

基于龙格-库塔法，matlab求微分方程数值解函数：[x,y]=ode23(‘函数’,‘求解区间’,‘初始值’)其它求解命令：ode45,ode15s,ode23s等

（采用2,3阶龙格-库塔法）

x:自变量值，y:函数值注意：matlab求微分方程数值解时，高阶微分方程必须等价的变换成一阶微分方程组。例4.5求解微分方程初值问题：f=@(x,y)1./(1+x.^2)-2*y.^2ode23(f,[0,6],0)figure(2)ode45(f,[0,6],0)例4.6求解微分方程初值问题：解:设y1=y,问题化为：建立函数文件：function

dy=li406(x,y)dy=zeros(2,1);dy(1)=y(2);dy(2)=1000*(1-y(1)^2)*y(2)-y(1);

[x,y]=ode15s('li406',[0,3000],[2,0])

plot(x,y(:,1))求微分方程(组)数的值解时，需先定义表示方程(组)的函数文件，如：dy=ffff(x,y)然后调用：[x,y]=ode23(‘ffff’,[a,b],x0)

在求解n个方程的微分方程组时，y与x0均为n维向量，函数文件ffff中的待解方程组应以y的分量形式写出。

例4.7导弹追踪问题设位于坐标原点的甲舰向位于x轴上点(1,0)处的乙舰发射导弹，导弹头始终对准乙舰，若乙舰以最大速度v0（常数）沿平行于y轴的直线行驶，导弹速度是5v0

，求导弹运行的曲线；又乙舰行驶多远时会被击中？解：设在t时刻导弹位于P(x,y)处，乙舰位于Q(1,v0t)处，导弹运动曲线：y=f(x)由导弹头始终对准乙舰：由导弹速度是5v0：弧OP长=5|AQ|解:设y1=y,问题化为：

定义函数文件

function

dy=fun407(x,y)

dy=zeros(2,1);dy(1)=y(2);dy(2)=1/5*sqrt(1+y(1)^2)/(1-x);[x,y]=ode15s('fun407',[0,0.99999],[0,0]);

yy=y(:,1);

plot(x,yy,'*')s=yy(end)得：s=0.2017例4.8放射性废料的处理有一个时期美国原子能机构处理核废料时，一直采用将它们装入密封的圆桶里，然后扔到水深约91米海底的方法。对此科学家们表示担心，怕圆桶下沉到海底时，与海底碰撞发生破裂而造成核污染。为此工程师们进行碰撞试验，发现当圆桶下沉速度超过12.2米/秒与海底碰撞时，圆桶就可能发生破裂。按此试验结论，分析这种处理核废料的方式是否安全。已知圆桶重量239.456千克，体积为0.208立方米，海水浮力1025.94千克/立方米。（大量试验表明圆桶下沉时的阻力与下沉速度成正比，比例系数为0.12）分析：求出当圆桶下沉到海底时，下沉速度是否大于12.2米/秒重量W=239.456千克，体积V=0.208立方米，海水浮力B=1025.94V=213.396千克，阻力系数k=0.12，y=dsolve('m*D2y=W-B-k*Dy','y(0)=0','Dy(0)=0')y=(B*m-W*m-B*k*t+W*k*t)/k^2-(B*m-W*m)/(k^2*exp((k*t)/m))得：求解方程：B=1025.94*0.208;k=0.12;W=239.456;m=W/9.8;y=@(t)(B*m-W*m-B*k*t+W*k*t)/k^2-(B*m-W*m)./(k^2*exp((k*t)/m))-91;t=0:0.1:15;

yy=y(t);plot(t,yy,t,0,'r.-')t0=input('输入初值t0=')

tt=fzero(y,t0)方程的解在t=13附近，键盘输入t0=13得：tt=13.2042求圆桶沉到海底时的速度：B=1025.94*0.208;k=0.12;W=239.456;m=W/9.8;tt=13.2042;symstdyy=diff((B*m-W*m-B*k*t+W*k*t)/k^2-(B*m-W*m)./(k^2*exp((k*t)/m)))dd=subs(dyy,tt)dyy=916921625133564767/4222124650659840-1021030434653936875/(4701511720370176*exp((21*t)/4276))dd=13.6361大于12.2米/秒，故这种处理核废料的方式是不安全的。练习题：1.求螺旋线与平面的交点3.求解微分方程初值问题：2.求下列微分方程的通解4.设y1与y2分别为捕食者与被捕食者两种动物种群的数量，那么Lotka-Voltera捕食者与被捕食者模型为：对求解该模型，作出这两种动物种群的数量的变化曲线。5.有一地下油库，它的内部设计是顶部和底部的半径（1）求出油库的体积（2）求出油库中油量分别为10m3,20m3,30m3,40m3,50m3,60m3,…时油的高度。为2米的圆盘，侧面为抛物线绕z轴旋转所成的旋转曲面（如图），（3）若油的高度h与油量V之间的函数为h=h(V)，画出该函数的曲线图。

第五章拟合与插值

t(h)0.250.511.523468u(g/ml)19.2118.1515.3614.1012.899.327.455.243.01例5.1试验得快速静脉注射下的一组血药浓度数据求血药浓度随时间的变化规律u(t).(t=0时注射300mg，药物排除速率与血药浓度成正比)解：记血液容积v,t=0注射剂量d,血药浓度即为d/v.如何确定常数k和v？一、用最小二乘法确定函数中的参数问题提出：对函数(曲线)y=f(x,k)