MATLAB语言与现代科学计算(还有整套课件大家可以参考我的文库).ppt

MATLAB程序语言设计,余义斌 五邑大学信息工程学院,,第四讲 MATLAB语言与现代科学计算,4.1 数值线性代数问题 4.2 数值微积分 4.3 统计分析 4.4 数学问题的解析运算 4.5 数学问题的高精度运算,4.1 数值线性代数问题,一些 MATLAB 特殊矩阵,零矩阵： A=zeros(m, n); 其中 （m, n） 定义零矩阵维数大小,全1矩阵： A=ones(m, n); 其中 （m, n） 定义矩阵维数大小,单位矩阵： A=eye(m, n); 其中 （m, n） 定义零矩阵维数大小,随机元素矩阵： A=rand(m, n); 0，1上均匀分布 A=randn(m, n); 正态分布,对角矩阵： A=diag(v); V 为对角向量,4.1 数值线性代数问题,Hilbert 矩阵：A=hilb(n);,The elements of the Hilbert matrices are: H(i, j) = 1/(i+j-1),伴随矩阵：A=compan(p); 其中 p 为多项式系数向量,多项式 对应的向量为 p= 1 2 7 6 9 8 ，它的伴随矩阵为 compan(p) ans = -2 -7 -6 -9 -8 1 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 1 0,4.1 数值线性代数问题,Hankel 矩阵：A=hankel(c, r);,A Hankel 矩阵是反对角线对称且为常数的矩阵，其分量是 h(i, j) = p(i+j-1), 向量p = c r(2:end) 可以完全决定Hankel 矩阵.,c = 1:3; r = 7:10; h = hankel(c,r) h = 1 2 3 8 2 3 8 9 3 8 9 10 p = 1 2 3 8 9 10,4.1 数值线性代数问题,范德蒙德矩阵：A=vander(c, r);,A = vander(v) 返回一个范德蒙德矩阵，其列向量是向量V的乘方，即 A(i,j) = v(i)(n-j), 其中 n = length(v). c=1:5, v=vander(c) c = 1 2 3 4 5 v = 1 1 1 1 1 16 8 4 2 1 81 27 9 3 1 256 64 16 4 1 625 125 25 5 1,4.1 数值线性代数问题,矩阵行列式： det(A),矩阵的迹： trace(A) (即对角线元素之和）,矩阵的秩： rank(A) 即线性无关的列数或行数,矩阵的特征多项式： poly(A) 即,矩阵的特征根 roots(poly(A),4.1 数值线性代数问题,矩阵的特征值与特征向量,V, D= eig(A) 矩阵V的各列为特征向量，矩阵D的对角元素为特征值, a=1 2 3;4 5 6;7 8 0; v,d=eig(a) v = -0.2998 -0.7471 -0.2763 -0.7075 0.6582 -0.3884 -0.6400 -0.0931 0.8791 d = 12.1229 0 0 0 -0.3884 0 0 0 -5.7345,4.1 数值线性代数问题,矩阵求逆与线性方程求解,B= inv(A) 求矩阵A的逆矩阵, a=1 2 3;4 5 6; 7 8 0; b=inv(a); c=a*b;d= b*a; b c d ans = -1.7778 0.8889 -0.1111 1.0000 0 -0.0000 1.0000 0.0000 0 1.5556 -0.7778 0.2222 -0.0000 1.0000 0 -0.0000 1.0000 0 -0.1111 0.2222 -0.1111 0.0000 -0.0000 1.0000 0.0000 0.0000 1.0000,4.1 数值线性代数问题,面向矩阵各个元素的函数,B= 函数名(A),用命令 help elfun 可以查看这些命令列表，主要有,sin sinh asin asinh cos cosh acos acosh tan tanh atan atan2 atanh sec sech asec asech csc csch acsc acsch cot coth acot acoth exp log log10 log2 pow2 sqrt abs angle complex conj imag real unwrap isreal cplxpair fix floor ceil round mod rem sign,4.1 数值线性代数问题,对矩阵进行数值分析的函数,B= 函数名(A),用命令 help datafun 可以查看这些命令列表，主要有,max min mean median std var sort sortrows sum prod hist histc trapz cumsum cumprod cumtrapz diff gradient del2 corrcoef cov subspace filter filter2 conv conv2 convn deconv detrend fft fft2 fftn ifft ifft2 ifftn fftshift ifftshift,4.2 数值微积分,数值差分运算 dy=diff(y) (按列运算）, y=magic(6) y = 35 1 6 26 19 24 3 32 7 21 23 25 31 9 2 22 27 20 8 28 33 17 10 15 30 5 34 12 14 16 4 36 29 13 18 11 diff(y) ans = -32 31 1 -5 4 1 28 -23 -5 1 4 -5 -23 19 31 -5 -17 -5 22 -23 1 -5 4 1 -26 31 -5 1 4 -5,4.2 数值微积分,数值积分 y, n=quad( F, a, b, tol); y, n=quadl(F, a, b, tol);,其中F为被积函数， a, b为积分上下限， tol 为误差限 y 为积分结果，n 为被积分函数的调用次数。,tol 的缺省值对 quad 函数为 1e-3 quadl 函数使用的算法比 quad 高得多,4.2 数值微积分,例子：求无穷定积分, f=inline(1/sqrt(2*pi)*exp(-x.2/2),x); y,kk=quad(f,-8,8) y = 1.00000197533430 kk = 81, y,kk=quad(f,-15,15) y = 0.99999920879563 kk = 89,该无穷定积分的理论值为 1,4.2 数值微积分,例子：双重积分,用函数 dblquad(F, xmin, xmax, ymin, ymax, tol) 求双重积分, f=inline(exp(-x.2/2).*sin(x.2+y),x,y); y=dblquad(f,-2,2,-1,1) y = 1.57450259485169,4.3 统计分析,数据统计分析,x, i= max(X), 或 x , i= min(X) 返回矩阵 X 各列的最大值或最小值， i 为各列最大值或最小值所在位置的行号构成的向量,x, i= sort(X) 对矩阵X各列排序,y= mean(X) 求矩阵X中各列向量的均值,y= std(X) 求矩阵X中各列向量的标准差,y= median(X) 求矩阵X中各列向量的中间值,y= cov(X) 求协方差,y= corrcoef(X) 求相关系数矩阵,4.3 统计分析,例子, s=randn(10000, 5); M=mean(s), D=std(s), MD=median(s) M = 0.0148 -0.0149 -0.0026 -0.0067 -0.0037 D = 0.9940 1.0024 1.0098 1.0044 1.0043 MD = 0.0250 -0.0104 0.0064 -0.0090 -0.0185 V=corrcoef(s) V = 1.0000 0.0047 0.0100 0.0197 -0.0079 0.0047 1.0000 0.0219 -0.0062 0.0009 0.0100 0.0219 1.0000 0.0021 0.0081 0.0197 -0.0062 0.0021 1.0000 -0.0028 -0.0079 0.0009 0.0081 -0.0028 1.0000,4.4 数学问题的解析运算,MATLAB本身没有提供解析运算，即符号运算的能力，但借助于著名的符号运算软件 MAPLE， MATLAB提供了一个符号运算工具箱，使得MATAB可以进行符号运算及进行高精度运算.,syms 变量列表 定义符号变量,s=sym(A) 从A中定义符号对象 s。,观察下例, a=hilb(5)，b=sym(hilb(5)，whos a = 1.0000 0.5000 0.3333 0.2500 0.2000 0.5000 0.3333 0.2500 0.2000 0.1667 0.3333 0.2500 0.2000 0.1667 0.1429 0.2500 0.2000 0.1667 0.1429 0.1250 0.2000 0.1667 0.1429 0.1250 0.1111 b = 1, 1/2, 1/3, 1/4, 1/5 1/2, 1/3, 1/4, 1/5, 1/6 1/3, 1/4, 1/5, 1/6, 1/7 1/4, 1/5, 1/6, 1/7, 1/8 1/5, 1/6, 1/7, 1/8, 1/9 Name Size Bytes Class a 5x5 200 double array b 5x5 2478 sym object,4.4 数学问题的解析运算,可以用命令 help symbolic 查看MATLAB提供的符号运算 命令，可以看出有很多，我们只介绍少数几个命令。, a=sym(hilb(4); b=inv(a); a b ans = 1, 1/2, 1/3, 1/4, 16, -120, 240, -140 1/2, 1/3, 1/4, 1/5, -120, 1200, -2700, 1680 1/3, 1/4, 1/5, 1/6, 240, -2700, 6480, -4200 1/4, 1/5, 1/6, 1/7, -140, 1680, -4200, 2800 c=a*b c = 1, 0, 0, 0 0, 1, 0, 0 0, 0, 1, 0 0, 0, 0, 1,矩阵求逆 inv(A),4.4 数学问题的解析运算,4.4 数学问题的解析运算,极限、微分、积分和 Taylor级数展开 limit, diff, int, taylor, syms x a t h; limit(sin(x)/x) ans =1 limit(x-2)/(x2-4),2) ans =1/4 limit(1+2*t/x)(3*x),x,inf) ans =exp(6*t) limit(1/x,x,0,right) ans =inf limit(1/x,x,0,left) ans =-inf limit(sin(x+h)-sin(x)/h,h,0) ans =cos(x) v = (1 + a/x)x, exp(-x);limit(v,x,inf,left) ans = exp(a), 0,4.4 数学问题的解析运算,极限、微分、积分和 Taylor级数展开 limit, diff, int, taylor, syms x t; diff(sin(x2) ans =2*cos(x2)*x diff(t6,1), diff(t6, 2), diff(t6, 3),diff(t6, 4),diff(t6, 5),diff(t6, 6) ans = 6*t5, 30*t4, 120*t3, 360*t2, 720*t, 720, syms x x1 alpha u t; int(1/(1+x2) ans =atan(x) int(sin(alpha*u),alpha) ans =-1/u*cos(alpha*u) int(x1*log(1+x1),0,1) ans =1/4,4.4 数学问题的解析运算,极限、微分、积分和 Taylor级数展开 limit, diff, int, taylor, taylor(exp(-x) ans =1-x+1/2*x2-1/6*x3+1/24*x4-1/120*x5 taylor(log(x),6,1) ans =x-1-1/2*(x-1)2+1/3*(x-1)3-1/4*(x-1)4+1/5*(x-1)5 taylor(sin(x),pi/2,6) ans =1-1/2*(x-1/2*pi)2+1/24*(x-1/2*pi)4 taylor(xt,3,t) ans =1+log(x)*t+1/2*log(x)2*t2,4.4 数学问题的解析运算,无穷级数求和 symsum, syms n k; sym