数学实验与 Matlab

数学实验与Matlab 制作人周晓阳华中科技大学数学系 zxyhust@netease.com http://math.hust.edu.cn/~zhou/index.html

实验一 矩阵运算与Matlab命令

日常矩阵及其运算 • 矩阵应用实例: • 榄球防护用品的生产管理

应用问题 • 一个工厂生产三种橄榄球用品：防护帽、垫肩、臀垫。 • 需要不同数量的：硬塑料、泡沫塑料尼龙线、劳动力。 • 为监控生产，管理者对它们之间的关系十分关心。 • 为把握这些量的关系，他列出下面的表

原料产品关系表

订单 • 管理者接到四份订单如上表所示。 • 问应该如何计算每份订单所需的原材料，以便组织生产？

将表格写成矩阵形式

计算 输入下面Matlab指令 • A=[4 2 3;1 3 2;1 3 3;3 2 2], • B=[35 20 60 45;10 15 50 40;20 12 45 20] • C=A*B • 请自行计算观看结果

Matlab基本指令 向量的创建和运算

1. 直接输入向量 》x1=[1 2 4],x2=[1,2,1],x3=x1’ • 运行结果 • x1 = 1 2 4 x2 = 1 2 1 x3 = 1 2 4

2.冒号创建向量 • x1=3.4:6.7, x2=3.4:2:6.7 x3=2.6:-0.8:0 • 运算结果 x1 = 3.4000 4.4000 5.4000 6.4000 x2 = 3.4000 5.4000 x3 = 2.6000 1.8000 1.0000 0.2000

3.生成线性等分向量 指令x=linspace(a,b,n) 在[a,b]区间产生 n 个等分点(包括端点) • x=linspace(0,1,5) • 结果 • x = • 0 0.2500 0.5000 0.7500 1.0000

工作空间 • 在Matlab窗口创建向量后并运行后，向量就存在于工作空间，可以被调用。

向量的运算 • 设x=[x1 x2 x3]; y=[y1 y2 y3];为两个三维向量，a,b为标量。 • 向量的数乘：a*x=[a*x1 a*x2 a*x3] • 向量的平移： x+b=[x1+b x2+b x3+b] • 向量和： x+y=[x1+y1 x2+y2 x3+y3] • 向量差： x-y=[x1-y1 x2-y2 x3-y3] • 数的乘幂：如 a^2

