matlab第四章数值计算(编辑修改稿)

matlab第四章数值计算(编辑修改稿)内容摘要：

（ 4） xx=0::。 ff=exp(xx.^2)。 pp=spline(xx,ff)。 int_pp=fnint(pp)。 Ssp=ppval(int_pp,[0,1])*[1。 1] Ssp = （ 5） 图 开型数值积分 [] function g = gauss10(fun,a,b) %GAUSS10(fun,a,b) % fun %====================================================== x = [。 ...。 ]。 w = [。 ...。 ]。 t = .5*(b+a)+.5*(ba)*[flipud(x)。 x]。 W = [flipud(w)。 w]。 g = sum(W.*feval(fun,t))*(ba)/2。 【例 】当 )cos()( xxf  时，比较解析积分和近似积分。 （ 1） syms x。 F=int(39。 cos(x)39。 ,39。 x39。 ,1,1),vpa(F) F = 2*sin(1) ans = （ 2） aF=cos(1/sqrt(3))+cos(1/sqrt(3)) aF = 【例 】求  101ln dxxI ，准确结果是 .886226922 。 （ 1） 15 syms x。 IS=vpa(int(39。 sqrt(log(1/x))39。 ,39。 x39。 ,0,1)) Warning: Explicit integral could not be found. In D:\MATLAB6P5\toolbox\symbolic\@sym\ at line 58 In D:\MATLAB6P5\toolbox\symbolic\@char\ at line 9 IS = .88622692545275801364908374167057 （ 2）用 quad指令求积 ff=inline(39。 sqrt(log(1./x))39。 ,39。 x39。 )。 Isim=quad(ff,0,1) Warning: Divide by zero. In D:\MATLAB6P5\toolbox\matlab\funfun\ at line 13 In D:\MATLAB6P5\toolbox\matlab\funfun\@inline\ at line 34 In D:\MATLAB6P5\toolbox\matlab\funfun\ at line 62 Isim = （ 3） Ig=gauss10(ff,0,1) Ig = 多重数值积 分 积分限为常数的二重积分指令 【例 】计算    211001 dydxxS yx和    102112 dydxxS yx。 （ 1） syms x y ssx01=vpa(int(int(x^y,x,0,1),y,1,2)) ssx12=vpa(int(int(x^y,y,0,1),x,1,2)) Warning: Explicit integral could not be found. In D:\MATLAB6P5\toolbox\symbolic\@sym\ at line 58 ssx01 = .40546510810816438197801311546435 ssx12 = （ 2） zz=inline(39。 x.^y39。 ,39。 x39。 ,39。 y39。 )。 nsx01=dblquad(zz,0,1,1,2) nsx12=dblquad(zz,1,2,0,1) nsx01 = nsx12 = 内积分限为函数的二重积分 [] function SS=double_int(fun,innlow,innhi,outlow,outhi) %double_int 16 % %fun %innlow,innhi %outlow,outhi y1=outlow。 y2=outhi。 x1=innlow。 x2=innhi。 f_p=fun。 SS=quad(@G_yi,y1,y2,[],[],x1,x2,f_p)。 [] function f=G_yi(y,x1,x2,f_p) %G_yi %y %x1,x2 % %f_p y=y(:)。 n=length(y)。 if isnumeric(x1)==1。 xx1=x1*ones(size(y))。 else xx1=feval(x1,y)。 end if isnumeric(x2)==1。 xx2=x2*ones(size(y))。 else xx2=feval(x2,y)。 end for i=1:n f(i)=quad(f_p,xx1(i),xx2(i),[],[],y(i))。 end f=f(:)。 【例 】计算 dydxyxIy   412 22 )(。 （ 1） [] function f=x_low(y) f=sqrt(y)。 （ 2） （ 3） ff=inline(39。 x.^2+y.^239。 ,39。 x39。 ,39。 y39。 )。 SS=double_int(ff,@x_low,2,1,4) Warning: Minimum step size reached。 singularity possible. In D:\MATLAB6P5\toolbox\matlab\funfun\ at line 88 In D:\MATLAB6p5\work\ at line 11 In D:\MATLAB6P5\toolbox\matlab\funfun\@inline\ at line 20 In D:\MATLAB6P5\toolbox\matlab\funfun\ at line 62 In D:\MATLAB6p5\work\ at line 8 SS = （ 4） Ssym=vpa(int(int(39。 x^2+y^239。 ,39。 x39。 ,39。 sqrt(y)39。 ,2),39。 y39。 ,1,4)) Ssym = 17 卷积 “完整”离散序列的数值卷积 “截尾”离散序列的数值卷积 多项式乘法与离散卷积的算法同构 连续时间函数的数值卷积 卷积的 MATLAB 实现 【例 】有序列  e lsennA 12,4,301)(  和   e lsennB 9,3,201)( 。 （ A） % a=ones(1,10)。 n1=3。 n2=12。 b=ones(1,8)。 n3=2。 n4=9。 c=conv(a,b)。 nc1=n1+n3。 nc2=n2+n4。 kc=nc1:nc2。 % aa=a(1:6)。 nn1=3。 nn2=8。 cc=conv(aa,b)。 ncc1=nn1+n3。 nx=nn2+n4。 ncc2=min(nn1+n4,nn2+n3)。 kx=(ncc2+1):nx。 kcc=ncc1:ncc2。 N=length(kcc)。 stem(kcc,cc(1:N),39。 r39。 ,39。 filled39。 ) axis([nc12,nc2+2,0,10]),grid,hold on stem(kc,c,39。 b39。 ),stem(kx,cc(N+1:end),39。 g39。 ),hold off 4 6 8 10 12 14 16 18 20 22012345678910 图 【例 】求函数 )()( tUetu t 和 )()( 2/ tUteth t 的卷积。 本例展示：（ A）符号 Laplace 变换求卷积的理论表示；（ B） SIMULINK 卷积法的执行过程和它的快速精确性。 （ C）从理论符号解产生相应的理论数值序列。 （ 1） syms tao。 t=sym(39。 t39。 ,39。 positive39。 )。 US1=laplace(exp(t))。 HS1=laplace(t*exp(t/2)) yt1=simple(ilaplace(US1*HS1)) HS1 = 18 1/(1/2+s)^2 yt1 = 4*exp(t)+(2*t4)*exp(1/2*t) （ 2） 图 （ 3） t=yt2(:,1)。 yyt1=eval(vectorize(char(yt1)))。 [dy,kd]=max(abs(yyt1yt2(:,2)))。 dy12=dy/abs(yyt1(kd)) dy12 = 【例 】用“零阶” 近似法求 )()( tUetu t 和 )()( 2/ tUteth t 的卷积。 本例演示：（ A）连续函数的有限长度采样。 （ B）卷积数值计算三个误差（“截尾”误差、“零阶”近似误差、计算机字长误差）的影响。 （ C）卷积“无截尾误差”区间、“非平凡”区间端点的确定。 （ D）离散卷积和连续卷积之间的关系。 （ E）指令 conv 的使用。 （ F）绘图分格线的运用。 （ 1） （ 2） % t2=3。 t4=11。 T=。 tu=0:T:t2。 N2=length(tu)。 th=0:T:t4。 N4=length(th)。 u=exp(tu)。 h=th.*exp(th/2)。 tx=0:T:(t2+t4)。 Nx=length(tx)。 yt3=T*conv(u,h)。 % t=tx。 yyt1=eval(vectorize(char(yt1)))。 [dy,kd]=max(abs(yyt1(1:N2)yt3(1:N2)))。 dy13(1)=dy/abs(yyt1(kd))。 [dy,kd]=max(abs(yyt1(N2+1:N4)yt3(N2+1:N4)))。 dy13(2)=dy/abs(yyt1(N2+kd))。 [dy,kd]=max(abs(yyt1(N4+1:Nx)yt3(N4+1:Nx)))。 dy13(3)=dy/abs(yyt1(N4+kd))。 （ 3） disp(39。 与理论结果的相对误差 39。 ) disp([blanks(4),39。 [0,3]段 [3,11]段 [11,14]段 39。 ]),disp(dy13) plot(t,yyt1,39。 :b39。 ,tx,yt3,39。 r39。 ) set(gca,39。 Xtick39。 ,[0,3,11,14]),grid 与理论结果的相对误差 [0,3]段 [3,11]段 [11,14]段 19 0 3 11 140 图 随机数据的统计描述 统计分布的 数字特征 【例 】样本统计特征量计算示例。 % X(:,1)=ones(10,1)。 X(1,1)=100。 X(10,1)=。 rand(39。 state39。 ,1)。 randn(39。 state39。 ,1) X(:,2)=rand(10,1)。 X(:,3)=randn(10,1)。 X(:,3)=2*abs(min(X(:,3)))+X(:,3)。 % =mean(X)。 =median(X)。 =geomean(X)。 =harmmean(X)。 % =std(X)。 =var(X)。 =mad(X)。 =range(X)。 % X,Moment1,Moment2 X = Moment1 = arithmetic: [ ] median: [1 ] geometric: [1 ] harmmonic: [ ] 20 Moment2 = Standard: [ ] variance: [ ] absolute: [ ] range: [ ] 样本分布的频数直方图描述 【例 】 hist 指令的使用示例。 randn(39。 state39。 ,1),rand(39。 state39。 ,31) x=randn(100,1)。 y=rand(100,1)。 % subplot(1,2,1),hist(x,7) subplot(1,2,2),histfit(x,20) 5 0 50510152025305。