日本免费高清视频-国产福利视频导航-黄色在线播放国产-天天操天天操天天操天天操|www.shdianci.com

學無先后,達者為師

網站首頁 編程語言 正文

利用Python實現數值積分的方法_python

作者:趙卓不凡? ? 更新時間: 2022-04-15 編程語言

1. 栗子

為了加深大家的印象,首先我們來看個例子:

圖示如下:

2. 矩形計算面積

我們知道,在數學中,積分運算表示上述曲線和x軸圍成的封閉區域的面積,為此,我們在數值預算中,來近似計算上述區域的面積,最直觀的想法就是拆成一個個小的矩形來計算對應的面積。

2.1 左側邊長計算面積

為了計算每個小矩形的面積,設計到邊長高的選擇,這里我么以左側函數取值作為對應矩形的高來計算相應的小矩形的面積,圖示如下:

對應的代碼如下:

import numpy as np
x = np.linspace(0, 3, 1001)
f = lambda x: x**3 - 4*x**2 + 4*x + 2
a = 0.5
b = 2.5
Ax = np.linspace(a, b, 101)
Ay = f(Ax)
def defInt_left(f, a, b, N):
? ? # left-hand point
? ? result = 0; FX = []; Xn = []
? ? dx = abs(b - a)/N
? ? while a < b:
? ? ? ? result += f(a)*dx
? ? ? ? FX += [f(a)]
? ? ? ? Xn += [a]
? ? ? ? a += dx
? ? return result, FX, Xn, dx
N = 4
I_left, FX, Xn, dx = defInt_left(f, a, b, N)
print(I_left)

上述代碼中,我們將橫坐標拆分為4小份,也就是拆分成4個小矩形,然后使用函數左側的點坐位小矩形的高,上述代碼的運行結果如下:

5.25

2.2 右側邊長計算面積

這里和上述原理類似,只不過每個小矩形的高采用右側邊長函數取值來近似計算,圖例如下:

樣例代碼如下:

def defInt_right(f, a, b, N):
? ? # right-hand point
? ? result = 0; FX = []; Xn = []
? ? dx = abs(b - a)/N
? ? while a < b:
? ? ? ? result += f(a + dx)*dx
? ? ? ? FX += [f(a + dx)]
? ? ? ? Xn += [a]
? ? ? ? a += dx
? ? return result, FX, Xn, dx

N = 4
I_right, FX, Xn, dx = defInt_right(f, a, b, N)
print(I_right)

運行結果如下:

5.0

2.3 中值邊長計算面積

看了上述兩種近似計算方式,有同學就說有取左側點算出來面積大的,有取右側點算出來面積小的,那干脆折中一下,我們來以中值坐位矩形的高來計算對應的面積。圖例如下:

代碼實現如下:

def defInt_middle(f, a, b, N):
? ? # middle point
? ? result = 0; FX = []; Xn = []
? ? dx = abs(b - a)/N
? ? while a < b:
? ? ? ? result += f(a + dx/2)*dx
? ? ? ? FX += [f(a + dx/2)]
? ? ? ? Xn += [a]
? ? ? ? a += dx
? ? return result, FX, Xn, dx

N = 4
I_mid, FX, Xn, dx = defInt_middle(f, a, b, N)
print(I_mid)

運行結果如下:

5.0625

3. 梯形計算面積

讀到這里的同學可能會思考,既然可以將封閉區域劃分成一個個的小矩形,那當然也可以將其劃分成梯形來近似計算相應的面積,圖例如下:

樣例代碼如下:

def defInt_trapezoid(f, a, b, N):
? ? # trapezoidal rule
? ? result = 0; FXa, FXb = [], []; Xn = []
? ? dx = abs(b - a)/N
? ? while a < b:
? ? ? ? result += (f(a) + f(a + dx))*dx/2
? ? ? ? FXa += [f(a)]; FXb += [f(a + dx)]
? ? ? ? Xn += [a]
? ? ? ? a += dx
? ? return result, FXa, FXb, Xn, dx

N = 4
I_trap, FXa, FXb, Xn, dx = defInt_trapezoid(f, a, b, N)
print(I_trap)

運行結果如下:

5.125

4. 真值比對

最后,我們來針對不同的N來講封閉區域劃分成對應的小份,分別針對性的計算上述四種方式的積分值,樣例代碼如下:

Nx = range(1, 11)
I1, I2, I3, I4 = [], [], [], []
for Ni in Nx:
? ? i1, *_ = defInt_left(f, a, b, Ni); I1 += [i1];
? ? i2, *_ = defInt_right(f, a, b, Ni); I2 += [i2];
? ? i3, *_ = defInt_middle(f, a, b, Ni); I3 += [i3];
? ? i4, *_ = defInt_trapezoid(f, a, b, Ni); I4 += [i4];

最后將其與真值進行對比,如下:

可以看出,隨著劃分區域的增多,采用梯形計算面積方式最逼近真值。

5. 總結

本文重點介紹了使用不同面積劃分方法來近似計算積分取值的原理和相應的代碼實現,其中采用梯形計算面積的方式隨著劃分子區域數目的增加最接近真值,推薦大家使用。

原文鏈接:https://blog.csdn.net/sgzqc/article/details/122817967

欄目分類
最近更新