数值分析实验之数值积分法（Python 代码）

Posted 2020-12-19 ぺあ紫泪°冰封ヤ

详细实验指导见上一篇，此处只写内容啦

 

实验内容

     选择 y=arctan(x) 在0-1上的积分

 

• 复化simpson算法

 1 from sympy import *
 2 import math
 3 
 4 def func(x):
 5     return math.atan(x)
 6 def sum_fun_xk(xk, func):
 7     return sum([func(each) for each in xk])
 8 
 9 def integral(a, b, n, func):
10     h = (b - a)/float(n)
11     xk = [a + i*h for i in range(1, n)]
12     return h/2 * (func(a) + 2 * sum_fun_xk(xk, func) + func(b))
13 
14 if __name__ == "__main__":
15      
16     a, b = 0, 1
17     n = 50
18     y=integral(a, b, n, func)
19     print ("复化simpson公式：",y)

运行结果：

    将区间分成50份所得结果：

      

    将区间分成100份所得结果：

     

 

• 复化梯形算法

 1 import numpy 
 2 from scipy import integrate
 3 import matplotlib.pyplot as plt
 4 import math
 5 
 6 #x [0,2pi]
 7 def f(x):
 8     return math.atan(x)
 9 
10 def T(a,b,n=50):
11     h = (b - a) / n
12     temp = 0
13     for i in range(1,n):
14         x = a + i * h
15         temp += 2 * f(x)
16     return (b - a) / (2 * n) * (f(a) + temp + f(b))
17 
18 #50 100 200 500 1000
19 def S(a,b,n=50):
20     h = (b - a) / n
21     temp1 = 0
22     temp2 = 0
23     for i in range(1,n):
24         xk1 = a + h * i
25         xk2 = a + h * (i + 1)
26         xk12 = (xk1 + xk2) / 2
27         temp1 += f(xk1)
28         temp2 += f(xk12)
29     temp2 += f((a + a + h) / 2)
30     return (b - a) / (6 * n) * (f(a) + 4 * temp2 + 2 * temp1 + f(b))
31 
32 if __name__ == \'__main__\':
33     n =1000 #50 100 200 500 1000
34     a = 0
35     b = 1
36 
37     print (\'Truth-value:\',integrate.quad(f,a,b)[0])
38     print (\'T-Estimated-value:\',T(a,b,n))
39     print (\'S-Estimated-value:\',S(a,b,n))

运行结果：

      

   

体会

  首先同一方法的比较，我们将区间分成了50份和100份进行计算，得出结果：将区间分成100份的计算结果精度较高。其次，不同方法之间的比较，我们选用的复化梯形算法和复化辛普森算法对该积分进行计算，发现在将区间分成50份时，所得精度相差不大。总而言之，通过本次实验，我对数值分析中相关的知识更加熟悉。其次，对各种编程软件的使用也进一步了解，虽然java是上学期学的，但是通过本次课进行熟悉，收获也颇多。