快捷導(dǎo)航

利用Python實(shí)現(xiàn)數(shù)值積分的方法

更新時(shí)間：2022年02月10日 08:43:40 作者：趙卓不凡?

這篇文章主要介紹了利用Python實(shí)現(xiàn)數(shù)值積分。本文主要用于對(duì)比使用Python來實(shí)現(xiàn)數(shù)學(xué)中積分的幾種計(jì)算方式，并和真值進(jìn)行對(duì)比，加深大家對(duì)積分運(yùn)算實(shí)現(xiàn)方式的理解

1. 栗子

為了加深大家的印象，首先我們來看個(gè)例子：

圖示如下：

2. 矩形計(jì)算面積

我們知道，在數(shù)學(xué)中，積分運(yùn)算表示上述曲線和x軸圍成的封閉區(qū)域的面積，為此，我們?cè)跀?shù)值預(yù)算中，來近似計(jì)算上述區(qū)域的面積，最直觀的想法就是拆成一個(gè)個(gè)小的矩形來計(jì)算對(duì)應(yīng)的面積。

2.1 左側(cè)邊長(zhǎng)計(jì)算面積

為了計(jì)算每個(gè)小矩形的面積，設(shè)計(jì)到邊長(zhǎng)高的選擇，這里我么以左側(cè)函數(shù)取值作為對(duì)應(yīng)矩形的高來計(jì)算相應(yīng)的小矩形的面積，圖示如下：

對(duì)應(yīng)的代碼如下：

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)

上述代碼中，我們將橫坐標(biāo)拆分為4小份，也就是拆分成4個(gè)小矩形，然后使用函數(shù)左側(cè)的點(diǎn)坐位小矩形的高，上述代碼的運(yùn)行結(jié)果如下：

5.25

2.2 右側(cè)邊長(zhǎng)計(jì)算面積

這里和上述原理類似，只不過每個(gè)小矩形的高采用右側(cè)邊長(zhǎng)函數(shù)取值來近似計(jì)算，圖例如下：

樣例代碼如下：

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)

運(yùn)行結(jié)果如下：

5.0

2.3 中值邊長(zhǎng)計(jì)算面積

看了上述兩種近似計(jì)算方式，有同學(xué)就說有取左側(cè)點(diǎn)算出來面積大的，有取右側(cè)點(diǎn)算出來面積小的，那干脆折中一下，我們來以中值坐位矩形的高來計(jì)算對(duì)應(yīng)的面積。圖例如下：

代碼實(shí)現(xiàn)如下：

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)

運(yùn)行結(jié)果如下：

5.0625

3. 梯形計(jì)算面積

讀到這里的同學(xué)可能會(huì)思考，既然可以將封閉區(qū)域劃分成一個(gè)個(gè)的小矩形，那當(dāng)然也可以將其劃分成梯形來近似計(jì)算相應(yīng)的面積，圖例如下：

樣例代碼如下：

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)

運(yùn)行結(jié)果如下：

5.125

4. 真值比對(duì)

最后，我們來針對(duì)不同的N來講封閉區(qū)域劃分成對(duì)應(yīng)的小份，分別針對(duì)性的計(jì)算上述四種方式的積分值，樣例代碼如下：

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];

最后將其與真值進(jìn)行對(duì)比，如下：

可以看出，隨著劃分區(qū)域的增多，采用梯形計(jì)算面積方式最逼近真值。

5. 總結(jié)

本文重點(diǎn)介紹了使用不同面積劃分方法來近似計(jì)算積分取值的原理和相應(yīng)的代碼實(shí)現(xiàn),其中采用梯形計(jì)算面積的方式隨著劃分子區(qū)域數(shù)目的增加最接近真值，推薦大家使用。

到此這篇關(guān)于利用Python實(shí)現(xiàn)數(shù)值積分的方法的文章就介紹到這了,更多相關(guān)Python實(shí)現(xiàn)數(shù)值積分內(nèi)容請(qǐng)搜索腳本之家以前的文章或繼續(xù)瀏覽下面的相關(guān)文章希望大家以后多多支持腳本之家！

您可能感興趣的文章: