数値解析シラバス Ｃ言語環境と数値解析の概要（１ 回） Ｃ言語環境と数値解析の概要（１ 回） 方程式の根（１回） 方程式の根（１回） 連立一次方程式（２回） 連立一次方程式（２回） 補間と近似（２回） 補間と近似（２回） 数値積分（１回） 数値積分（１回） 中間試験（１回） 中間試験（１回） 常微分方程式（１回）

Presentation on theme: "数値解析シラバス Ｃ言語環境と数値解析の概要（１ 回） Ｃ言語環境と数値解析の概要（１ 回） 方程式の根（１回） 方程式の根（１回） 連立一次方程式（２回） 連立一次方程式（２回） 補間と近似（２回） 補間と近似（２回） 数値積分（１回） 数値積分（１回） 中間試験（１回） 中間試験（１回） 常微分方程式（１回）"— Presentation transcript:

1 数値解析シラバス Ｃ言語環境と数値解析の概要（１ 回） Ｃ言語環境と数値解析の概要（１ 回） 方程式の根（１回） 方程式の根（１回） 連立一次方程式（２回） 連立一次方程式（２回） 補間と近似（２回） 補間と近似（２回） 数値積分（１回） 数値積分（１回） 中間試験（１回） 中間試験（１回） 常微分方程式（１回） 常微分方程式（１回） 偏微分方程式（２回） 偏微分方程式（２回） 数値計算の応用（３回） 数値計算の応用（３回） 期末試験（１回） 期末試験（１回） ６月１２日（水）５限 今週 来週

2 数値積分 解析的に積分が困難な関数を数値的に積分する

3 数値積分 台形公式 台形公式 関数を直線で近似し、台形の集合体とす る。 関数を直線で近似し、台形の集合体とす る。 x y x0x0 x1x1 x2x2 x3x3 y0y0 y1y1 y2y2 y3y3

4 台形 I = 0.5×h {(y 0 +y 1 )+(y 1 +y 2 )+.. + (y n-1 +y n )} I = 0.5×h {(y 0 +y 1 )+(y 1 +y 2 )+.. + (y n-1 +y n )} = 0.5×h(y 0 +2y 1 + 2y 2 +.. + y n ) = 0.5×h(y 0 +2y 1 + 2y 2 +.. + y n ) x y x0x0 x1x1 x2x2 x3x3 y0y0 y1y1 y2y2 y3y3 h

5 台形公式計算例 x y 012 25 50 75 100 I = x 4 +2x dx = I = x 4 +2x dx = 3 0 57.6 3 I = 0.5×1×(0+2 ･ 3+2 ･ 20+87)=66.5 I = 0.5×1×(0+2 ･ 3+2 ･ 20+87)=66.5

6 Ｅｘｃｅｌで確認 ABC 1XY 20.0=A2^4+2*A2=B2 30.5=2*B3 41.0 51.5 62.0 72.5 83.0=B8 =(0.5/2)*sum(C2:C8) 区間の幅 数値の合計 I = x 4 +2x dx = I = x 4 +2x dx = 3 0 57.6 N = 6

7 Ｅｘｃｅｌで確認 ABC 1XY 20.0=A2^4+2*A2=B2 30.1=2*B3 40.2 302.8 312.9 323.0=B8 =(0.1/2)*sum(C2:C32) 区間の幅 数値の合計 I = x 4 +2x dx = I = x 4 +2x dx = 3 0 57.6 N = 30

8 台形公式による積分 #include #defineN30 doublefunc_y( double ); int main( ) { doubley[N+1]; doublexa = 0.0, xb = 3.0; doublez = 0.0, h, x, s; inti;

9 h = (xb - xa)/(double)N; for(i=0;i<=N;i++){ x = xa + h*(double)i; y[i] = func_y(x); } for(i=1;i<N;i++)z += 2.0*y[i]; s = (h/2.0)*(y[0]+z+y[N]); printf("ANS = %8.4lf\n",s); return 0; } doublefunc_y(double x){ returnpow(x, 4.0) + 2.0*x; } I =0.5×h(y 0 +2y 1 + 2y 2 +.. + y n ) I =0.5×h(y 0 +2y 1 + 2y 2 +.. + y n ) y = x 4 + 2x y = x 4 + 2x

