请教高手解决下面这题
用(梯形)龙贝格求积公式近似计算I=(积分上限1,下限0)2e(-x)/根号(x)dx 注:由于这里输不了公式,我用文字代替了,使其精度达到1e-5 我运行时出错,找不到原因,请各位高手帮助我,非常感激!
程序代码:
function [val,M]=Romberg(f,ab,dalt) % Romberg算法计 算积分 % [val,M]=Romberg(f,ab,dalt) % val 返回积分值 % M T数表 % f 被积函数,函数文件名 % ab 积分区间 % dalt 精度 if nargin<3 dalt=1e-7; %设置默认精度1e-7 end i=1; h=ab(2)-ab(1); t=0; j=0; T(1,1)=(h/2)*(feval(f,ab)*ones(2,1)); while i<50 %最大跌?0次 i=i+1; h=h/2; T=[T zeros(i-1,1);zeros(1,i)]; x=ab(1)+h:2*h:ab(2)-h; y=feval(f,x)*ones(size(x')); % feval(f,x)求f在x点值 T(i,1)=T(i-1,1)/2+h*y; % 细化区间,求得梯形值 for t=2:i j=t-1; T(i,t)=(4^j*T(i,j)-T(i-1,j))/(4^j-1); %外推加速 end if abs(T(i,i)-T(i-1,i-1))<=dalt %控制精度 break end end if nargout==2 M=T; %按需求返回T数表 end val=T(i,i); %积分值 脚本文件: y=inline('2*exp(-x).*(x.^(-1/2))'); [I,M]=Romberg(y,[0 1],1.0e-5)