Агуулгыг алгасах

Симпсоны томьёогоор интегралчлах

Бид тодорхой интегралын утгыг тооцоолох гэж байна

$$\int_a ^ b f (x) dx$$

Энд тайлбарласан шийдлийг Томас Симпсон 1743 онд диссертацийнхаа нэгэнд нийтэлсэн.

Симпсоны томьёо

$n$ нь ямар нэг натурал тоо байг. Бид интегралчлалын $[a, b]$ хэрчмийг $2n$ тэнцүү хэсэгт хуваана:

$$x_i = a + i h, ~~ i = 0 \ldots 2n,$$
$$h = \frac {b-a} {2n}.$$

Одоо бид $[x_ {2i-2}, x_ {2i}]$, $i = 1 \ldots n$ хэрчим тус бүр дээр интегралыг тусад нь тооцоолж, дараа нь бүх утгыг нэмнэ.

Тэгэхээр $[x_ {2i-2}, x_ {2i}], i = 1 \ldots n$ дараагийн хэрчмийг авч үзье гэж бодъё. Түүн дээрх $f(x)$ функцийг $(x_ {2i-2}, x_ {2i-1}, x_ {2i})$ гэсэн 3 цэгээр дайрч буй $P(x)$ параболоор сольё. Ийм парабол үргэлж оршин байх ба цор ганц; түүнийг аналитик байдлаар олж болно. Жишээ нь бид Лагранжийн олон гишүүнтийн интерполяц ашиглан түүнийг байгуулж болно. Үлдсэн цорын ганц зүйл бол энэ олон гишүүнтийг интегралчлах явдал юм. Хэрэв та үүнийг ерөнхий $f$ функцийн хувьд хийвэл гайхалтай энгийн илэрхийллийг олж авна:

$$\int_{x_ {2i-2}} ^ {x_ {2i}} f (x) ~dx \approx \int_{x_ {2i-2}} ^ {x_ {2i}} P (x) ~dx = \left(f(x_{2i-2}) + 4f(x_{2i-1})+(f(x_{2i})\right)\frac {h} {3} $$

Эдгээр утгыг бүх хэрчмээр нэмбэл бид эцсийн Симпсоны томьёо-г олж авна:

$$\int_a ^ b f (x) dx \approx \left(f (x_0) + 4 f (x_1) + 2 f (x_2) + 4f(x_3) + 2 f(x_4) + \ldots + 4 f(x_{2N-1}) + f(x_{2N}) \right)\frac {h} {3} $$

Алдаа

Симпсоны томьёогоор интегралыг ойролцоолоход гарах алдаа нь

$$ -\tfrac{1}{90} \left(\tfrac{b-a}{2}\right)^5 f^{(4)}(\xi)$$

энд $\xi$ нь $a$ ба $b$ хоёрын хооронд орших ямар нэг тоо.

Алдаа нь асимптотоор $(b-a)^5$-д пропорциональ. Гэвч дээрх гаргалгаанууд нь $(b-a)^4$-д пропорциональ алдааг санал болгож байгаа. Симпсоны дүрэм нэмэлт эрэмбэ олж авдаг, учир нь интеграл доорх функцийг тооцоолох цэгүүд $[a, b]$ интервалд тэгш хэмтэйгээр тархсан байдаг.

Implementation

Here, $f(x)$ is some user-defined function.

const int N = 1000 * 1000; // number of steps (already multiplied by 2)

double simpson_integration(double a, double b){
    double h = (b - a) / N;
    double s = f(a) + f(b); // a = x_0 and b = x_2n
    for (int i = 1; i <= N - 1; ++i) { // Refer to final Simpson's formula
        double x = a + h * i;
        s += f(x) * ((i & 1) ? 4 : 2);
    }
    s *= h / 3;
    return s;
}

Дасгал бодлогууд