10 数値積分での誤差縮小 分割数を多くする． 分割数を多くする． 曲線で近似 → シンプソン 曲線で近似 → シンプソン

11 厳密解の推定 nanswer 462.6309 858.8636 1657.9163 3257.6791 6457.6198 ∞57.6

12 台形公式（第４章１）

13 台形公式（第４章１） 1 64 1 32 1 16 1818 1414

14 曲線による近似 シンプソン公式 シンプソン公式 ２区間毎を二次関数で近似し、積分する。 x y x0x0 x1x1 x2x2 x3x3 y0y0 y1y1 y3y3 y4y4 y2y2 x4x4 ２次関数

15 ２次関数の積分 x y x0x0 x1x1 x2x2 y0y0 h y1y1 y2y2 h I = h(y 0 + 4y 1 + y 2 ) I = h(y 0 + 4y 1 + y 2 )13

16 ２次関数の積分 x012 1 4 y y = x 2 I = x 2 dx =[ ] = I = x 2 dx =[ ] = 2 0 x3x333x3x3333 2 0 83 解析解 I = h(y 0 + 4y 1 + y 2 ) I = h(y 0 + 4y 1 + y 2 ) = 1 × (0 + 4 × 1 + 4)= = 1 × (0 + 4 × 1 + 4)= 13 13 83 数値解

17 シンプソン公式 x y x0x0 x1x1 x2x2 x3x3 y0y0 y1y1 y3y3 y4y4 x4x4 I = h{(y 0 + 4y 1 + y 2 )+ I = h{(y 0 + 4y 1 + y 2 )+ (y 2 + 4y 3 + y 4 ) + ……} (y 2 + 4y 3 + y 4 ) + ……}13 y2y2

18 シンプソン公式 x y x0x0 x1x1 x2x2 x3x3 y0y0 y1y1 y3y3 y4y4 x4x4 I = h{y 0 + 4(y 1 + y 3 +.. ) I = h{y 0 + 4(y 1 + y 3 +.. ) 2(y 2 + y 4 +.. ) + y n } 2(y 2 + y 4 +.. ) + y n }13 奇数番目を４ 倍 偶数番目を２ 倍 偶数番目 奇数番目 y2y2

19 Ｅｘｃｅｌで確認 ABC 1XY 20.0=A2^4+2*A2=B2 30.5=?*B3 41.0=?*B4 51.5=?*B5 62.0=?*B6 72.5=?*B7 83.0=B8 =(?/?)*sum(C2:C8) 区間の幅数値の合計 I = x 4 +2x dx = I = x 4 +2x dx = 3 0 57.6

20 Ｅｘｃｅｌで確認 ABC 1XY 20.0=A2^4+2*A2=B2 30.5=4*B3 41.0=2*B4 51.5=4*B5 62.0=2*B6 72.5=4*B7 83.0=B8 =(0.5/3)*sum(C2:C8) 初め 数値の合計 奇数番目 偶数番目 最後 I = x 4 +2x dx = I = x 4 +2x dx = 3 0 57.6

21 シンプソン計算例 x y 012 8 16 24 32 I = x 4 +2x dx = I = x 4 +2x dx = 3 0 57.6 3

22 nclude #include #defineN30 doublefunc_y( double ); int main( ){ doubley[N+1]; doublexa = 0.0, xb = 3.0; doublez1 = 0.0, z2 = 0.0, h, x, s; inti; h = (xb - xa)/(double)N;

23 for(i=0;i<=N;i++){ x = xa + h*(double)i; y[i] = func_y(x); } for(i=1;i<=N-1;i+=2)z1 += 4.0*y[i]; for(i=2;i<=N-2;i+=2)z2 += 2.0*y[i]; s = (h/3.0)*(y[0]+z1+z2+y[N]); printf("ANS = %8.4lf\n",s); return 0; } doublefunc_y(double x) { returnpow(x, 4.0) + 2.0*x; } 奇数番目を４ 倍 偶数番目を２ 倍 最初と最後は１ 倍

24 一次（台形）と二次（シンプソン） の精度 n台形 シンプソ ン 462.630957.7266 858.863657.6079 1657.916357.6005 3257.679157.600 6457.619857.600 ∞57.6

25 一次（台形）と二次（シンプソン） の精度 1 64 1 32 1 16 1818 1414