元素群运算(四则运算） • x.*y=[x1*y1 x2*y2 x3*y3] (元素群乘积) • x./y=[x1/y1 x2/y2 x3/y3] (元素群右除，右边的y做分母) • x.\y=[y1/x1 y2/x2 y3/x3] (元素群左除，左边的x做分母) • x.^5=[x1^5 x2^5 x3^5] (元素群乘幂) • 2.^x=[2^x1 2^x2 2^x3] (元素群乘幂) • x.^y=[x1^y1 x2^y2 x3^y3] (元素群乘幂)

元素群运算（函数计算） • Matlab有许多内部函数，可直接作用于向量产生一个同维的函数向量。 • x=linspace(0,4*pi,100);（产生100维向量x） • y=sin(x); (y也自动为100维向量) • y1=sin(x).^2; • y2=exp(-x).*sin(x); • 观察结果

创建矩阵（数值矩阵的创建） • 直接输入法创建简单矩阵。 • A=[1 2 3 4; 5 6 7 8; 9 10 11 12] • B=[-1.3,sqrt(3);(1+2)*4/5,sin(5);exp(2),6] • 观察运行结果

创建矩阵（符号矩阵的创建） • 用指令“syms”说明符号变量。 • syms a11 a12 a13 a14 a21 a22 a23 a24 a31 a32 a33 a34 b11 b12 b13 b14 b21 b22 b23 b24 b31 b32 b33 b34 • A1=[a11 a12 a13 a14 ;a21 a22 a23 a24; a31 a32 a33 a34], • B1=[b11 b12 b13 b14 ;b21 b22 b23 b24; b31 b32 b33 b34] • 运行

矩阵的运算(矩阵的加减、数乘、乘积) • C=A1+B1 • D=A1-B1 • syms c, cA=c*A1 • A2=A1(:,1:3), B1 • G=A2*B1

矩阵的运算(矩阵的加减、数乘、乘积) • A, A_trans=A’ • H=[1 2 3 ; 2 1 0 ; 1 2 3 ], • K=[1 2 3 ; 2 1 0 ; 2 3 1] • h_det=det(H), k_det=det(K), • H_inv=inv(H),K_inv=K^-1

矩阵的运算(左除和右除) • 左除“ \ ”：求矩阵方程AX=B的解；（A 、B的行要保持一致）解为 X=A\B；当A为方阵且可逆时有X=A\B=inv(A)*B； • 右除“ / ”：求矩阵方程XA=B的解（A 、B的列要保持一致）解为 X=B/A ，当A为方阵且可逆时有X=B/A=B*inv(A)

矩阵的运算(左除和右除) • 求矩阵方程：设A、B满足关系式：AB＝2B+A,求B。其中A=[3 0 1; 1 1 0; 0 1 4]。 • 解：有(A-2I)B＝A • 程序： A=[3 0 1; 1 1 0;0 1 4]; B=inv(A-2*eye(3))*A, B=(A-2*eye(3))\A • 观察结果：

分块矩阵（矩阵的标识） • 1.矩阵元素的标识： A(i,j)表示矩阵A 的第 i 行 j 列的元素； • 2.向量标识方式 A(vr,vc)： vr=[i1,i2,…,ik]、vc=[j1,j2,…,ju]分别是含有矩阵A的行号和列号的单调向量。 A(vr,vc)是取出矩阵A的第i1,i2,…,ik行与j1,j2,…,ju列交叉处的元素所构成新矩阵。

分块矩阵（矩阵的标识） • 取出A的1、3行和1、3列的交叉处元素构成新矩阵A1 • 程序 • A=[1 0 1 1 2;0 1 -1 2 3; • 3 0 5 1 0;2 3 1 2 1], • vr=[1, 3];vc=[1, 3]; • A1=A(vr, vc) • 观察结果

分块矩阵（矩阵的标识） • 将A分为四块，并把它们赋值到矩阵B中，观察运行后的结果。 • 程序 • A11=A(1:2,1:2),A12=A(1:2,3:5), • A21=A(3:4,1:2),A22=A(3:4,3:5) • B=[A11 A12;A21 A22] • 结果

分块矩阵（矩阵的修改和提取） • 修改矩阵A，将它的第1行变为0。 • 程序： • A=[1 0 1 1 2;0 1 -1 2 3; • 3 0 5 1 0;2 3 1 2 1] • A(1,:)=[0 0 0 0 0]; A • 删除上面矩阵A的第1、3行。 • 程序： • A([1,3],:)=[ ] • 结果

生成特殊矩阵 • 全1阵 • ones(n), ones(m,n), ones(size(A)) • 全零阵： • zeros(n) ，zeros(m,n), zeros(size(A)) • 常常用于对某个矩阵或向量赋0初值 • 单位阵： • eye(n)，eye(m,n) • 随机阵： • rand(m,n)， rand(n)=rand(n,n)用于随机模拟，常和rand('seed',k)配合使用。

生成特殊矩阵 • 将rand指令运行多次，观察结果。 • 程序： • y1=rand(1,5), y2=rand(1,5), • rand('seed',3), x1=rand(1,5), • rand('seed',3), x2=rand(1,5) • 结果

常用矩阵函数 • det(A) ：方阵的行列式； • rank(A)：矩阵的秩； • eig(A)：方阵的特征值和特征向量； • trace(A)：矩阵的迹； • rref(A)：初等变换阶梯化矩阵A • svd(A)：矩阵奇异值分解。 • cond(A)：矩阵的条件数；

数据的简单分析 • 1.当数据为行向量或列向量时，函数对整个向量进行计算. • 2.当数据为矩阵时，命令对列进行计算，即把每一列数据当成同一变量的不同观察值。 • max(求最大)、min(求最小)、mean(求平均值)、sum(求和)、std(求标准差)、cumsum(求累积和)、median(求中值)、diff(差分)、sort(升序排列)、sortrows(行升序排列)等等。

数据的简单分析 • 观察：生成一个3×6的随机数矩阵，并将其各列排序、求各列的最大值与各列元素之和。 • 程序 • rand('seed',1);A=rand(3,6), • Asort=sort(A), Amax=max(A), Asum=sum(A) • 结果

实验二 函数可视化与Matlab作图

函数的可视化 • 令 • f (x), g (x)是周期函数吗？观察它们的图象。 • 程序 • clf, x=linspace(0,8*pi,100); • F=inline('sin(x+cos(x+sin(x)))'); • y1=sin(x+cos(x+sin(x))); y2=0.2*x+sin(x+cos(x+sin(x))); • plot(x,y1,'k:',x,y2,'k-') legend('sin(x+cos(x+sin(x))','0.2x+sin(x+cos(x+sin(x)))',2)

绘制平面曲线(plot指令) • plot(x,y)： • 以x为横坐标、y为纵坐标绘制二维图形 • x,y是同维数的向量； • plot(y)： • 相当于x=[1,2,…,length(y)]时情形。

绘制平面曲线（绘制多个图形) • 1. plot(x,[y1;y2;…]), • x是横坐标向量，[y1;y2;…]是由若干函数的纵坐标拼成的矩阵 • 2. plot(x,y1), hold on, plot(x,y2), hold off • 3. plot(x,y1,x,y2,…) • 4.plotyy • 两个坐标系，用于绘制不同尺度的函数。

