Tích hợp vùng dưới đường cong python

Tích hợp chung ()

Hàm được cung cấp để tích phân hàm một biến giữa hai điểm. Điểm có thể là \(\pm\infty\) ( \(\pm\)

>>> print(abs(result[0]-I))
1.03761443881e-11

>>> print(abs(result[0]-I))
1.03761443881e-11

\[I=\int_{0}^{4. 5}J_{2. 5}\trái(x\phải)\, dx. \]

Điều này có thể được tính toán bằng cách sử dụng

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

>>> print(abs(result[0]-I))
1.03761443881e-11

Đối số đầu tiên của quad là một đối tượng Python “có thể gọi được” (i. e. , một hàm, phương thức hoặc thể hiện của lớp). Lưu ý việc sử dụng hàm lambda- trong trường hợp này làm đối số. Hai đối số tiếp theo là các giới hạn của tích hợp. Giá trị trả về là một bộ, với phần tử đầu tiên giữ giá trị ước tính của tích phân và phần tử thứ hai giữ giới hạn trên của lỗi. Lưu ý rằng trong trường hợp này, giá trị thực của tích phân này là

\[I=\sqrt{\frac{2}{\pi}}\left(\frac{18}{27}\sqrt{2}\cos\left(4. 5\right)-\frac{4}{27}\sqrt{2}\sin\left(4. 5\right)+\sqrt{2\pi}\textrm{Si}\left(\frac{3}{\sqrt{\pi}}\right)\right),\]

ở đâu

\[\textrm{Si}\left(x\right)=\int_{0}^{x}\sin\left(\frac{\pi}{2}t^{2}\right)\, dt. \]

là tích phân sin Fresnel. Lưu ý rằng tích phân được tính toán bằng số nằm trong \(1. 04\times10^{-11}\) của kết quả chính xác — thấp hơn rất nhiều so với giới hạn lỗi được báo cáo.

Nếu chức năng tích hợp có các tham số bổ sung, chúng có thể được cung cấp trong đối số args. Giả sử rằng tích phân sau sẽ được tính

\[I(a,b)=\int_{0}^{1} ax^2+b \, dx. \]

Tích phân này có thể được đánh giá bằng cách sử dụng đoạn mã sau

________số 8

Số lượng đầu vào vô hạn cũng được cho phép bằng cách sử dụng \(\pm\)

>>> print(abs(result[0]-I))
1.03761443881e-11

\[E_{n}\left(x\right)=\int_{1}^{\infty}\frac{e^{-xt}}{t^{n}}\, dt. \]

là mong muốn (và thực tế là tích phân này có thể được tính là

>>> from scipy.integrate import quad
>>> def integrand(x, a, b):
..     return a*x**2 + b
...
>>> a = 2
>>> b = 1
>>> I = quad(integrand, 0, 1, args=(a,b))
>>> I
(1.6666666666666667, 1.8503717077085944e-14)

>>> from scipy.integrate import quad
>>> def integrand(x, a, b):
..     return a*x**2 + b
...
>>> a = 2
>>> b = 1
>>> I = quad(integrand, 0, 1, args=(a,b))
>>> I
(1.6666666666666667, 1.8503717077085944e-14)

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

Hàm được tích hợp thậm chí có thể sử dụng đối số quad (mặc dù giới hạn lỗi có thể đánh giá thấp lỗi do lỗi số có thể xảy ra trong tích phân từ việc sử dụng ). Tích phân trong trường hợp này là

\[I_{n}=\int_{0}^{\infty}\int_{1}^{\infty}\frac{e^{-xt}}{t^{n}}\, dt\, dx . \]

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

Ví dụ cuối cùng này cho thấy rằng nhiều tích hợp có thể được xử lý bằng các lệnh gọi lặp lại tới

Tích hợp bội chung ( , , )

Cơ chế tích hợp kép và ba đã được gói gọn trong các chức năng và. Các hàm này lấy hàm để tích phân và bốn hoặc sáu đối số tương ứng. Các giới hạn của tất cả các tích phân bên trong cần được xác định là các hàm

Một ví dụ về việc sử dụng tích phân kép để tính toán một số giá trị của \(I_{n}\) được hiển thị bên dưới.

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

Như ví dụ cho các giới hạn không đổi, xét tích phân

\[I=\int_{y=0}^{1/2}\int_{x=0}^{1-2y} x y \, dx\, dy=\frac{1}{96}. \]

Tích phân này có thể được ước tính bằng cách sử dụng biểu thức bên dưới (Lưu ý việc sử dụng các hàm lambda không cố định cho giới hạn trên của tích phân bên trong)

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

Để tích hợp n-fold, scipy cung cấp chức năng. Các giới hạn tích hợp là một đối tượng có thể lặp lại. danh sách các giới hạn không đổi hoặc danh sách các hàm cho giới hạn tích hợp không đổi. Thứ tự của tích phân (và do đó là các giới hạn) là từ tích phân trong cùng đến tích phân ngoài cùng

Tích phân trên

\[I_{n}=\int_{0}^{\infty}\int_{1}^{\infty}\frac{e^{-xt}}{t^{n}}\, dt\, dx

có thể được tính như

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

Lưu ý rằng thứ tự của các đối số cho f phải khớp với thứ tự của các giới hạn tích phân; . e. , tích phân bên trong của \(t\) nằm trên khoảng \([1, \infty . and the outer integral with respect to \(x\) is on the interval \([0, \infty]\).

Giới hạn tích hợp không cố định có thể được xử lý theo cách tương tự;

\[I=\int_{y=0}^{1/2}\int_{x=0}^{1-2y} x y \, dx\, dy=\frac{1}{96}. \]

có thể được đánh giá bằng

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

đó là kết quả tương tự như trước đây

cầu phương Gaussian

Một vài chức năng cũng được cung cấp để thực hiện cầu phương Gaussian đơn giản trong một khoảng thời gian cố định. Đầu tiên là , thực hiện cầu phương Gaussian có thứ tự cố định. Hàm thứ hai là , thực hiện phép cầu phương Gaussian của nhiều bậc cho đến khi chênh lệch trong ước tính tích phân nằm dưới một số dung sai do người dùng cung cấp. Cả hai hàm này đều sử dụng mô-đun

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

Tích hợp Romberg

Phương pháp của Romberg là một phương pháp khác để đánh giá bằng số một tích phân. Xem chức năng trợ giúp để biết thêm chi tiết

Tích hợp sử dụng mẫu

Nếu các mẫu cách đều nhau và số lượng mẫu có sẵn là \(2^{k}+1\) đối với một số . Tích hợp Romberg sử dụng quy tắc hình thang ở các kích thước bước liên quan đến lũy thừa của hai và sau đó thực hiện phép ngoại suy Richardson trên các ước tính này để tính gần đúng tích phân với độ chính xác cao hơn. \(k\), then Romberg integration can be used to obtain high-precision estimates of the integral using the available samples. Romberg integration uses the trapezoid rule at step-sizes related by a power of two and then performs Richardson extrapolation on these estimates to approximate the integral with a higher degree of accuracy.

Trong trường hợp các mẫu có khoảng cách tùy ý, hai chức năng và có sẵn. Họ đang sử dụng các công thức Newton-Coates bậc 1 và bậc 2 tương ứng để thực hiện tích phân. Quy tắc hình thang xấp xỉ hàm dưới dạng một đường thẳng giữa các điểm liền kề, trong khi quy tắc Simpson xấp xỉ hàm giữa ba điểm liền kề dưới dạng parabol

Đối với một số lẻ các mẫu cách đều nhau, quy tắc Simpson là chính xác nếu hàm là đa thức bậc 3 trở xuống. Nếu các mẫu không cách đều nhau thì kết quả chỉ chính xác nếu hàm là đa thức bậc 2 trở xuống

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

Điều này tương ứng chính xác với

\[\int_{1}^{4} x^2 \, dx = 21,\]

trong khi tích hợp chức năng thứ hai

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

không tương ứng với

\[\int_{1}^{4} x^3 \, dx = 63. 75\]

vì bậc của đa thức ở f2 lớn hơn hai

Tích hợp nhanh hơn bằng cách sử dụng các chức năng gọi lại cấp thấp

Người dùng muốn giảm thời gian tích hợp có thể chuyển một con trỏ hàm C đến , hoặc và nó sẽ được tích hợp và trả về kết quả bằng Python. Sự gia tăng hiệu suất ở đây phát sinh từ hai yếu tố. Cải tiến chính là đánh giá chức năng nhanh hơn, được cung cấp bằng cách biên dịch chính chức năng đó. Ngoài ra, chúng tôi có một sự tăng tốc được cung cấp bằng cách loại bỏ các lệnh gọi hàm giữa C và Python trong. Phương pháp này có thể cải thiện tốc độ ~2 lần đối với các hàm tầm thường như sin nhưng có thể tạo ra các cải tiến đáng chú ý hơn nhiều (10x+) đối với các hàm phức tạp hơn. Sau đó, tính năng này hướng đến người dùng có tích hợp số chuyên sâu sẵn sàng viết một chút C để giảm đáng kể thời gian tính toán

Cách tiếp cận có thể được sử dụng, ví dụ, thông qua một vài bước đơn giản

1. ) Viết một hàm tích phân trong C với chữ ký hàm là

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> import scipy.integrate as integrate
>>> import scipy.special as special
>>> result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)
>>> result
(1.1178179380783249, 7.8663172481899801e-09)

2. ) Bây giờ hãy biên dịch tệp này thành thư viện dùng chung/động (tìm kiếm nhanh sẽ giúp ích cho việc này vì nó phụ thuộc vào hệ điều hành). Người dùng phải liên kết bất kỳ thư viện toán học nào, v.v. , được sử dụng. Trên linux, điều này trông giống như

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

Thư viện đầu ra sẽ được gọi là

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

3. ) Tải thư viện dùng chung vào Python bằng cách sử dụng và đặt

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

4. ) Bây giờ hãy tích hợp chức năng thư viện như bình thường, ở đây sử dụng

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

