梯形公式的数值积分的Python程序

本文通过实例展示了如何使用梯形积分算法计算函数f(t)的积分,并探讨了如何根据精度需求调整积分区域块数,以确保截断误差小于特定阈值。从10%误差到10^-6误差的精度提升过程中,逐步揭示了划分步长对积分结果的影响。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

在这里插入图片描述
在这里插入图片描述

代码

例子

下面的代码的意图是改变积分区域的梯形块数,得到不同精度的积分值。
以对函数f(t) = 2000×ln(140000/(140000-2100×t))-9.8×t 进行积分为例,积分区域为[8,30]

from sympy import *

def f(t):
    f = 2000*log(140000/(140000-2100*t))-9.8*t
    return f

x = symbols('x')
truth = integrate(f(x),(x,8,30)).evalf()  
print(truth) #真值

n = 10 #步长,就是将(a,b)区间分为多少个块
a = 8
b = 30
h = (b-a)/n
tra_result = 0
for i in range(n):
    tra_result += 1/2*h*(f(a+i*h)+f(a+(i+1)*h))   #梯形积分算法
print(tra_result) #梯形积分值

结果:

11061.3355350810
11069.5835420617

下面来探讨,积分区域内,需要划分多少块,梯形积分算法计算出来的积分值才能与真值之间的截断误差小于0.1。

from sympy import *

def f(t):
    f = 2000*log(140000/(140000-2100*t))-9.8*t
    return f

x = symbols('x')
truth = integrate(f(x),(x,8,30)).evalf()  #真值

a = 8
b = 30
n = 0 #步长,就是将(a,b)区间分为多少个块
while True:
    n += 1
    h = (b-a)/n
    tra_result = 0
    for i in range(n):
        tra_result += 1/2*h*(f(a+i*h)+f(a+(i+1)*h))   #梯形积分算法
    R = abs(truth - tra_result)
    if R <= 10**(-1):
        print(f'需要将积分区域划分{n}块,截断误差<0.1')
        break

结果:

需要将积分区域划分91,截断误差<0.1

最后我们来看看真值和梯形积分算法的积分值误差小于10^(-6),需要将积分区域划分多少块?

from sympy import *

def f(t):
    f = 2000*log(140000/(140000-2100*t))-9.8*t
    return f

x = symbols('x')
truth = integrate(f(x),(x,8,30)).evalf()  #真值

a = 8
b = 30
n = 0 #步长,就是将(a,b)区间分为多少个块
while True:
    n += 1
    h = (b-a)/n
    tra_result = 0
    for i in range(n):
        tra_result += 1/2*h*(f(a+i*h)+f(a+(i+1)*h))   #梯形积分算法
    error = abs((truth - tra_result)/truth)
    if error <= 10**(-6):
        print(f'需要将积分区域划分{n}块')
        break

结果:

需要将积分区域划分274

可以对比一下Simpson积分算法:抛物线公式即辛卜生(Simpson)公式的数值积分的Python程序

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值