实验六-数值积分与数值微分.ppt

数值实验六 数值积分与数值微分 1 注意 严禁在试验课上打游戏、登 录 外网、QQ等与上课无关行为 。 一经发现，从重扣除平时成 绩 2 上课时间 n10-16周 , 周三下午 n14：00 16：30 n(课间不休) 3 实验报告要求 n实验报告在课后一周内上交； n算法使用编程语言不限； n实验报告可以提交书面试验报告，也可以 提交电子版试验报告； n电子版试验报告与程序代码发送至服务器 如下文件夹： 张昆学生作业实验x n或发送至邮箱： n文件名称：实验三 08119000 张三 4 实验四题目 n参考题目(选做二题） n 数值试验七 （1） （p280） n 数值试验七（2） （p281） 5 Matlab基础（实验一） n基本命令 n基本数学运算 n符号解 n函数文件格式 n函数调用 n文件保存位置 6 基本命令 1. help 帮助命令 help format %查看format命令的帮助； 2. clear 清除变量命令 clear x %清除变量x(无参数时清理工作空间 ） 3. 清理命令窗口 clc %清屏 4. format 输出格式 format long %15位有效数字格式输出 format short %4-7位有效数字格式输出 format compact %紧凑格式(输出不加空行） 7 基本数学运算 n乘法：* n除法：/ n乘方： n根号：sqrt(x) n正弦：sin(x) n余弦：cos(x) n自然对数：log(x) n自然指数：exp(x) 8 1. vpa 控制运算精度 vpa(表达式,运算精度) vpa(pi,20) %显示至100位有效数字 2. solve 方程的符号解 solve( 表达式 , x) 符号解 9 符号解 n例 解方程 n解 在MATLAB工作窗口输入命令: ny=solve(x3-sin(x)-12*x+1=0,x) nvpa(y,5） ny=roots(2,0,-1,-1) 10 函数文件格式 函数文件由function语句引导，其格式为： function 输出形参表=函数名(输入形参表) % 注释说明部分 函数体： 例： function k,xk,yk,p=jhnewtonqx(x0,ddmax) % 牛顿切线法求非线性方程的根 注：其中函数名的命名规则与变量名相同。输入形参为函 数的输入参数，输出形参为函数的输出参数。当输出形 参多于1个时，则应该用方括号括起来。 11 函数调用 函数调用的一般格式是： 输出实参1,实参2,.=函数名(输入实参1,实参 2,.） 例： k,xk,yk,p=jhnewtonqx(1.5,20) 12 文件保存位置 设置当前目录位置 cd C:Matlab7work %设置Matlab当前目录为： C:Matlab7work 13 Matlab基础（实验二） n1. 矩阵输入 n2. 特殊矩阵 n3. 误差分析 n4. 内置函数 14 矩阵输入 n a= %建立空矩阵a ,可在workspace中 编辑 n a=1,2,3;4 5 6% 分行“,”或空格 分列 “;” n a=1:5% 建立15矩阵(向量) a=(1,2,3,4,5) n a=0:pi/5:2*pi n% (pi= ) 建立111矩阵（向量）a ：初始元0，步长/5 , 终止元2 n 【注意 “ ：”的用法】 15 常用矩阵 n1、零矩阵 n z=zeros(3,4)%产生 34零矩阵z； n zeros(5)%产生 5阶零矩阵； n2、单位矩阵 n E=eye(6)%产生6阶单位矩阵E； n2、幺矩阵 n ones(3,2);%产生3行2列元素为1的矩阵； n3、随机矩阵 n r=rand(7)%产生01间分布的随机矩阵r n s=round(rand(7)*30) %产生7阶030间均匀分部的随机矩 阵s 16 特殊矩阵 1、Hilbert矩阵 h=hilb(5)%产生 5阶Hilbert矩阵h h=sym(hilb(5)%产生 5阶Hilbert符号矩阵h 2、Vandermonde矩阵 v=vander(1,2,3,4,5)%产生5阶Vandermonde矩阵v 3、魔方矩阵 m=magic(3)%产生3阶魔方矩阵 3、Toeplitz矩阵 t=toeplitz(0:-1:-5,0:5) 17 矩阵运算 a=magic(3);b=round(10*rand(3);r=(1:3); n a+b%矩阵加法 n a-b%矩阵减法 n a*b%矩阵乘法 n n=inv(a)% inv(a) =a(-1)= a的逆 n a/b% a/b=a*inv(b)=a*b(-1) n x=ar,a*x% x=ar=inv(a)*r=a(-1)*r n c=a% a= a的（共轭）转置 18 元素运算 a=magic(3);b=round(10*rand(3);r=(1:3); n a.*b%矩阵对应元素乘积 n a b a.*b,a*b%比较 a.*b,a*b n a./b%矩阵对应元素右除 n a.b%矩阵对应元素左除 n a.(1/3)%矩阵对应元素的立方根 19 矩阵索引 n a=1:30%产生120行向量 n a=reshape(a,5,6)% 变更a的结构为35的矩阵 n a(3,2)% 取元素a（3,2） n a(3,:)% 取a第3行，取a第二 列 n a(:,2)% 取a第2列 n b=a(1,3,3:5)%取a的1,2行，2,4,5列元 n m,n=size(a)%输出a的行列数m、n n b=a(end,:) % b取a的最后一行 20 矩阵操作 n a=1:36; a=reshape(a,6,6); t=toeplitz(0:-1:-5,0:5); n c=diag(a)%提取a的列向量 n b=diag(c)% b为以c为对角元的对角阵 n d=diag(diag(a) n d=diag(diag(t,-1),-1) n l=tril(a)%l为a的下三角阵 n u=triu(a,1) n% u为a从第1条对角线的开始取的上三角阵 n l=tril(a,-2) n% l为a从第-2条对角线的开始取的下三角阵 21 矩阵分析 n rank(a)% a的秩 n norm(a)% a的2-范数 n norm(a,inf)% a的无穷范数 n cond(a)% a的谱（2-）条件数 n cond(a,1)% a的1-条件数 22 矩阵分解与特征值 n L,U = lu(a)%矩阵a的LU分解 n Q,R = qr(a)%矩阵a的QR分解 n eig(a) %矩阵a的特征值向量 n V,D = eig(a) n%D主对角线元素为A的全部近似特征值 n%V第k列元素为对应于A的特征值D(k,k)的特征向量 23 最大元搜索 n m=magic(5) n x,q=max(max(m) n%矩阵m每列最大元素向量构成x，x向量最大元下标 n C,p=max(m(:,q) n%向量m(:,q)最大元素C， m(:,q)最大元下标p n%例：搜索A的非对角线上绝对值最大元素A(p,q) A=abs(m- triu(m);x,q=max(max(A);. M,p=max(A(:,q);p,q,A(p,q) 24 条件搜索 n m=magic(5),n,n=size(m); n s=find(m24,1) n%搜索m中大于24的元素的前1个，返回此元素下标s n q=ceil(s/n)%计算第s个元素所在列 n p=s-n*(q-1)%计算第s个元素所在行 25 二维绘图 （一）plot 最基本的二维图形指令 plot的功能： n plot命令自动打开一个图形窗口Figure n 用直线连接相邻两数据点来绘制图形 n根据图形坐标大小自动缩扩坐标轴，将 数据标尺及单位标注自动加到两个坐标 轴上，可自定坐标轴，可把x, y 轴用对 数坐标表示 26 n如果已经存在一个图形窗口，plot命 令则清除当前图形，绘制新图形 n可单窗口单曲线绘图；可单窗口多 曲线绘图；可单窗口多曲线分图绘 图；可多窗口绘图 n可任意设定曲线颜色和线型 n可给图形加坐标网线和图形加注功 能 27 plot的调用格式 n plot(x) 缺省自变量绘图格式 ，x为向量, 以x元素值为纵坐标， 以相应元素下标为横坐标绘图 n plot(x,y) 基本格式，以y(x)的 函数关系作出直角坐标图，如果y 为nm的矩阵，则以x 为自变量， 作出m条曲线 nplot(x1,y1,x2,y2) 多条曲线绘 图格式 28 nplot(x,y,s) 开关格式，开关 量字符串s设定曲线颜色和绘图方 式，使用颜色字符串的前13个字 母，如 yellowyel表示等。 或plot(x1,y1,s1,x2,y2,s2，) 29 S的标准设定值 字母 颜色 标点 线型 y 黄色 点线 m 粉红 圈线 c 亮蓝 线 r 大红 字线 g 绿色 实线 b 蓝色 星形线 w 白色 ： 虚线 k 黑色 (-) 点划线 30 +、 o、*、 . 、 x、 、 v、 、 、 x=0:pi/100:2*pi; y=sin(x);y1=cos(x); plot(x,y,r:,x,y1,gp) 32 子图绘制 x=-3:0.1:3; y1=x;y2=x.2;y3=x.3;y4=x.4; subplot(2,2,1); plot(x,y1);title(y1=x) subplot(2,2,2) plot(x,y2);title(y2=x2) subplot(2,2,3) plot(x,y3);title(y3=x3) subplot(2,2,4) plot(x,y4);title(y4=x4) 33 图形标注 标题函数：tiltle x=0:0.1:5; y=exp(-0.2*x).*sin(x); plot(x,y) title(it e0.2xsin(x),FontWeight,Bold) 34 图形标注 坐标轴标注：xlabel , ylabel x=1990:2:2000; y=1.25 0.81 2.16 2.73 0.06 0.55; xin=1990:0.2:2000; yin=spline(x,y,xin); plot(x,y,ob,xin,yin,-.r) title(1990年到2000年某地区年平均降水量图) xlabel(it 年份,FontSize,15) ylabel(降雨量,FontSize,8) 35 图形标注 曲线标注：legend x=0:0.02*pi:2*pi; y1=sin(x); y2=cos(x); y3=sin(3*x).*cos(x); plot(x,y1;y2;y3) axis(0 2*pi -1.5 1.5) legend(sin(x),cos(x),sin(3x)cos(x) ,-1) 36 TeX字符 bf %粗体 it %斜体 rm %正常字体 fontname %指定字体 fontsize %指定字号 x %上标 _x %下标 pi ; lambda alpha ; beta ; 37 文本标注 text x=linspace(-1,1,500); f=exp(x); l2=0.9963+1.1036*x+0.5367*x.*x; l3=0.9963+0.9979*x+0.5367*x.*x+0.1761*x.*x.*x; plot(x,f,-b,x,l3,-r); xlabel(0 leq x leq1 );ylabel(y); title(bf ex的3次最佳平方逼近多项式,fontname,隶书,fontsize,12) legend(it f(x)=ex,it l_2=0.9963+0.9979*x+0.5367*x2+0.1761*x3,2) text(0.96,exp(0.96),it f(x) = ex ,rightarrow ,. HorizontalAlignment,right,fontsize,12,color,0,0,1,fontname,times new roman); text(0.6,0.9963+0.9979*0.6+0.5367*0.36+0.1761*0.36*0.6,it l_3=0.9963+0.9979*x+0.5367*x2+0.1761*x3,rightarrow,. HorizontalAlignment,right,fontsize,12,color,1,0,0,fontname,times new roman); 38 多项式函数 p=1:5;v=1:3; poly2sym(p) r=roots(p) poly(v) y=polyval(p,1) y=polyval(p,p) conv(p,v) q,r=deconv(p,v) 39 多项式插值 p=polyfit(x,y,n)%将数据点x，y拟合为n次多项 式p x=linspace(-5,5,11);y=5./(1+x.*x); p=polyfit(x,y,10); %10次插值多项式 t=linspace(-5,5); z=polyval(p,t);plot(t,z) hold on;x=linspace(-5,5,5);y=5./(1+x.*x); p=polyfit(x,y,4);%4次插值多项式 t=linspace(-5,5); z=polyval(p,t);plot(t,z,r) 40 多项式插值 p=polyfit(x,y,n)%将数据点x，y拟合为n次多项 式p x=linspace(-5,5,11);y=5./(1+x.*x); p=polyfit(x,y,10); t=linspace(-5,5); z=polyval(p,t);plot(t,z) yi = interp1(x,y,xi,method) %计算数据点(x,y)按method指定的插值函数在xi点的值，method参数 ： linear分段线性插值spline默认三次样条插值;cubic三次Hermite插值 ； x = 0:10; y = sin(x); xi = 0:0.25:10; yi = interp1(x,y,xi); %默认为一维线性插值 plot(x,y,o,xi,yi) 41 多项式插值 x = 0:10;y = sin(x);xi = 0:.25:10;yi=sin(xi); grid on y2