计算与仿真技术详解

计算与仿真技术详解演示文稿目前一页\总数七十七页\编于十一点优选计算与仿真技术目前二页\总数七十七页\编于十一点什么是MATLAB?MATLAB代表MATrixLABoratory它的首创者是美国新墨西哥大学计算机系的系主任CleveMoler博士，他在教授线性代数课程发现其他语言很不方便，篇构思开发了MATLAB。最初采用FORTRAN语言编写，20世纪80年代后出现了MATLAB的第二版，全部采用C语言编写.1984年Moler博士和一批数学家及软件专家创建了MathWorks公司，专门开发MATLAB。目前三页\总数七十七页\编于十一点MATLAB的发展历史1984年，MATLAB第1版(DOS版)

1992年，MATLAB4.0版

1994年，MATLAB4.2版

1997年，MATLAB5.0版

1999年，MATLAB5.3版

2000年，MATLAB6.0版

2001年，MATLAB6.1版

2002年，MATLAB6.5版

2004年，MATLAB7.0版2008年，MATLAB7.6版（支持4线程计算）2009年，MATLAB7.8版（支持8线程计算）是一种解释性语言，一种演草纸式的科学计算语言目前四页\总数七十七页\编于十一点MATLAB能干什么？MATLAB可以进行：数学计算、算法开发、数据采集建模、仿真、原型数据分析、开发和可视化科学和工程图形应用程序的开发，包括图形用户界面的创建。MATLAB广泛应用于：数值计算、图形处理、符号运算、数学建模、系统辨识、小波分析、实时控制、动态仿真等领域。目前五页\总数七十七页\编于十一点掌握MATLAB……MATLAB的构成：MATLAB开发环境：进行应用研究开发的交互式平台MATLAB数学与运算函数库：用于科学计算的函数MATLAB语言：进行应用开发的编程工具图形化开发：二维、三维图形开发的工具应用程序接口(API)：用于与其他预言混编面向专门领域的工具箱：小波工具箱、神经网络工具箱、信号处理工具箱、图像处理工具箱、模糊逻辑工具箱、优化工具箱、鲁棒控制工具箱等几十个不同应用的工具箱。目前六页\总数七十七页\编于十一点MATLAB界面：命令窗口、图形窗口、编辑窗口、帮助窗口。目前七页\总数七十七页\编于十一点常用MATLAB系统命令help在线帮助who显示当前变量whos显示当前变量的详细信息clear清空工作间的变量和函数pack整理工作间的内存load把文件调入变量到工作间save把变量存入文件中echo命令回显what显示指定的matlab文件dir显示目录内容目前八页\总数七十七页\编于十一点基本数学运算在MATLAB下进行基本数学运算，只需将运算式直接打在提示号>>后面，并按Enter键。MATLAB将计算的结果以ans显示。【例】求的算术运算结果。（1）用键盘在MATLAB指令窗中输入以下内容>>(12+2*(7-4))/3^2

（2）在上述表达式输入完成后，按【Enter】键，该就指令被执行。（3）在指令执行后，MATLAB指令窗中将显示以下结果。ans=

2

我们也可给运算式的结果设定一个变量x：x=(5*2+1.3-0.8)*10^2/25x=42变量x的值可以在下个语句中调用：y=2*x+1y=85

目前九页\总数七十七页\编于十一点例：计算在命令窗口中输入sqrt（18）+（7*sin（pi/6））/（5+cos（pi/6））注意：（1）“》”为命令提示符；（2）字母要区分大小写，所有命令和函数名都必须用小写。目前十页\总数七十七页\编于十一点变量命名规则1.变量名的大小写敏感。2.变量的第一个字符必须为英文字母，而且不能超过31个字符。3.变量名可以包含下连字符、数字，但不能为空格符、标点。注：数值型变量和符号型变量在MATLAB中是不相同的，它们之间不能直接进行转化。MATLAB提供了一个将数值型转化成符号型的命令，即sym（syms）。目前十一页\总数七十七页\编于十一点MATLAB提供基本的算术运算有：加(+)、减(-)、乘(*)、除(/)、幂次方(^)，范例为：5+3,5-3,5*3,5/3,5^3MATLAB书写表达式的规则与“手写算式”差不多相同，但要求所有表达式都是以纯文本形式输入。如果一个指令过长可以在结尾加上...（代表此行指令与下一行连续），例如：>>

