數值分析實驗之數值積分法（Python 代碼）

詳細實驗指導見上一篇，此處只寫內容啦

 

實驗內容

     選擇 y=arctan(x) 在0-1上的積分

 

• 復化simpson算法

from sympy import *
import math

def func(x):
    return math.atan(x)
def sum_fun_xk(xk, func):
    return sum([func(each) for each in xk])

def integral(a, b, n, func):
    h = (b - a)/float(n)
    xk = [a + i*h for i in range(1, n)]
    return h/2 * (func(a) + 2 * sum_fun_xk(xk, func) + func(b))

if __name__ == "__main__":
     
    a, b = 0, 1
    n = 50
    y=integral(a, b, n, func)
    print ("復化simpson公式：",y)

運行結果：

    將區間分成50份所得結果：

      

    將區間分成100份所得結果：

     

 

• 復化梯形算法

import numpy 
from scipy import integrate
import matplotlib.pyplot as plt
import math

#x [0,2pi]
def f(x):
    return math.atan(x)

def T(a,b,n=50):
    h = (b - a) / n
    temp = 0
    for i in range(1,n):
        x = a + i * h
        temp += 2 * f(x)
    return (b - a) / (2 * n) * (f(a) + temp + f(b))

#50 100 200 500 1000
def S(a,b,n=50):
    h = (b - a) / n
    temp1 = 0
    temp2 = 0
    for i in range(1,n):
        xk1 = a + h * i
        xk2 = a + h * (i + 1)
        xk12 = (xk1 + xk2) / 2
        temp1 += f(xk1)
        temp2 += f(xk12)
    temp2 += f((a + a + h) / 2)
    return (b - a) / (6 * n) * (f(a) + 4 * temp2 + 2 * temp1 + f(b))

if __name__ == '__main__':
    n =1000 #50 100 200 500 1000
    a = 0
    b = 1

    print ('Truth-value:',integrate.quad(f,a,b)[0])
    print ('T-Estimated-value:',T(a,b,n))
    print ('S-Estimated-value:',S(a,b,n))

運行結果：

      

   

體會

  首先同一方法的比較，我們將區間分成了50份和100份進行計算，得出結果：將區間分成100份的計算結果精度較高。其次，不同方法之間的比較，我們選用的復化梯形算法和復化辛普森算法對該積分進行計算，發現在將區間分成50份時，所得精度相差不大。總而言之，通過本次實驗，我對數值分析中相關的知識更加熟悉。其次，對各種編程軟件的使用也進一步了解，雖然java是上學期學的，但是通過本次課進行熟悉，收獲也頗多。