Numerikus integrálás
Contents
Trapézformula
A közelítő integrált kiszámító függvény:
function kvadr=trapez(a,b,n,fv)
h=(b-a)/n;
x=(a+h):h:(b-h);
kvadr=h*sum(fv(x))+h*fv(a)/2+h*fv(b)/2;
A Matlab beépített intergálófüggvénye a integral(), ezzel ellenőrizhetjük a konvergencia sebességét. Legyen a függvény f(x)=x*sin(x), ezt egy külön fájlba mentjük fv.m néven. Ekkor:
for i=1:4 trapez(1,2,10^i,@fv)-integral(@fv,1,2) end
ans = -0.0011 ans = -1.0873e-05 ans = -1.0873e-07 ans = -1.0873e-09
A hiba jól láthatóan O(h^2), hiszen egy lépsében 1/10-edére vettem az osztásközt, ettől a hiba kb 1/100-adára változott.
Simpson formula
Órai feladat: Valósítsuk meg Matlabban a Simpson-kvadratúrát!
function kvadr=simpsonkv(a,b,n,fv)
h=(b-a)/n;
x=a+h : h: b-h;
x2=a+h/2 : h : b-h/2;
kvadr=h*(2/6*sum(fv(x))+4/6*sum(fv(x2))+fv(a)/6+fv(b)/6);
Innen:
for i=1:4 simpsonkv(1,2,10^i,@fv)-integral(@fv,1,2) end
ans =
4.0605e-08
ans =
4.0594e-12
ans =
0
ans =
0
A hiba O(h^4). Persze létezik jobb, ha a 3-adfokú (4-edfokú,...) Lagrange-interpolációs polinommal közelítünk.