Komputasi Modern Matematika Integral Numerik

Integral Numerik

✍ Objektif :

⊲ Mengenalkan metode Trapezoida ⊲ Mengenalkan metode Simpson

⊲ Mengenalkan metode Composite-Simpson ⊲ Mengenalkan metode Adaptive Quardrature ⊲ Mengenalkan metode Gaussian Quadrature

9.1        Metode Trapezoida

Integral terhadap suatu fungsi, f(x), yang dievaluasi dari a hingga b dapat dinyatakan oleh rumus berikut ini

Z _b

f (x)dx

(9.1)

a

Pendekatan numerik yang paling dasar dalam memecahkan masalah integral adalah metode Trapezoida, yang dirumuskan sebagai berikut

b		h		h3
Za	f (x)dx =		[f (x0) + f (x1)] −		f ′′(ξ)	(9.2)
		2		12

dimana x₀ = a, x₁ = b dan h = b − a. Akan tetapi, suku yang terakhir pada ruas kanan dimana terdapat faktor turunan ke-2, f ^′′, seringkali diabaikan dengan tujuan agar persamaan (9.2) menjadi lebih sederhana.

Za	b	h
	f (x)dx =		[f (x0) + f (x1)]	(9.3)
		2

Akibatnya pendekatan Trapezoida hanya bekerja efektif pada fungsi-fungsi yang turunan kedua-nya bernilai nol (f ^′′ = 0). Gambar (9.1) memperlihatkan prinsip metode trapezoida dalam bentuk grafik. Sementara, script berikut ini dibuat berdasarkan persamaan (9.3).

189

f(x)			f(x)
	f(x1)
	f(x0)

x₀=a                                                             x₁=b                                              x₀=a                           x₁=b

Gambar 9.1: Metode Trapezoida. Gambar sebelah kiri menunjukkan kurva fungsi f (x) dengan batas bawah integral adalah a dan batas atas b. Gambar sebelah kanan menunjukan cara metode Trapesoi-da menghitung integral dengan cara menghitung luas area integrasi, dimana luas area integrasi sama dengan luas trapesium di bawah kurva f (x) dalam batas-batas a dan b. Jika anda perhatikan dengan teliti, ada area kecil dibawah garis kurva dan diatas garis miring yang berada diluar bidang trapesium. Metode Trapesoida tidak menghitung luas area kecil tersebut. Disinilah letak kelemahan metode trape-zoida.

1        clear all

2        clc

3

4        a = ...  %batas bawah integral;

5        b = ...  %batas atas integral;

6

7        x0 = a;

8        x1 = b;

9        h = b-a;

10

11        % -- metode trapezoida --

12        Int_trapezoida = h/2_*(f(x0)+f(x1))

Dengan fungsi eksternal fungsi f(x) adalah

1        function y = f(x)

2        y = ... % rumus fungsi yang di-integralkan;

9.2        Metode Simpson

Metode pendekatan yang lebih baik dibanding metode Trapezoida dalam integral numerik adalah metode Simpson yang diformulasikan sebagai berikut

Z	b		h		h5
		f (x)dx =		[f (x0) + 4f (x1) + f (x2)] −		f 4(ξ)	(9.4)
	a		3		90

dengan x₀ = a, x₂ = b, dan x₁ = a + h dimana h = (b − a)/2. Jika suku terakhir diabaikan, maka

Za	b	h
	f (x)dx =		[f (x0) + 4f (x1) + f (x2)]	(9.5)
		3

Gambar (9.2) memperlihatkan prinsip metode trapezoida dalam bentuk grafik. Sementara, script berikut ini dibuat berdasarkan persamaan (9.5).

9.2.  METODE SIMPSON

191

	f(x)				f(x)
		f(x 2)
		f(x 1)
		f(x0)
					h
				h
x0=a	x1=b		x0=a	x1	x2=b

Gambar 9.2: Metode Simpson. Gambar sebelah kiri menunjukkan kurva fungsi f (x) dengan batas bawah integral adalah a dan batas atas b. Gambar sebelah kanan menunjukan cara metode Simpson menghitung luas area integrasi, dimana area integrasi di bawah kurva f (x) dibagi 2 dalam batas interval a − x₁ dan x₁ − b dengan lebar masing-masing adalah h