1*2+3*4+5*6+7*8+9*10+11*12+...13*14+15*16ans=744目前十二页\总数七十七页\编于十一点MATLAB常用数学函数三角函数和双曲函数名称含义名称含义名称含义sin正弦csc余割atanh反双曲正切cos余弦asec反正割acoth反双曲余切tan正切acsc反余割sech双曲正割cot余切sinh双曲正弦csch双曲余割asin反正弦cosh双曲余弦asech反双曲正割acos反余弦tanh双曲正切acsch反双曲余割atan反正切coth双曲余切atan2四象限反正切acot反余切asinh反双曲正弦sec正割acosh反双曲余弦目前十三页\总数七十七页\编于十一点MATLAB常用数学函数指数函数名称含义名称含义名称含义expE为底的指数log1010为底的对数pow22的幂log自然对数log22为底的对数sqrt平方根复数函数名称含义名称含义名称含义abs绝对值conj复数共轭real复数实部angle相角imag复数虚部目前十四页\总数七十七页\编于十一点其他函数名称含义名称含义min最小值max最大值mean平均值median中位数std标准差diff相邻元素的差sort排序length个数norm欧氏（Euclidean）长度sum总和prod总乘积dot内积cumsum累计元素总和cumprod累计元素总乘积cross外积目前十五页\总数七十七页\编于十一点矩阵的MATLAB表示直接输入

>>A=[12,3;456;7,89]冒号操作符>>a=0:1:10>>a=linspace(0,1,10)>>a=logspace(1,2,10)

目前十六页\总数七十七页\编于十一点MATLAB下矩阵的运算矩阵的代数运算转置B=A’加减乘A+BA-BA*B左除A\B即AX=B的解X=A-1B右除A/B即XB=A的解X=AB-1翻转fliplrflipudrot90乘方A^B点运算A.*BA./BA.\BA.^BA.’函数目前十七页\总数七十七页\编于十一点矩阵的逻辑运算运算符&(与)|(或)~(非)结果是一个0-1矩阵。当逻辑表达式的值为真时，赋值1，否则为0

逻辑函数all：当某列的元素都为真时，返回值为1，否则返回0。最终运算结果为一个0-1行向量any：当向量中至少有一个元素为真时，返回值为1，否则返回0。最终运算结果为一个0-1行向量find：用于查找向量中的真元素的下标，返回由所有真元素下标构成的列向量。目前十八页\总数七十七页\编于十一点矩阵的比较关系在MATLAB里共有六个关系运算符<小于 >大于<=小于等于>=大于等于==等于 ~=不等于关系运算符将生成一个0-1矩阵，当运算数相应元素为真时，对应位置上生成1，否则为0。

目前十九页\总数七十七页\编于十一点矩阵元素的数据变换floor(A)向－∞取整ceil(A)向＋∞取整round(A)四舍五入取整fix(A)按离0近的方向取整目前二十页\总数七十七页\编于十一点下标操作矩阵下标（按列优先排列）用两个下标来表示。第i行j列的元素用A(i,j)来表示用一个下标来表示。对于矩阵按列操作，可用单下标引用A(i)。在下标的表达式中使用冒号表示矩阵的一部分。下标引用从序号1开始。目前二十一页\总数七十七页\编于十一点下标操作A=123456789>>A(4)+A(8)目前二十二页\总数七十七页\编于十一点行列删除>>B=[123;456;789];>>B(2,:)=[]>>B=123789目前二十三页\总数七十七页\编于十一点矩阵的尺寸为了判断矩阵的尺寸，Matlab提供了以下函数：[m,n]=size(A) %返回矩阵的行数m和列数nM=size(A) %在矩阵M中返回矩阵的行数、列数

len=length(A) %返回矩阵A行数和列数中的最大值

目前二十四页\总数七十七页\编于十一点特殊矩阵单位矩阵

