import numpy as np import math from scipy import integrate def f(x): return math.sin(x)*math.sin(x) #复化梯形法 def func(a,b,n,f): x = np.linspace(a,b,n+1) sum1 = 0 h =(b-a)/n for i in range(补充代码1): 补充代码2 return sum1 #复化辛普森法 def func1(a,b,n,f): x = np.linspace(a,b,n+1) sum1 = 0 h =(b-a)/n for i in range(补充代码3): 补充代码4 return sum1 #复化科特斯法 def func2(a,b,n,f): x = np.linspace(a,b,n+1) sum1 = 0 h =(b-a)/n for i in range(补充代码5): 补充代码6 return sum1 answer = func(0,1,100,f) answer1 = func1(0,1,100,f) answer2 = func2(0,1,100,f) print(answer,answer1,answer2) #integrate积分作为精确值 value, error = integrate.quad(f,0, 1) print(value,error) print("error: ", abs(answer-value), abs(answer1-value), abs(answer2-value))
时间: 2023-08-18 15:03:47 浏览: 166
代码中缺少部分内容,需要进行补充。以下是三种数值积分方法的完整代码:
```python
import numpy as np
import math
from scipy import integrate
def f(x):
return math.sin(x)*math.sin(x)
# 复化梯形法
def func(a, b, n, f):
x = np.linspace(a, b, n+1)
sum1 = 0
h = (b-a)/n
for i in range(1, n):
sum1 += f(x[i])
return h*(0.5*f(a) + sum1 + 0.5*f(b))
# 复化辛普森法
def func1(a, b, n, f):
x = np.linspace(a, b, n+1)
sum1 = 0
h = (b-a)/n
for i in range(1, n):
if i%2 == 0:
sum1 += 2*f(x[i])
else:
sum1 += 4*f(x[i])
return h*(f(a) + sum1 + f(b))/3
# 复化科特斯法
def func2(a, b, n, f):
x = np.linspace(a, b, n+1)
sum1 = 0
h = (b-a)/n
for i in range(1, n):
if i%3 == 0:
sum1 += 2*f(x[i])
else:
sum1 += 3*f(x[i])
return h*(f(a) + sum1 + f(b))*3/8
answer = func(0,1,100,f)
answer1 = func1(0,1,100,f)
answer2 = func2(0,1,100,f)
print(answer,answer1,answer2)
# integrate积分作为精确值
value, error = integrate.quad(f, 0, 1)
print(value,error)
print("error: ", abs(answer-value), abs(answer1-value), abs(answer2-value))
```
其中,补充代码1,补充代码3,补充代码5 分别表示循环的次数,可设置为 `for i in range(1, n):`。补充代码2,补充代码4,补充代码6 分别表示不同积分方法中的积分公式,具体可参考课程中的讲解。
阅读全文