数值分析 辛普森公式论文

数值分析辛普森公式论文 辛普森公式是牛顿-科特斯公式当n=2时的情形，也称为三点公式。该公式利用区间二等分的三个点来进行积分插值。其科特斯系数分别为1/6，4/6，1/6。 一、辛普森公式简介 辛普森公式是一种数值分析方法，用于计算立体图形的体积。设拟柱体的高（两底面α、β间的距离）为H，如果用平行于底面的平面γ去截该图形，所得到的截面面积是平面γ与平面α之间距离h的不超过3次的函数，那么该拟柱体的体积V为： V = H (S_1 + 4S_0 + S_2) /6 式中，S_1和S_2是两底面的面积，S_0是中截面的面积（即平面γ与平面α之间距离h=H/2时得到的截面的面积）。 实际上，不光是拟柱体，其他符合条件（所有顶点都在两个平行平面上、用平行于底面的平面去截该图形时所得到的截面面积是该平面与一底之间距离的不超过3次的函数）的立体图形也可以利用该公式求体积。 二、计算实例 例1：计算底面积为S、高为h的柱体的体积。 解：此题中S_1 = S_0 = S_2 = S，H = h，所以V = H (S_1 + 4S_0 + S_2) /6 =h (S + 4S + S) /6 = S h。 例2：计算底面积为S、高为h的锥体的体积。 解：此题中S_1 = S，S_0 = S /4，S_2 = 0，H = h，所以V = H (S_1 + 4S_0 + S_2) /6 = h (S + 4S /4 + 0) /6 = S h /3。 例3：计算半径为r的球的体积。 解：此题中S_1 = S_2 = 0，S_0 = πr^2，H = 2r，所以V = H (S_1 + 4S_0 + S_2) /6 = 2r (0 + 4πr^2 + 0) /6 = 4πr^3 /3。 三、公式证明 只需要证明根据公式算出来的体积和用积分算出来的体积相等即可。 设截面面积是截面高h的不超过3次的函数：f(h)= ah^3 + bh^2 + ch + d。 那么，利用积分计算体积，可以得到： V = ∫ f(x) dx = ah^4 /4 + bh^3 /3 + ch^2 /2 +dh 利用公式计算体积，可以得到： V = H (S_1 + 4S_0 + S_2) /6 = h ( f(0) + 4f(h/2) + f(h) ) /6 = h [ d + 4 (ah^3 /8 + bh^2 /4 + ch /2 + d) + (ah^3 + bh^2 + ch + d) ]/6 = ah^4 /4 + bh^3 /3 + ch^2 /2 +dh 因此两式相等，公式得证。 Remark：当函数f(h)次数高于或等于4次时，公式一般不成立。这只需验证f(h)=h^4时公式不成立即可。 二、关于辛普森公式的编程 MATLAB用辛普森公式求重积分 function q=DblSimpson(f,a,A,b,B,m,n) if(m==1 && n==1) %辛普森公式 q=((B-b)*(A-a)/9)*(subs(sym(f),findsym(sym(f)),{a,b})+... subs(sym(f),findsym(sym(f)),{a,B})+... subs(sym(f),findsym(sym(f)),{A,b})+... subs(sym(f),findsym(sym(f)),{A,B})+... 4*subs(sym(f),findsym(sym(f)),{(A-a)/2,b})+... 4*subs(sym(f),findsym(sym(f)),{(A-a)/2,B})+... 4*subs(sym(f),findsym(sym(f)),{a,(B-b)/2})+... 4*subs(sym(f),findsym(sym(f)),{A,(B-b)/2})+... 16*subs(sym(f),findsym(sym(f)),{(A-a)/2,(B-b)/2})); else %复合辛普森公式 q=0; for i=0:n-1 for j=0:m-1 x=a+2*i*(A-a)/2/n; y=b+2*j*(B-b)/2/m; ... ... ... q=q+... end end end