eye(n)

%生成n维单位阵

eye(m,n)

%生成m×n的单位阵

eye([m,n])%生成m×n的单位阵

eye(size(A))%生成与A矩阵同样大小的单位矩阵oneszeros

目前二十五页\总数七十七页\编于十一点特殊矩阵全1矩阵

ones(n)ones(m,n)ones([m,n]) ones(size(A))全零矩阵

zeros(n)zeros(m,n)zeros([m,n]) zeros(size(A))随机矩阵

rand(n)rand(m,n)rand([m,n])rand(size(A))生成的随机矩阵的各个元素值在0和1之间。对角阵diag(n)对角阵魔术矩阵

magic(n)魔术矩阵的所有行元素和、列元素和、对角线元素和均相等目前二十六页\总数七十七页\编于十一点阵列与矩阵MATLAB的运算事实上是以阵列(array)及矩阵(matrix)方式在做运算.阵列强调元素对元素的运算，而矩阵则采用线性代数的运算方式.宣告一变量为阵列或是矩阵时，须用中括号[]将元素置于其中。阵列为一维元素所构成，而矩阵为多维元素所组成.例如:»x=[12345678];%一维1x8阵列【例】简单矩阵的输入步骤。（1）在键盘上输入下列内容:(以;区隔各列的元素)

A=[1,2,3;4,5,6;7,8,9]（2）按【Enter】键，指令被执行。（3）在指令执行后，MATLAB指令窗中将显示以下结果：A=123456789目前二十七页\总数七十七页\编于十一点阵列的运算符号

除了加减符号外其余的阵列运算符号均须多加.符号。阵列运算功能(注意:一定要多加.符号)+加-减.*乘./左除.\右除.^次方.’转置矩阵的几种基本变换操作（1）

通过在矩阵变量后加’的方法来表示转置运算>>a=[10,2,12;34,2,4;98,34,6];>>a'ans=10349822341246目前二十八页\总数七十七页\编于十一点（2）

矩阵求逆>>inv(a)ans=-0.01160.0372-0.00150.0176-0.10470.03450.0901-0.0135-0.0045（3）

矩阵的特征值>>[u,v]=eig(a)u=-0.29600.3635-0.3600-0.2925-0.41280.7886-0.9093-0.83520.4985v=48.8395000-19.8451000-10.9943目前二十九页\总数七十七页\编于十一点（4）

取出上三角和下三角>>triu(a)ans=10212024006>>tril(a)ans=1000342098346（5）LU分解>>[l,u]=lu(a)l=0.10200.15001.00000.34691.000001.000000u=98.000034.00006.00000-9.79591.91840011.1000目前三十页\总数七十七页\编于十一点MatLab图形绘制功能基本绘图函数命令含义plot建立向量或矩阵各队队向量的图形loglogx、y轴都取对数标度建立图形semilogxx轴用于对数标度，y轴线性标度绘制图形semilogyy轴用于对数标度，x轴线性标度绘制图形title给图形加标题xlabel给x轴加标记ylabel给y轴加标记text在图形指定的位置上加文本字符串gtext在鼠标的位置上加文本字符串grid打开网格线目前三十一页\总数七十七页\编于十一点x=0:0.001:10;%0到10的1000个点的x座标y=sin(x);%对应的y座标plot(x,y);%绘图（可以设定显示范围）目前三十二页\总数七十七页\编于十一点画椭圆

a=0:pi/50:2*pi; X=cos(a)*3;