绘制平面曲线（线型、点形和颜色的控制） • plot(x,y,‘颜色＋线型＋点形’) • plot(x,y,‘颜色＋线型＋点形’,x,y,‘颜色＋线型＋点形’,… ) • 句柄图形和set命令改变属性值，可套用： • h=plot(x,y), • set(h,‘属性’,属性值,‘属性’,属性值,…) • 也可用plot(x,y,'属性','属性值')设置图形对象的属性。

绘制平面曲线（属性变量和属性值） • 线宽：LineWidth • 点的大小： MarkerSize • 线型：LineStyle • 颜色：color

绘制平面曲线（例） • 观察： • 改变绘图的线型和颜色。 • 用grid on 指令为图形窗口加上网格线，并改变网格的线型和字体的大小。 • 程序 • h=plot([0:0.1:2*pi],sin([0:0.1:2*pi])); • set(h,'LineWidth',5,'color','red'); grid on • set(gca,'GridLineStyle','-','fontsize',16) • 观察结果

绘制平面曲线（坐标轴的控制） • axis指令 • axis([xmin xmax ymin ymax])：设定二维图形的x和y坐标的范围； • axis([xmin xmax ymin ymax zmin ymax])：设定三维图形的坐标范围； • 其中xmin<x<xmax， ymin<y<ymax ，zmin<z<zmax。

绘制平面曲线（gca属性控制） • 改变当前轴对象句柄gca属性 • 用set(gca,‘属性’,属性值,…)可改变字体大小、坐标刻度等轴对象的内容。例如： • set(gca,'ytick',[-1 -0.5 0 0.5 1]) • 将 y 坐标按向量[-1 -0.5 0 0.5 1]将刻度分成4格； • set(gca,'yticklabel','a|b|c|d|e') • 改变y坐标刻度的说明。

绘制平面曲线（gca属性控制,例） • 设置y坐标的刻度并加以说明，并改变字体的大小。 • 程序 • plot([0:0.1:2*pi],sin([0:0.1:2*pi]),'k.-',);grid on, • axis([0 6.3 -1.1 1.1])， • set(gca,'ytick',[-1 -0.5 0 0.5 1]), • set(gca,'yticklabel','a|b|c|d|e'), • set(gca,'fontsize',20) • get(gca) • 运行结果