1        clc

2        clear all

3

4        a = ... %batas bawah integrasi ;

5        b = ... %batas atas integrasi ;

6

7        x0 = a;

8        x2 = b;

9        h = (b-a)/2;

10     x1 = a + h;

11

12        % -- metode simpson --

13        Int_simpson = h/3_*(f(x0)+4_*f(x1)+f(x2))

Contoh

Metode Trapezoida untuk fungsi f pada interval [0,2] adalah

Z ₂

f (x)dx ≈ f (0) + f (2)

0

dimana x₀ = 0, x₁ = 2 dan h = 2 − 0 = 2. Sedangkan metode Simpson untuk fungsi f pada interval [0,2] adalah

Z0	2	1
	f (x)dx ≈		[f (0) + 4f (1) + f (2)]
		3

dengan x₀ = 0, x₂ = 2, dan x₁ = a + h = 1 dimana h = (b − a)/2 = 1.

Tabel berikut ini memperlihatkan evaluasi integral numerik terhadap beberapa fungsi dalam interval [0,2] beserta solusi exact-nya. Jelas terlihat, metode Simpson lebih baik dibanding Trapezoida. Karena hasil intergral numerik metode Simpson lebih mendekati nilai exact

	x2	x4		√			ex
f (x)			1/(x + 1)		1 + x2	sin x
Nilai exact	2,667	6,400	1,099	2,958		1,416	6,389
Trapezoida	4,000	16,000	1,333	3,326		0,909	8,389
Simpson	2,667	6,667	1,111	2,964		1,425	6,421

9.3        Peran faktor pembagi, n

Kalau diamati lebih teliti, akan kita dapatkan bahwa interval [0,2] telah dibagi 2 pada metode Simpson, sementara pada metode Trapesoida tidak dibagi sama sekali. Sebenarnya dengan membagi interval lebih kecil lagi, maka error -nya akan semakin kecil. Misalnya, banyaknya pembagian interval dinyatakan dengan n

ketika n = 1: Trapesioda

x1

h

h3

Zx0

f (x)dx =

[f (x0) + f (x1)] −

f ′′(ξ)

(9.6)

2

12

ketika n = 2: Simpson

x2

h

h5

Zx0

f (x)dx =

[f (x0) + 4f (x1) + f (x2)] −

f 4(ξ)

(9.7)

3

90

ketika n = 3: Simpson tiga-per-delapan

x3

3h

3h5

Zx0

f (x)dx =

[f (x0) + 3f (x1) + 3f (x2) + f (x3)] −

f 4(ξ)

(9.8)

8

80

ketika n = 4:

x4

2h

8h7

Zx0

f (x)dx =

[7f (x0) + 32f (x1) + 12f (x2) + 32f (x3) + 7f (x4)] −

f 6(ξ)

(9.9)

45

945

9.3.1        Source code metode integrasi

Source code untuk persamaan (9.8) disajikan sebagai berikut

1        clc

2        clear all

3

4        % -- batas integrasi --

5        a = 0;

6        b = 2;

7

8        x0 = a;

9        x3 = b;

10        h = (b-a)/3;

11        x1 = a + h;

12        x2 = a + 2_*h;

13        % ---------------------

14

15        % -- metode simpson 3/8 --

16        Int_38 = 3_*h/8_*(f(x0)+3_*f(x1)+3_*f(x2)+f(x3))

Sementara, source code untuk persamaan (9.9) disajikan sebagai berikut

1        clc

2        clear all

3

4        % -- batas integrasi --

5        a = 0;

6        b = 2;

7

8        x0 = a;

9        x4 = b;

10        h = (b-a)/4;

11        x1 = a + h;

12        x2 = a + 2_*h;

13         x3 = a + 3_*h;

14        % ---------------------

15

16        % -- metode simpson n=4 --

17        Int_n4 = 2_*h/45_*(7_*f(x0)+32_*f(x1)+12_*f(x2)+32_*f(x3)+7_*f(x4))

Perbandingan tingkat akurasi hasil perhitungan seluruh metode integral numerik yang su-

dah dibahas adalah sebagai berikut

	x2	x4		√			ex