%参数方程Y=sin(a)*2;plot(X,Y);xlabel('x'),ylabel('y');title('椭圆')目前三十三页\总数七十七页\编于十一点线条和颜色控制符:.点y黄。小圆圈m棕色xx标记c青色+加号r红色-实线g绿色*星号b蓝色：虚线w白色-.点划线h黑色目前三十四页\总数七十七页\编于十一点条形图当资料点数量不多时，长条图是很适合的表示方式：closeall;%关闭所有的图形视窗x=1:10;y=rand(size(x));bar(x,y)目前三十五页\总数七十七页\编于十一点误差图x=linspace(0,2*pi,30);y=sin(x);e=std(y)*ones(size(x));%std:标准差errorbar(x,y,e)目前三十六页\总数七十七页\编于十一点极坐标图形theta=linspace(0,2*pi);r=cos(4*theta);polar(theta,r);目前三十七页\总数七十七页\编于十一点x=linspace(0,10,50);y=sin(x).*exp(-x/3);stairs(x,y);%阶梯图x=linspace(0,10,50);y=sin(x).*exp(-x/3);stem(x,y);%针状图目前三十八页\总数七十七页\编于十一点符号表达式画图在许多的场合，将表达式可视化是有利的。MATLAB提供了函数ezplot来完成该任务。>>y='16*x^2+64*x+96'y= 16*x^2+64*x+96>>ezplot(y)目前三十九页\总数七十七页\编于十一点三维曲线和曲面

1．空间曲线绘制函数plot3()格式一：plot3(x,y,z)功能：若x、y、z为相同长度的向量，则根据向量x、y、z绘制空间三维曲线。如果x、y、z为同阶矩阵，则绘制对应列的多条曲线。格式二：plot3(x,y,z,’s’)功能：按字符串s设置的线型、颜色、标记符号绘制三维空间曲线。字符串s的设置与plot命令的设置相同。2．三维网线图函数mesh、meshc和meshzMATAL提供的函数mesh用来绘制三维网线图。格式一：mesh(X,Y,Z)功能：根据矩阵X、Y和Z绘制彩色的空间三维网线图。X、Y和Z中对应的元素为三维 空间上的点，点与点之间用线连接。其中网线的颜色随着网点高度的改变而改 变。格式二：mesh(x,y,Z)功能：n维向量x、m维向量y和m×n矩阵Z绘制网线图，节点的坐标为 (x(j),y(i),Z(i,j))，网线的颜色随着网点高度的改变而改变。目前四十页\总数七十七页\编于十一点三维曲面图函数surf及surfc

格式一：surf(X,Y,Z,C)功能：绘制由四个矩阵所指定的带色参数的网状表面图。视角是由view所指定。轴的刻度决定于X，Y及Z的范围，或当前对轴的设定；颜色范围由C指定。格式二：surf(X，Y，Z)功能：将C设为与Z相等，则颜色与网的高度成正比。目前四十一页\总数七十七页\编于十一点三维螺旋线t=0:pi/50:10*pi;plot3(sin(t),cos(t),t)grid目前四十二页\总数七十七页\编于十一点空间中的曲线t=linspace(0,20*pi,501);plot3(t.*sin(t),t.*cos(t),t);%注意用点乘.*t=linspace(0,10*pi,501);plot3(t.*sin(t),t.*cos(t),t,t.*sin(t),t.*cos(t),-t);目前四十三页\总数七十七页\编于十一点立体网状图x=linspace(-2,2,25);%在x轴上取25点y=linspace(-2,2,25);%在y轴上取25点[xx,yy]=meshgrid(x,y);%xx和yy都是21x21的矩阵zz=xx.*exp(-xx.^2-yy.^2);%计算函数值，zz也是21x21的矩阵mesh(xx,yy,zz);%画出立体网状图

目前四十四页\总数七十七页\编于十一点网状图填补x=linspace(-2,2,25);%在x轴上取25点y=linspace(-2,2,25);%在y轴上取25点[xx,yy]=meshgrid(x,y);%xx和yy都是21x21的矩阵zz=xx.*exp(-xx.^2-yy.^2);%计算函数值，zz也是21x21的 矩阵surf(xx,yy,zz);%画出立体曲面图目前四十五页\总数七十七页\编于十一点peaks函数peaks函数，可产生一个凹凸有致的曲面，包含了三个局部极大点及三个局部极小点，其方程式为：[x,y,z]=peaks;waterfall(x,y,z);目前四十六页\总数七十七页\编于十一点等高线图meshc同时画出网状图与等高线：[x,y,z]=peaks;meshc(x,y,z);surfc同时画出曲面图与等高线：[x,y,z]=peaks;surfc(x,y,z);目前四十七页\总数七十七页\编于十一点剔透玲珑球[X0,Y0,Z0]=sphere(30);%产生单位球面的三维坐标X=2*X0;Y=2*Y0;Z=2*Z0;%产生半径为2的球面的三维坐标surf(X0,Y0,Z0);%画单位球面shadinginterp%采用插补明暗处理holdon;mesh(X,Y,Z);holdoff%画外球面hiddenoff%产生透视效果axisoff%不显示坐标轴目前四十八页\总数七十七页\编于十一点彗星状轨迹图

