プログラミング論 数値積分 http://www.ns.kogakuin.ac.jp/~ct13140/Prog/

Presentation on theme: "プログラミング論 数値積分 http://www.ns.kogakuin.ac.jp/~ct13140/Prog/"— Presentation transcript:

1 プログラミング論 数値積分

2 区分求積法と台形公式のみ解説 多くの変数が登場する．しっかり整理すること
積分 区分求積法と台形公式のみ解説 多くの変数が登場する．しっかり整理すること

3 積分 不定積分は，「微分」の逆操作．すなわち，f(x)を不定積分する/の不定積分を求めるとは，「微分したらf(x)になる関数を探す」こと．
定積分は，関数f(x)とx軸に挟まれる部分の面積を求めること． 本講義では，定積分のみを扱う．

4 定積分 とは？ y=f(x)と，x軸と x=aと，x=bで 囲まれた領域の面積を 求めること x=a x=b

5 この短冊の面積は，f(xi)(xi+1－xi)
区分求積法 細かい区間毎の長方形の面積の合計で近似． 全短冊の面積の合計と， はほぼ等しい f(xi+1) f(xi) y=f(x) この短冊の面積は，f(xi)(xi+1－xi) a xi xi+1 b

6 区分求積法 a=x0<x1<…<xn-1<xn=bとして， f (x0)，f (x1)，…，f (xn)が既知なら

7 区分求積法 当然， 長方形の高さは 左の値を使用しても， 右の値を使用しても良い． f(xi) f(xi+1) 左の値を使用 xi xi+1
右側の値を使用 xi xi+1

8 区分求積法の誤差 f(xi+1) f(xi) この部分が近似の誤差 xi xi+1

9 区分求積法 より細かい区間に分割するほど誤差は少なくなる． 長方形でなく台形で近似した方が，より正確な近似になると期待される．
台形公式 (これは1次近似である) 近傍3点を用い，2次式で近似した方がより正確な近似になると期待される． シンプソンの公式 3次式，4次式で近似した方が更に正確． ニュートン・コーツの積分公式

10 区分求積法台形公式 長方形ではなく，台形で近似する． より少ない誤差で計算が可能 f(xi) f(xi+1) xi xi+1

11 区分求積法台形公式 a=x0<x1<…<xn-1<xn=bとして，
f (x0)，f (x1)，…，f (xn)が既知なら 区間[a,b]をn等分した場合， として

12 C言語プログラミング を求めるには？ h = (b-a)/N; N等分すると仮定 f(x1) f(x2) h a x1 x2 x3 x4 b
xn

13 C言語プログラミング を求めるには？ x0～xnを表示してみる． for(i=0; i<=N; i++){
h = (b-a)/N; x0 = a; x1 = a + h*1; x2 = a + h*2; : xi = a + h*i; なので， for(i=0; i<=N; i++){ printf("%lf\n", a+h*i); } N等分すると仮定 x1 x2 f(x1) f(x2) a b h x0 x3 x4 xn

14 C言語プログラミング を求めるには？ 各短冊の面積を表示してみる． for(i=0; i<N; i++){ xi = a+h*i;
xi～xi+1の短冊の横の長さは h， 縦の長さは f (xi) なので， h = (b-a)/N; for(i=0; i<N; i++){ 　xi = a+h*i; area = f(xi)*h; printf("%lf\n",area); } x1 x2 f(x1) f(x2) a b h x0 x3 x4 xn 関数 double f(double d)が 存在していると仮定.

15 C言語プログラミング を求めるには？ 全短冊の面積の合計を出す． すなわち，定積分． for(i=0; i<N; i++){
h = (b-a)/N; sum = 0.0; for(i=0; i<N; i++){ 　xi = a+h*i; area = f(xi)*h; sum += area; } printf("%lf\n", sum); x1 x2 f(x1) f(x2) a b h x0 x3 x4 xn

16 例(伊) y=x2 と x軸と x=1 と x=3 で囲まれた領域の面積は? y=x2 2 4 6 8 10 12 0.5 1 1.5
2 4 6 8 10 12 0.5 1 1.5 2.5 3 3.5 x y y=x2

17 例(伊) 100分割，非台形 #define NUM 100 #include <stdio.h> void main(){
int i; double x, y, sum, h; h = 2.0/NUM; sum = 0.0; for(i=0; i<NUM; i++){ x = h*i; y = x*x; sum += y*h; } printf("%lf\n", sum); 例(伊) 100分割，非台形

18 例(伊) 1000分割，台形 #define NUM 1000 #include <stdio.h> void main(){
int i; double xa, xb, ya, yb, sum, h; h = 2.0/NUM; sum = 0.0; for(i=0; i<NUM; i++){ xa = h*i; xb = h*(i+1); ya = xa*xa; yb = xb*xb; sum += (ya+yb)*h/2.0; } printf("%lf\n", sum); 例(伊) 1000分割，台形

19 例(伊) 100分割 1000分割 面積 誤差 左値を 使用 8.5868 右値を 8.7468 台形 8.6668

20 例(呂) 下の扇形(1/4円)の面積は? 1.0 0.8 0.6 0.4 0.2 0.2 0.4 0.6 0.8 1.0

21 例(呂) 1000分割，非台形 面積=0.785889 PI=3.143555 #define NUM 1000
#include <stdio.h> #include <math.h> double f(double x){ return sqrt(1-x*x); } void main(){ int i; double x, y, sum, h; h = 1.0/NUM; sum = 0.0; for(i=0; i<NUM; i++){ x = h*i; y = f(x); sum += y*h; printf("面積=%lf\nPI=%lf\n", sum, sum*4); 例(呂) 1000分割，非台形 面積= PI=