f (x)			1/(x + 1)		1 + x2	sin x
Nilai exact	2,667	6,400	1,099	2,958		1,416	6,389
Trapezoida	4,000	16,000	1,333	3,326		0,909	8,389
Simpson n=2	2,667	6,667	1,111	2,964		1,425	6,421
Simpson n=3	2,667	6,519	1,105	2,960		1,420	6,403
Simpson n=4	2,667	6,400	1,099	2,958		1,416	6,389

Keempat bentuk persamaan integral numerik di atas dikenal dengan closed Newton-Cotes formulas. Keterbatasan metode Newton-Cotes terlihat dari jumlah pembagian interval. Di

atas tadi pembagian interval baru sampai pada n = 4.               Bagaimana bila interval evaluasinya

dipersempit supaya solusi numeriknya lebih mendekati solusi exact?  Atau dengan kata lain

n > 4.

9.4        Metode Composite-Simpson

Persamaan (9.9) terlihat lebih rumit dibandingkan persamaan-persamaan sebelumnya. Bisakah

anda bayangkan bentuk formulasi untuk n = 5 atau n = 6 dan seterusnya?  Pasti akan lebih

kompleks dibandingkan persamaan (9.9).

Metode Composite Simpson menawarkan cara mudah menghitung intergal numerik ketika nilai n > 4. Perhatikan contoh berikut, tentukan solusi numerik dari ^R₀⁴ e^xdx. Metode Simpson dengan h = 2 (atau interval evaluasi integral dibagi 2 , n = 2) memberikan hasil

Z		4	2
		exdx ≈		e0 + 4e2 + e4	= 56, 76958
	0		3
Padahal solusi exact dari integral tersebut adalah e4					− e0  =  53, 59815, artinya terdapat er-

ror sebesar 3,17143 yang dinilai masih terlampau besar untuk ditolerir.  Bandingkan dengan

                                                                                                                                     

f(x)

h

x₀=a                                                     x₁      x₂     x₃      x₄      x₅     x₆      x₇   x_n=b

Gambar 9.3: Metode Composite Simpson. Kurva fungsi f (x) dengan batas bawah integral adalah a dan batas atas b. Luas area integrasi dipecah menjadi 8 area kecil dengan lebar masing-masing adalah h.

metode yang sama namun dengan h = 1 (atau interval evaluasi integral dibagi 4 , n = 4)

Z0 4 exdx   =	Z0 2 exdx + Z2 4 exdx
	1			e0 + 4e + e2		1		e2 + 4e3 + e4
≈						+
	3						3
	1			0
					2			3	4
=				e  + 4e + 2e  + 4e  + e
	3
=	53, 86385

Hasil ini memperlihatkan error yang makin kecil, yaitu menjadi 0,26570. Jadi dengan mem-perkecil h, error menjadi semakin kecil dan itu artinya solusi integral numerik semakin mendekati solusi exact. Sekarang kita coba kecilkan lagi nilai h menjadi h = ¹₂ (atau interval evaluasi in-tegral dibagi 8 , n = 8),

Z ₄                                             Z ₁                 Z ₂                 Z ₃                 Z ₄

e^xdx                                 =             e^xdx +         e^xdx +         e^xdx +         e^xdx

0                                                    0                      1                     2                      3

_≈        ¹_{6   e}0 _{+ 4e}1/2 _{+ e  +} ¹_{6   e + 4e}3/2 _{+ e}2   ₊

1	e2 + 4e5/2 + e3   +	1	e3 + 4e7/2 + e4
6		6

₌         ¹ _e0 _{+ 4e}1/2 _{+ 2e + 4e}3/2 _{+ 2e}2 _{+ 4e}5/2 _{+ 2e}3 _{+ 4e}7/2 _{+ e}4 6

=     53, 61622

dan seperti yang sudah kita duga, error -nya semakin kecil menjadi 0,01807.

Prosedur ini dapat digeneralisir menjadi suatu formula sebagai berikut

Z _b

f (x)dx            =

a

=

n/2		x2j
j=1
	Zx2j−2 f (x)dx
X
n/2	3 [f (x2j−2) + 4f (x2j−1) + f (x2j )] −
j=1
X		h

h5
90 f (4)(ξj )		(9.10)

dimana h = (b −a)/n dan x_j = a + jh, untuk j = 1, ..., n/2, dengan x₀ = a dan x_n = b. Formula ini dapat direduksi menjadi

