震荡的瑕积分的数值算法wangtianyi.docVIP

武汉大学数学与统计学院 数值分析实验报告 实验名称 对于震荡的瑕积分的数值算法的讨论 实验时间 2006年 12 27 日 姓名 王天怡 班级 数类一班 学号 200531000227 成绩 一、实验目的，内容 二、相关背景知识介绍 三、代码 四、数值结果 五、计算结果的分析 六、计算中出现的问题，解决方法及体会 一、实验目的，内容 通过对一个具体的震荡的函数在瑕点的积分（），讨论对于震荡的瑕积分的数值算法。 二、相关背景知识介绍 在瑕积分中，如果自变量在趋向瑕点时原函数值趋向于某个正常极限，在求积分时可以用补点的方法构造一个按段光滑的函数（可积），再对构造出的按段光滑的函数积分即可。如果自变量在趋向瑕点时原函数值趋向于某个反常极限，在求积分时可以用局部的方法构造改变积分区间使之变成一个正常的数值积分（对上一种情况也可以使用类似的处理手法）。但对于当自变量在趋向瑕点时原函数值震荡的瑕积分，则只对积分进行一些变化。其中比较有效的是：级数法，换元法和复积分。 Ⅰ级数法： 则： 然后可以近似： Ⅱ换元法 其余的部分和方法Ⅰ类似 注：方法Ⅰ的收敛速度比方法Ⅱ快，但方法Ⅱ的原函数振幅为1可以进行误差分析。 Ⅲ复积分 参考文献：Cleve B.Moler ,Numerical Computing with MATLAB F.Bornemann,D.Laurie,S.Wagon,and ,J.Waldogel,The SIAM 100-DigitChallenge 三、代码（Matlab） Ⅰ级数法 1.在积分的局部使用Simposon公式积分的代码： x1=2; x2=1; s=0; k=0; sq=0; while k1000 psq=sq; sq=0; x1=x2; k=k+1; w=(k-1/2)*pi; x2=lambertw(w)/w; c=x2; h=(x1-x2)/1000; for i=1:1000 sq=sq+(fun(c)+4*fun(c+h)+fun(c+2*h))*h/6; c=c+2*h; end s=s+sq; end st=s-sq^2/(sq-psq) 2.在积分中fun的代码 function y=fun(x); y=(cos(log(x)/x))/x; end 3.在积分的局部使用自适应积分的代码： x1=2; x2=1; s=0; k=0; sq=0; while k500 psq=sq; sq=0; x1=x2; k=k+1; w=(k-1/2)*pi; x2=lambertw(w)/w; sq=quadl(fun1,x2,x1,1.e-12); s=s+sq; end st=s-sq^2/(sq-psq) 4.在积分中fun1的代码 function y=fun1(x); y=(cos(log(x)./x))./x; end Ⅱ换元法 主程序： x1=2; x2=0; s=0; k=0; sq=0; pst=0; while k10000 psq=sq; sq=0; x1=x2; k=k+1; w=(k-1/2)*pi; x2=-lambertw(w); sq=quadl(fu,x2,x1,1.e-12); s=s+sq; if k==1 st=s; else st=s-sq^2/(sq-psq); end pst=st; if (x1-x2)1.e-11 k break; end end st 在积分中fu的代码： function y=fu(x) y=cos(x.*exp(-x)); end Ⅲ复积分 function y=func(t) z=t+i.*t.*(1-t); y=real(z.^(i./z-1).*(1+i*(1-2*t))); end 在求值时使用命令： quad(func,0,1,1.e-14) 四、数值结果 Ⅰ级数法 在积分的局部使用Simposon公式积分方法的结果： st = 0.41280001534037 2. 在积分的局部使用自适应积分方法的结果： st = 0.32336743114790 Ⅱ换元法 st = 0.323367431677