龙贝格算法求数值积分的Python程序

最新推荐文章于 2022-10-27 17:10:51 发布
最新推荐文章于 2022-10-27 17:10:51 发布
阅读量2.4k 收藏 31
点赞数 6
分类专栏： 算法 文章标签： python
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/weixin_48615832/article/details/111496112
版权

Ronberg Integration

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

分成n等份时辛普森公式的值

在这里插入图片描述

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) #真值

def T(n,a,b):
    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))  # 梯形积分算法
    return tra_result

a = 8
b = 30
list = []
for k in range(4):
    m = 2**k
    list.append(T(m,a,b))
print(list)

while True:
    new_list = []
    for i in range(len(list)-1):
        R = list[i+1] + (list[i+1]-list[i])/3
        new_list.append(R)
    print(new_list)
    if len(new_list) == 1:
        print(new_list[0])
        break
    else:
        list = new_list

结果：

11061.3355350810
[11868.3481898411, 11266.3742932594, 11112.8206763693, 11074.2212976601]
[11065.7163277322, 11061.6361374059, 11061.3548380903]
[11060.2760739638, 11061.2610716518]
[11061.5894042144]
11061.5894042144

辛普森法的积分值Sn

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) #真值

def T(n,a,b):
    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))  # 梯形积分算法
    return tra_result

n = 0 #步长,就是将(a,b)区间分为多少个块
a = 8
b = 30
while True:
    n += 1
    S = T(2*n, a, b) + (T(2*n, a, b) - T(n, a, b)) / 3
    error = abs(truth - S) / truth
    if error <= 10 ** (-6):
        print('S为：', S)
        print('n为：', n)
        break

结果：

11061.3355350810
S为： 11061.3434684075
n为： 5

柯特斯法的积分值Cn

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) #真值

def T(n,a,b):
    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))  # 梯形积分算法
    return tra_result

def S(n,a,b):
    I = T(2*n, a, b) + (T(2*n, a, b) - T(n, a, b)) / (4**1-1)
    return I

n = 1 #步长,就是将(a,b)区间分为多少个块
a = 8
b = 30
C = S(2*n, a, b) + (S(2*n, a, b) - S(n, a, b)) / (4**2-1)
print(C)

结果：

11061.3355350810
11061.3641247175

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) #真值

def T(n,a,b):
    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))  # 梯形积分算法
    return tra_result

def S(n,a,b):
    I = T(2*n, a, b) + (T(2*n, a, b) - T(n, a, b)) / (4**1-1)
    return I

n = 0 #步长,就是将(a,b)区间分为多少个块
a = 8
b = 30
while True:
    n += 1
    C = S(2*n, a, b) + (S(2*n, a, b) - S(n, a, b)) / (4**2-1)
    error = abs(truth - C) / truth
    if error <= 10 ** (-6):
        print('C为：', C)
        print('n为：', n)
        break

结果：

11061.3355350810
C为： 11061.3360848026
n为： 2

龙贝格（Romberg）公式Rn

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) #真值

def T(n,a,b):
    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))  # 梯形积分算法
    return tra_result

def S(n,a,b):
    I = T(2*n, a, b) + (T(2*n, a, b) - T(n, a, b)) / (4**1-1)
    return I

def C(n,a,b):
    C = S(2*n, a, b) + (S(2*n, a, b) - S(n, a, b)) / (4**2-1)
    return C

n = 1 #步长,就是将(a,b)区间分为多少个块
a = 8
b = 30
R = C(2*n, a, b) + (C(2*n, a, b) - C(n, a, b)) / (4**3-1)
print(R)

结果：

11061.3355350810
11061.3356397246

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) #真值

def T(n,a,b):
    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))  # 梯形积分算法
    return tra_result

def S(n,a,b):
    I = T(2*n, a, b) + (T(2*n, a, b) - T(n, a, b)) / (4**1-1)
    return I

def C(n,a,b):
    C = S(2*n, a, b) + (S(2*n, a, b) - S(n, a, b)) / (4**2-1)
    return C

n = 0 #步长,就是将(a,b)区间分为多少个块
a = 8
b = 30
while True:
    n += 1
    R = C(2*n, a, b) + (C(2*n, a, b) - C(n, a, b)) / (4**3-1)
    error = abs(truth - R) / truth
    if error <= 10 ** (-6):
        print('R为：', R)
        print('n为：', n)
        break

结果：