Za

f (x)dx =

3

f (x0) + 2

j=1

f (x2j ) + 4 j=1 f (x2j−1) + f (xn)

−

90 j=1 f (4)

(ξj )    (9.11)

b

h

(n/2)−1

n/2

h5

n/2

X

X

X

Formula ini dikenal sebagai metode Composite Simpson.

9.5        Adaptive Quardrature

Metode composite mensyaratkan luas area integrasi dibagi menjadi sejumlah region dengan jarak interval yang seragam yaitu sebesar nilai h. Akibatnya, bila metode composite diterapkan pada fungsi yang memiliki variasi yang tinggi dan rendah sekaligus, maka interval h yang kecil menjadi kurang efektif, sementara interval h yang besar mengundang error yang besar pula. Metode Adaptive Quadrature muncul untuk mendapatkan langkah yang paling efektif dimana nilai interval h tidak dibuat seragam, melainkan mampu beradaptasi sesuai dengan tingkat variasi kurva fungsinya.

Misalnya kita bermaksud mencari solusi numerik dari integral ^R_a^b f (x)dx dengan toleransi

ǫ > 0. Sebagai langkah awal adalah menerapkan metode Simpson dimana step size h = (b −

a)/2

Z	b		h5
		f (x)dx = S(a, b) −			f (4)(µ)
	a			90

dengan

S(a, b) = ^h₃ [f (a) + 4f (a + h) + f (b)]

Langkah berikutnya adalah men

Za

b

6

f (a) + 4f

a + 2

+ 2f (a + h) + 4f   a +

2

+ f (b)

f (x)dx   =

h

h

3h

4

−

h

(b − a)

f (4)(µ˜)

2

180

(9.12)

(9.13)

9.6        Gaussian Quadrature

Suatu integral dapat ditransformasi kedalam bentuk Gaussian quadrature melalui formulasi berikut

Z	b	Z	1		2		2
	a		−1
	f (x)dx =			f	(b − a)t + (b + a)		(b − a)	dt	(9.14)

dimana perubahan variabel memenuhi

t =

2x − a − b

⇔

x =

1

[(b

−

a)t + a + b]

(9.15)

2

b

−

a

Berikut adalah table polinomial Legendre untuk penyelesaian Gaussian quadrature

Tabel 9.1: Polinomial Legendre untuk n=2,3,4 dan 5
	n	Akar rn,i	Koefisien cn,i
2		0,5773502692	1,0000000000
		-0,5773502692	1,0000000000
3		0,7745966692	0,5555555556
		0,0000000000	0,8888888889
		-0,7745966692	0,5555555556
4		0,8611363116	0,3478548451
		0,3399810436	0,6521451549
		-0,3399810436	0,6521451549
		-0,8611363116	0,3478548451
5		0,9061798459	0,2369268850
		0,5384693101	0,4786286705
		0,0000000000	0,5688888889
		-0,5384693101	0,4786286705
		-0,9061798459	0,2369268850

9.6.1    Contoh

Selesaikan integrasi berikut ini

^Z                                                                                            1,5

e^−x2 dx                                                                                                                                                                             (9.16)

1

(Solusi exact integral diatas adalah: 0.1093643)

jawab:

Pertama, integral tersebut ditransformasikan kedalam Gaussian quadrature melalui persamaan (9.14)

Z	1,5	e−x	2		1	Z	1		−(t+5)2
				dx =				e		dt	(9.17)
									16
	1				4		−1

Kedua, Gaussian quadrature dihitung menggunakan konstanta-konstanta yang tercantum pa-da tabel polinomial Legendre. Untuk n = 2

Z1

1,5

2

1

he(−(0,5773502692+5)

2

2

/16)i = 0, 1094003

e−x

dx ≈

/16) + e(−(−0,5773502692+5)

4

Untuk n = 3

1,5

2

1

2

2

Z1

[(0, 5555555556)e(−(0,7745966692+5)

/16) + (0, 8888888889)e(−(5)

/16)

e−x

dx

≈

4

+

(−(−0,7745966692+5)2

/16)

] = 0, 1093642

(0, 5555555556)e

Sumber : http://supriyanto.fisika.ui.ac.id/laci04/komputasi_matlab_3.pdf