第五章 MATLAB数值积分.pptVIP

定积分的MATLAB符号计算： f=int(s,a,b) 输入量：s可以是被积函数，也可以是几个被积 函数组成的矩阵，积分上下限可以是数也可以是符 号变量，积分变量为x。 f=int(s,x,a,b) 例:计算 解：syms t f=5/((t-1)*(t-2)*(t-3));F=int(f,t,4,5),y=double(F) 例： 解： syms x a b c d n f=(a*x+b)/((c*x+d)^(1/n)); F=int(f,0,1),y=simple(F) Matlab里提供了复合梯形公式的计算程序— trapz.m，其调用格式为： Z=trapz(y） y为向量，则z是y的积分，若y为矩阵，则行向量 z是y的每一列的积分； Z=trapz(x,y) 输入x、y为同长度的数组，计算y对x积分。 调用trapez_g函数：第一步建立函数文件： function y=f(x) y=sqrt(1 + exp(x)) 再调用： I = trapez_g(@f, 0, 2, 16) X坐标下的节点及被积函数： 注：函数quad (‘f_name’,a,b,tol,trace)用于求被积函 数f(x)在[a,b]上的定积分，tol是计算精度，调用 quad函数时，先要建立一个描述被积函数f(x)的函 数文件或语句函数。当被积函数f_name含有一个以 上的变量时，quad函数的调用格式为： quad (‘f_name’, a,b,tol,trace,g1,g2) 其中f_name,a,b,tol,trace等参数的含义同前。 例 用两种不同的方法求积分。 先建立一个函数文件ex.m： function ex=ex(x) ex=exp(-x.^2); %注意应用点运算 return 然后，在MATLAB命令窗口，输入命令： quad(ex,0,1,1e-6) %注意函数名应加字符引号 quad8(ex,0,1,1e-6) %用另一函数求积分 实际上，MATLAB提供了计算二重积分的函数： dblquad(f,a,b,c,d,tol,trace) 该函数求f(x,y)在[a,b]×[c,d]区域上的二重积分。参 数tol，trace的用法与函数quad完全相同。 命令如下： g=inline(exp(-x.^2-y.^2)); dblquad(g,0,1,0,1) %直接调用二重积分函数求解 其幂级数形式为： n 阶Legendre多项式的 n 个根全部位于-1与1之间，n＝7时的根如下，在端点处没有根，因此向两端点方向区间长度逐渐减少。 （系数可由Legen_pw(n)确定） 对积分 用 n 个Legendre节点的插值多项式拟合被积函数，然后积分： 其中， 是 n 阶Legendre多项式的根. 当积分限不是-1和1时，通过坐标变换可将这些根映 射到积分上、下限 a，b之间，如： Legendre求积公式： Legendre多项式的根可由roots(Legen_pw(n))给出， 该求积公式由函数Gauss_q实现。 function I = Gauss_q(f_name, a, b, n) p=legen_pw(n);x = roots(p);x = sort(x); for j=1:n y = zeros(1,n); y(j)=1; p = polyfit(x,y,n-1); P = poly_itg(p); w(j) = polyval(P,1) - polyval(P,-1); end x = 0.5*((b-a)*x + a + b);y=feval(f_name, x); I = sum(w.*y)*(b-a)/2; fprintf(\n x y w \n) for j=1:n，fprintf(%e %e %e\n, x(j),y(j), w(j))；end function py = poly_itg(p) n=length(p); py = [p.*[n:-1:1].^(-1),0]; 例4：用6个节点的Gauss公式计算积分 I=Gauss_q(fex,0,pi,6) function f=fex(x) f=1./(2+cos(x)); 结果： I＝1.8138 五、积分限无界及被积函数有奇点的数值积分法 考虑如下积分 无穷区间的积分 无界函数的积分 一个函数在无穷区间或半无穷区间可积，当且