11061.3355350810
R为： 11061.3356397246
n为： 1
Fo*(Bi)
关注
专栏目录
龙贝格求积分算法例题_龙贝格求积公式（Python实现）
weixin_29716207的博客
12-22 2125
#Author：glm233 #这个程序就是一个黑箱接口，把要进行龙贝格求积的函数放在func函数里，可以自行修改，然后范围就是在Romberg里改 import math ''' 给定一个函数，如：f(x)= x^(3/2），和积分上下限a,b，用机械求积Romberg公式求积分。 ''' import numpy as np def func(x): return x**2+...
python实现数值分析之龙贝格求积公式
weixin_42920648的博客
10-24 6118
复合梯形公式的提出： 1.首先，什么是梯形公式： 梯形公式表明：f(x)在[a,b]两点之间的积分（面积），近似地可以用一个梯形的面积表示。 2.显然，这个梯形公式对于不同的f(x)而言，其代数精度不同。为了能适合更多的f(x)，我们一般使用牛顿-科特斯公式其中比较高次的公式来进行数值求积。但高次的缺陷是当次数大于8次，求积公式就会不稳定。因此，我们用于数值积分的牛顿-科特斯公式通常是一次的梯形...
Python龙贝格法求积分实例
09-17
今天小编就为大家分享一篇Python龙贝格法求积分实例，具有很好的参考价值，希望对大家有所帮助。一起跟随小编过来看看吧
变步长梯形公式数值积分的Python程序
Fo*(Bi)的博客
12-21 2280
基本代码： 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) #真值 def T(n,a,b): h = (b-a)/n tra_result = 0 for i in range(n): ..
变步长复化梯形法的Python实现
Jaryn62的博客
06-09 1238
变步长复化梯形法的Python
python实现龙贝格积分算法
12-22
用python简洁的实现龙贝格积分算法，可供基本参考
[计算方法作业]利用python中matplotlib实现绘制龙贝格公式和变步长梯形法图像
12-25
[计算方法作业]利用python中matplotlib实现绘制龙贝格公式和变步长梯形法图像,利用龙贝格公式及变步长积分两种算法计算定积分，用python中的matplotlib库绘制图像
龙贝格求积分算法python
最新发布
05-25
龙贝格求积分算法是一种数值积分方法，可以用于计算定积分的近似值。下面是用Python实现的代码： ```python def func(x): return x ** 2 + 1 # 定义被积函数 def romberg(func, a, b, eps=1e-6, max_steps=10): ...
计算方法之用变步长梯形求积公式求定积分
he91-com
06-11 1722
/* * * 用变步长梯形求积公式求定积分 * 1 * ∫ （x/(4+x^2)）dx * 0 */ #include &lt;stdio.h&gt; #include &lt;string.h&gt; #include &lt;math.h&gt; #include &lt;conio.h&gt; #include &lt;stdlib.h&gt; #defin...
复化梯形求积分实例——用Python进行数值计算
09-18
今天小编就为大家分享一篇复化梯形求积分实例——用Python进行数值计算，具有很好的参考价值，希望对大家有所帮助。一起跟随小编过来看看吧
龙贝格算法计算积分sinx/x积分
06-18
数值积分 用龙贝格算法计算下面积分sinx/x,积分限0~1，精度要求为10-6.
变步长梯形算法
06-28
变步长梯形算法，本科生作业，CPP文件，可直接运行。
数值分析作业：梯形公式
02-05
C语言数值分析上机实践作业：用梯形公式求积分1/(1+x*x)
Romberg(龙贝格)数值积分算法较高效的python实现
K1052176873的博客
04-08 6268
1. 原理与公式 2. Python代码实现 # coding=utf-8 import numpy as np import math # 二分法梯形公式 # func:需要积分的函数 # x_min: 积分下限 # x_max: 积分上限 # epoch: 二分次数 def compute_Tn(func, x_min=0, x_max=1, epoch=10): Tn_list = [] Tn = 0 h0 = x_max - x_min # 积分区间的长度,即
python辛普森积分_Python龙贝格法求积分实例
weixin_39939303的博客
12-02 375
我就废话不多说了，直接上代码吧！# 龙贝格法求积分import matha=0 # 积分下限b=1 # 积分上限eps=10**-5 # 精度T=[] # 复化梯形序列S=[] # Simpson序列C=[] # Cotes序列R=[] # Romberg序列def func(x): # 被积函数y=math.exp(-x)return ydef Romberg(a,...
变步长梯形法与龙贝格算法
wujing1_1的博客
11-01 1830
文章目录1. 变步长梯形法算法描述流程图代码实现龙贝格算法算法描述例子代码实现 1. 变步长梯形法 提出背景： 复化求积公式虽然能提高精度，但需要给出步长，步长精度太大则精度低，步长太小则计算量大，难以找到一个合适的步长（划分成的小区间的个数） 算法描述 1.对所有已存在的子区间进行二分化，区间数由n变为2n 2.利用区间数为n时的积分值Tn以及新增的节点（即原来各子区间的中点）递推出区间数为2n时的积分值T2n 3.利用两次计算结果的差来估计误差，直到满足精度 公式如下： 流程图 代码实现 /**
变步长梯形积分算法求解函数定积分
Dear_Jia的博客
12-09 6630
算法基本原理：把原区间分为一系列小期间(n份)，在每个小区间上都用小的梯形面积来近似代替原函数的积分，当小区间足够小时，就可以得到原来积分的近似值。但是这个过程中有一个问题是步长的取值，步长太大精度难以保证，步长太小会导致计算量的增加，所以，实际计算中常常采用变步长的方法，使得在步长逐次减小的过程中，求得的积分结果满足要求的精度为止。        首先，给出两个计算公式           ...
用python实现龙贝格算法
m0_64027371的博客
10-27 1347
用python实现龙贝格算法
python梯形公式_復化梯形求積分——用Python進行數值計算
weixin_42399813的博客
02-10 416
用程序來求積分的方法有很多，這篇文章主要是有關牛頓-科特斯公式。學過插值算法的同學最容易想到的就是用插值函數代替被積分函數來求積分，但實際上在大部分場景下這是行不通的。插值函數一般是一個不超過n次的多項式，如果用插值函數來求積分的話，就會引進高次多項式求積分的問題。這樣會將原來的求積分問題帶到另一個求積分問題：如何求n次多項式的積分，而且當次數變高時，會出現龍悲歌現象，誤差反而可能會增大，並且高次...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

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

抵扣说明：

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

余额充值