Bộ dữ liệu Python được trả về như mong đợi trong một khoảng thời gian ngắn hơn. Tất cả các tham số tùy chọn có thể được sử dụng với phương pháp này bao gồm chỉ định điểm kỳ dị, giới hạn vô hạn, v.v.

Phương trình vi phân thường ()

Tích hợp một tập hợp các phương trình vi phân thông thường (ODE) với các điều kiện ban đầu là một ví dụ hữu ích khác. Hàm này có sẵn trong SciPy để tích hợp phương trình vi phân vectơ bậc nhất

\[\frac{d\mathbf{y}}{dt}=\mathbf{f}\left(\mathbf{y},t\right),\]

với các điều kiện ban đầu \(\mathbf{y}\left(0\right)=y_{0}\) , trong đó < . \) \(\mathbf{y}\) is a length \(N\) vector and \(\mathbf{f}\) is a mapping from \(\mathcal{R}^{N}\) to \(\mathcal{R}^{N}.\) Phương trình vi phân thường cấp cao hơn luôn có thể được rút gọn thành phương trình vi phân loại này bằng cách đưa đạo hàm trung gian vào \(\mathbf{y . vector.

Ví dụ: giả sử muốn tìm nghiệm của phương trình vi phân cấp hai sau đây

\[\frac{d^{2}w}{dz^{2}}-zw(z)=0\]

với điều kiện ban đầu \(w\left(0\right)=\frac{1}{\sqrt[3]{3^{2}}\Gamma\left . \frac{dw}{dz}\right. _{z=0}=-\frac{1}{\sqrt[3]{3}\Gamma\left(\frac{1}{3}\right)}. \) and \(\left.\frac{dw}{dz}\right|_{z=0}=-\frac{1}{\sqrt[3]{3}\Gamma\left(\frac{1}{3}\right)}.\) Được biết, nghiệm của phương trình vi phân với các điều kiện biên này là hàm Airy

\[w=\textrm{Ai}\left(z\right),\]

cung cấp một phương tiện để kiểm tra bộ tích hợp bằng cách sử dụng

Đầu tiên, chuyển đổi ODE này thành dạng chuẩn bằng cách đặt \(\mathbf{y}=\left[\frac{dw}{dz},w\right]\) . Do đó, phương trình vi phân trở thành and \(t=z\). Thus, the differential equation becomes

\[\begin{split}\frac{d\mathbf{y}}{dt}=\left[\begin{array}{c} ty_{1}\\ y_{0}\end{array}\right] . \end{split}\]

Nói cách khác,

\[\mathbf{f}\left(\mathbf{y},t\right)=\mathbf{A}\left(t\right)\mathbf{y}. \]

Xin lưu ý rằng nếu \(\mathbf{A}\left(t\right)\) đi làm bằng \(\int_{0}^{t}\mathbf{A}\left(\tau\right)\, d\tau\) under matrix multiplication, then this linear differential equation has an exact solution using the matrix exponential:

\[\mathbf{y}\left(t\right)=\exp\left(\int_{0}^{t}\mathbf{A}\left(\tau\right)d\tau\right)\mathbf

Tuy nhiên, trong trường hợp này, \(\mathbf{A}\left(t\right)\) và tích phân của nó thì không .

Phương trình vi phân này có thể được giải bằng hàm. Nó yêu cầu đạo hàm, fprime, khoảng thời gian [t_start, t_end] và vectơ điều kiện ban đầu, y0, làm đối số đầu vào và trả về một đối tượng có trường y là một mảng với các giá trị nghiệm liên tiếp dưới dạng cột. Do đó, các điều kiện ban đầu được đưa ra trong cột đầu ra đầu tiên

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

Như có thể thấy, nó tự động xác định các bước thời gian nếu không được chỉ định khác. Để so sánh nghiệm của hàm thoáng, vectơ thời gian được tạo bởi hàm thoáng

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

Giải pháp với các tham số tiêu chuẩn của nó cho thấy độ lệch lớn đối với chức năng thoáng khí. Để giảm thiểu sai lệch này, có thể sử dụng dung sai tương đối và tuyệt đối

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

Để chỉ định các mốc thời gian do người dùng xác định cho giải pháp của , cung cấp hai khả năng cũng có thể được sử dụng bổ sung. Bằng cách chuyển tùy chọn t_eval cho lệnh gọi hàm trả về nghiệm của các điểm thời gian này của t_eval trong đầu ra của nó

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

Nếu biết ma trận jacobian của hàm, nó có thể được chuyển đến để đạt được kết quả tốt hơn. Tuy nhiên, hãy lưu ý rằng phương thức tích hợp mặc định không hỗ trợ ma trận jacobian và do đó phải chọn một phương thức tích hợp khác. Ví dụ, một trong những phương pháp tích hợp hỗ trợ ma trận jacobian là phương pháp của ví dụ sau

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> np.trapz([1,2,3], x=[8,6,4])
-8.0

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701

>>> from numpy import sqrt, sin, cos, pi
>>> I = sqrt(2/pi)*(18.0/27*sqrt(2)*cos(4.5) - 4.0/27*sqrt(2)*sin(4.5) +
..                 sqrt(2*pi) * special.fresnel(3/sqrt(pi))[0])
>>> I
1.117817938088701