数値積分

台形公式とシンプソンの公式による数値積分を計算します。

台形公式

曲線の下方の台形の面積を求めて合計します。

関数 f(x) は定義されているものとし，引数を範囲と分割数にします。戻り値が台形の面積の合計です。

trapez(range,divn):=(

regional(sx,ex,dx,xn,yn,sum);

sx=range_1; //定義域左端

ex=range_2; //定義域右端

dx=(ex-sx)/divn; //刻み幅

xn=apply(0..divn,sx+#*dx); // sxからexまでをn等分した分点のリスト

yn=apply(xn,f(#));

sum=0;

repeat(divn,

drawpoly([[xn_#,0],[xn_#,yn_#],[xn_(#+1),yn_(#+1)],[xn_(#+1),0]],color->[0,0.8,0],alpha->0.7);

connect([[xn_#,0],[xn_#,yn_#],[xn_(#+1),yn_(#+1)],[xn_(#+1),0]]);

sum=sum+dx*(yn_#+yn_(#+1))/2;

);

);

シンプソンの公式

曲線の下方の図形を線分ではなく放物線で近似します。

sympson(range,divn):=(

regional(sx,ex,dx,xn,yn,sum,x0,x1,x2,y0,y1,y2);

sx=range_1; //定義域左端

ex=range_2; //定義域右端

dx=(ex-sx)/divn; //区間幅

xn=apply(0..n,sx+#*dx); // sxからexまでをn等分した分点のリスト

yn=apply(xn,f(#));

p(x,x0,y0,x1,y1,x2,y2):=(x-x1)*(x-x2)*y0/((x0-x1)*(x0-x2))+(x-x0)*(x-x2)*y1/((x1-x0)*(x1-x2))+(x-x0)*(x-x1)*y2/((x2-x0)*(x2-x1)); //3点を通る放物線

sum=0;

dx=dx/2;

repeat(n,i,

x0=xn_i;

x1=(xn_i+xn_(i+1))/2;

x2=xn_(i+1);

y0=yn_i;

y1=f(x1);

y2=yn_(i+1);

draw([x0,0],[x0,y0],color->[0,0.8,0],alpha->0.7);

draw([x1,0],[x1,y1],color->[0,0.8,0],alpha->0.7);

draw([x2,0],[x2,y2],color->[0,0.8,0],alpha->0.7);

plot(p(#,x0,y0,x1,y1,x2,y2),start->x0,stop->x2,color->[0,0.7,0]);

sum=sum+dx*(y0+4*y1+y2)/3;

);

sum;

);

たとえば，次のようにして使います。

f(x):=1/(x^2+1);

n=20;

plot(f(#),color->[0,0,0]);

ret=4*trapez([0,1],n); または ret=4*sympson([0,1],n);

drawtext([5,2],"分割の数="+n,size->16);

drawtext([5,1],"近似した面積の和="+ret,size->16);

< 戻る ＞