【*例】简单二维示例。（请指令窗中运行以下指令）n=10;t=n*pi*(0:0.0005:1);x=sin(t);y=cos(t);plot(x,y,'g');holdon;comet(x,y,0.02);holdoff目前四十九页\总数七十七页\编于十一点卫星返回地球的运动轨线示意R0=1; %以地球半径为一个单位a=12*R0;b=9*R0;T0=2*pi; %T0是轨道周期T=5*T0;dt=pi/100;t=[0:dt:T]';f=sqrt(a^2-b^2); %地球与另一焦点的距离th=12.5*pi/180; %卫星轨道与x-y平面的倾角E=exp(-t/20); %轨道收缩率

x=E.*(a*cos(t)-f);y=E.*(b*cos(th)*sin(t));z=E.*(b*sin(th)*sin(t));plot3(x,y,z,'g') %画全程轨线[X,Y,Z]=sphere(30);X=R0*X;Y=R0*Y;Z=R0*Z;%获得单位球坐标gridon,holdon,surf(X,Y,Z),shadinginterp %画地球x1=-18*R0;x2=6*R0;y1=-12*R0;y2=12*R0;z1=-6*R0;z2=6*R0;axis([x1x2y1y2z1z2]) %确定坐标范围view([11737]),comet3(x,y,z,0.02),holdoff %设视角、画运动轨线

目前五十页\总数七十七页\编于十一点求解线性方程组MATLAB语句实现：a=[2,3,-1；8,2,3；45,3,9]；b=[2；4；23]；x=inv(a)*b目前五十一页\总数七十七页\编于十一点求多项式方程的根roots(p)solve(p)例.求方程的所有根p=[1-49-10]r=roots(p)或s1=sym('x^3-4*x^2+9*x-10');solve(s1)目前五十二页\总数七十七页\编于十一点求方程组的根和超越方程解[x,y]=solve('x+y=1','x-11*y=5','x','y')求超越方程的根，其中p，r为常数ff=sym('p*sin(x)=r')solve(ff,'x')目前五十三页\总数七十七页\编于十一点求极限symsx;y1=(1+4*x)^(1/x);y2=(exp(x)-1)/x;limit(y1,x,0)limit(y2,x,0)目前五十四页\总数七十七页\编于十一点求导数与偏导数

symsx;f=3*x^3+5*x+1;diff(f,1)

symsxy;%对x偏导数z=x^2*sin(2*y);B=diff(z,x)目前五十五页\总数七十七页\编于十一点求解不定积分symsxb;f=cos(b*x);f1=int(f,x)f1=1/b*sin(b*x)ezplot(f1)目前五十六页\总数七十七页\编于十一点求定积分多重积分

symsx;y=(x^2+sin(x))/(1+x^2);int(y,x,-1,1)

symsxy;int(int(x^2+y^2+1,y,x,x+1),x,0,1)目前五十七页\总数七十七页\编于十一点求微分方程通解和特解dsolve('Dy=1+y^2','x')dsolve('Dy=1+y^2','y(0)=1','x')dsolve('x^2*D2y+x*Dy+(x^2-1/2)*y=0','y(pi/2)=2','Dy(pi/2)=-2/pi','x')目前五十八页\总数七十七页\编于十一点MATLAB语言编程介绍MATLAB中各种命令可以完成许多单一的任务，对于某些较为复杂的问题，仅靠现有的命令或函数来解决，往往是难以达到目的。为此，要运用MATLAB编程语言编制程序，形成M-文件。程序是使计算机完成各项运算的命令集，运行一个编制好的程序，计算机会从第一条命令行开始，一行接一行地执行相应的命令，直到终止。程序一般包含：数据输入，数据处理和结果输出三部分，其中数据处理是核心。程序编写调试完成后，需要存盘，形成永久性文件，可以随时对它进行调用或修改。文件名以字母开头，但不能用专用变量名，如pi等。MATLAB中每一个命令都是一个M-文件。目前五十九页\总数七十七页\编于十一点MATLAB程序类型脚本文件函数文件functiony=mean(x)%MEANAverageormeanValue%Forvectors,Mean(x)isthemeanvalueofX%Formatrices,Mean(x)isarowvectorcontaining%themeanvalueofeachcolumn[m,n]=size(x)ifm==1m=n;endy=sum(x)/m;函数定义行H1行函数帮助文件函数体目前六十页\总数七十七页\编于十一点M-文件中输入、输出命令

