作者Gwaewluin (神无月 孝臣)
看板MATLAB
标题Re: [请益] 如何写辛浦森积分
时间Tue Nov 21 10:48:13 2006
※ 引述《LESPAUL (大帅)》之铭言:
: 请教各位高手
: 如何不使用内建的积分quad积分指令来写一个辛浦森积分
: 我想写 W在[a,b]区间之内,z及Wn都为一个固定常数下
: 对(z.*(W/Wn).^2)/((1-(W/Wn).^2).^2+(2*z.*(W/Wn)).^2)做积分
: 如何使用辛浦森积分公式
: b
: S s(x) dx= h/3*{s + 4(s +s +s +...+s ) + 2(s +s +...s )+s }
: a 0 1 3 5 2n-1 2 4 2n-2 2n
: 其中a,b分别表示积分的下限与上限
: n指ab间成对的区间数目
: h=(b-a)/(2n) 代表每区间的半宽度
: 请教各位高手要怎麽写辛浦森的积分
: 感激不尽啊
新浦森法有一个限制
就是他拿来积分的点个数一定要是奇数个(观察一下你方程式的index就可以知道)
也就是使用时取的区间必须为偶数个
所以只要注意一下index的特色就很好写了
s = 5 % 任意值,你想切几等分一半的数字就是
N = 2 * s ; % 只是确定区间为偶数
xn = linspace( a , b , N + 1 ) ;
yn = f( xn ) ; % 你要的方程式给xn输出成yn,记得yn要输出成矩阵
Is = (xn(2)-xn(1))*(yn(1)+yn(end)+4*sum(yn(2:2:end))+2*sum(yn(3:2:end-1)))/3 ;
这样积分就出来了
这只是最简单的给几个区间就直接计算出答案的写法
如果要进行收敛测试的话
区间数量的变化必须为2 4 8 16 32 64 ......
这边也只需要小调整即可
不过要做到收敛测试的话辛浦森法其实比较没啥作用了
等於只是Romberg method的第二层而已
而第二层却可以使用梯形法计算的结果处理出来
也就是说
想要玩收敛测试的话只需要梯形法配合Romberg method
根本不必再写一个新浦森法
--
Deserves death! I daresay he does. Many that live deserve death. And some die
that deserve life. Can you give that to them? Then be not too eager to deal out
death in the name of justice, fearing for your own safty. Even the wise cannot
see all ends.
Gandalf to Frodo
--
※ 发信站: 批踢踢实业坊(ptt.cc)
◆ From: 140.120.25.235
※ 编辑: Gwaewluin 来自: 140.120.25.235 (11/21 10:53)