绘制平面曲线（文字标注） • title(‘图形标题’)； • xlabel（‘x轴名称’）；ylabel（‘y轴名称’）；zlabel（‘z轴名称’）； • text(‘说明文字’)：创建说明文字； • gtext（'说明文字'）：用鼠标在特定位置输入文字。 • 文字标注常用符号： • \pi （π）；\alpha （α）；\beta （β）； • \leftarrow （左箭头） \rightarrow （右箭头）； • \bullet （点号）

绘制平面曲线（程序讲解，exp2_1.m） • clf, t=0:0.1:3*pi;alpha=0:0.1:3*pi; • plot(t,sin(t),'r-');hold on; plot(alpha,3*exp(-0.5*alpha),'k:'); • set(gca,'fontsize',15,'fontname','times New Roman'), • xlabel('\it{t(deg)}'); • ylabel('\it{magnitude}'); • title(' \it{sine wave and {\it{Ae}}^{-\alpha{\itt}}wave}');

绘制平面曲线（程序讲解，exp2_1.m） • text(6,sin(6),'\fontsize{15}The Value \it{sin(t)} at {\itt}=6\rightarrow\bullet', 'HorizontalAlignment','right'), • text(2,3*exp(-0.5*2), • ['\fontsize{15}\bullet\leftarrow The Value of \it{3e}^{-0.5 \it{t}}=', • num2str(3*exp(-0.5*2)),' at \it{t} =2 ']); • legend('\itsin(t)','{\itAe}^{-\alphat}') • 注1： num2str： ['string1' ,num2str,'string2']，用方括号 • 注2：legend 请结合图形观察此命令的使用

图形窗口的创建和分割 • subplot(m,n,k)命令。 • 在图形区域中显示多个图形窗口。 • m为上下分割数，n为左右分割数，k为第k子图编号。 • 例：将一个图形分为9个子图，在第k个子图画sin(kx) 的图象. • 程序： • clf，b=2*pi;x=linspace(0,b,50); • for k =1:9 • y=sin(k *x); • subplot(3,3,k),plot(x,y),axis([0,2*pi,-1,1]) • end

若干有用的指令 • clf：清除图形窗口已有的内容. • shg：显示图形窗口。 • clear、 clear x：清除工作空间的已有变量。 • figure(n): 打开第n个图形窗口 • help: • … : 续行号

绘制二元函数 • 基本步骤： • 1.生成二维网格点 • 2. 计算函数在网格点上的值 • 3. 绘制函数图形

三维绘图（ meshgrid指令：生成网格点） • 观察meshgrid指令的效果。 • 程序： • a=-0.98;b=0.98;c=-1;d=1;n=10; • x=linspace(a,b,n); y=linspace(c,d,n); • [X,Y]=meshgrid(x,y); • plot(X,Y,'+') • 观察结果

三维绘图（计算函数值，定义域裁减） • 程序： • for i=1:n • for j=1:n • if (1-X(i,j))<eps1|X(i,j)-Y(i,j)<eps1 • z(i,j)=NaN; • else • z(i,j)=1000*sqrt(1-X(i,j))^-1.*log(X(i,j)-Y(i,j)); • end • end • end

三维绘图（绘图指令） • mesh(X,Y,z) ： • 在三维空间中绘出由(X,Y,z)表示的曲面； • meshz(X,Y,z)： • 除了具有mesh的功能外，还画出上下高度线， • meshc(X,Y,z)： • 除了具有mesh的功能外，还在曲面的下方画出函数z=f(x,y)的等值线图， • surf(X,Y,z)： • 也是三维绘图指令，与mesh的区别在于mesh绘出彩色的线，surf绘出彩色的面， • 运行exp2_1,观察效果

数学实验与 Matlab