直接赋值输入：提示对话输入(input命令)x=2;x=input('请输入参数x=');a=[1,2;3,4];a=input('请输入矩阵a=');s='anystring!';s=input('Pleaseinputs=');直接输出：格式控制输出(fprintf命令)xfprintf('x=%.0f,y=%.5f\n',pi,pi);disp([a,b]);fprintf('x=%5g,y=%10.5g\n',N,pi);目前六十一页\总数七十七页\编于十一点M-文件中循环控制命令(for命令)格式：fori=n1:(step):n2 commands;end作用：重复执行命令集commands.目前六十二页\总数七十七页\编于十一点例2：求奇数和：s=1+3+5+…+(2k-1)n=input('请输入n='); s=0;fori=1:2:n s=s+i;fprintf('i=%.0f,s=%.0f\n',i,s);end例1：求和：s=1+2+…+nn=input('请输入n='); s=0;fori=1:n s=s+i;fprintf('i=%.0f,s=%.0f\n',i,s);end目前六十三页\总数七十七页\编于十一点例3：求阶乘：p=1×2×3×…×n=n!n=input('请输入n='); p=1;fori=1:n p=p*i;fprintf('i=%.0f,p=%.0f\n',i,p);end例4：求e：e=1+1+1/2!+1/3!+…+1/n!n=input('请输入n='); p=1;e=1;fori=1:n p=p*i;p1=1/p;e=e+p1;fprintf('i=%.0f,p=%.0f,e=%.8f\n',i,p,e);end目前六十四页\总数七十七页\编于十一点例5：求pi：pi/4=1-1/3+1/5-1/7+1/9-1/11+……n=input('请输入n='); pi=0;k=-1;fori=1:2:n k=(-1)*k;pi=pi+k*1/i;fprintf('i=%.0f,pi=%.8f\n',i,4*pi);end目前六十五页\总数七十七页\编于十一点M-文件中条件循环命令(while命令)格式：while(conditionistrue) commands;end作用：重复执行命令集commands.目前六十六页\总数七十七页\编于十一点例6：求和：s=1+2+…+n<3000 s=0;k=1;whiles<3000 s=s+k;fprintf('k=%.0f,s=%.0f\n',k,s);k=k+1;end例7：求e：e=1+1+1/2!+1/3!+…(1/n!<1.0e-8) p=1;e=1;r=1;i=1;whiler>=1.0e-8 p=p*i;r=1/p;e=e+r;fprintf('i=%.0f,r=%.9f,e=%.9f\n',i,r,e);i=i+1;end目前六十七页\总数七十七页\编于十一点M-文件中选择控制命令(if命令)单项选择控制格式：if(conditionistrue) commands;end作用：若条件成立，则执行命令集commands.否则，不执行。目前六十八页\总数七十七页\编于十一点例8：求n个实数中最大的数M.a=input('请输入数组a[n]=');[m,n]=size(a);M=a(1); fori=2:n ifM<a(i)M=a(i);endfprintf('M=%.5f,a(%.0f)=%.5f\n',M,i,a(i));end目前六十九页\总数七十七页\编于十一点M-文件中选择控制命令(if命令)多项选择控制格式：if(conditionistrue) commands;elseif(conditionistrue)commands;elsecommands;end作用：若条件成立，则执行命令集commands.否则，不执行。目前