La Design 96: 2028 Kbytes of D:My DocumentsPublikasiMetoda NumerikMetoda Numerik.doc printed on Saturday 12/02/05 08:24
METODA NUMERIK
oleh Ir. Djoko Luknanto, M.Sc., Ph.D. November 2001
Bahan kuliah Metoda Numerik Jurusan Teknik Sipil FT UGM Yogyakarta
PRAKATA
Buku berjudul “Metoda Numerik” ini merupakan bahan kuliah di Jurusan Teknik Sipil FT UGM. Buku ini tidak menjelaskan secara rinci teori-teori numerik secara lengkap, namun hanya membahas teori-teori numerik yang sering digunakan di lapangan. Pembaca yang ingin mengetahui secara lengkap Metoda Numerik diharapkan mencari dari acuan-acuan di luar buku ini. Buku ini lebih merupakan petunjuk praktis bagi mahasiswa S1 maupun praktisi di lapangan. Dalam buku ini prinsip umum teori-teori numerik dijelaskan secara singkat, kemudian aplikasinya dijelaskan. Semoga buku kecil ini berguna, kritik membangun sangatlah diharapkan.
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Yogyakarta, November 2001 Dosen Jurusan Teknik Sipil FT UGM
Ir. Djoko Luknanto, M.Sc., Ph.D. Penyusun
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. ii Jack la Motta
DAFTAR ISI
halaman
PRAKATA .......................................................................................................................... ii DAFTAR ISI...................................................................................................................... iii DAFTAR GAMBAR .........................................................................................................vi
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
1. Error: Asal & Rambatannya....................................................................................... 7 1.1. Pendahuluan ....................................................................................................... 7 1.2. Bilangan Dalam Komputer ............................................................................. 11 1.2.1. Underflow and Overflow ................................................................. 12 1.3. Definisi dan Asal ‘error’ .................................................................................. 12 1.3.1. Angka signifikan (significant digits)............................................... 13 1.3.2. Asal dari ‘error’ .................................................................................. 13 1.4. Rambatan ‘Error’ .............................................................................................. 13 1.4.1. ‘Propagated Error’ pada Perkalian.................................................. 14 1.4.2. ‘Propagated Error’ pada Pembagian............................................... 14 1.4.3. ‘Propagated Error’ pada Penjumlahan dan Pengurangan .......... 14 2. Persamaan Non-Linier.............................................................................................. 16 2.1. Metode Bagi Paruh (Bisection) ....................................................................... 16 2.2. Metode Newton ................................................................................................ 17 2.3. Metode Sekan.................................................................................................... 19 2.4. Akar dari Persamaan Polinomial ................................................................... 20 3. Teori Interpolasi ........................................................................................................ 22 3.1. Metoda Beda Terbagi Newton........................................................................ 22 3.2. Interpolasi dengan tabel beda hingga ........................................................... 24 3.2.1. Beda Maju ........................................................................................... 24 3.2.2. Beda Mundur ..................................................................................... 26 3.3. Lagrange ............................................................................................................ 26 3.4. Beberapa fakta penting dari’beda terbagi’.................................................... 27
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. iii Jack la Motta
DAFTAR ISI
Buku kuliah
4. Integrasi Numeris...................................................................................................... 28 4.1. Rumus trapesium dan Simpson ..................................................................... 28 4.1.1. Rumus trapesium terkoreksi............................................................ 30 4.1.2. Rumus Simpson ................................................................................. 31 4.2. Rumus Newton–Cotes..................................................................................... 33 4.2.1. Rumus Newton-Cotes Tertutup ...................................................... 34 4.2.2. Rumus Newton–Cotes terbuka........................................................ 35 4.3. Kuadratur Gaussian......................................................................................... 36 4.3.1. Kuadratur Gauss-Legendre.............................................................. 37 4.4. Polinomial Orthogonal .................................................................................... 39 4.4.1. Kuadratur Gauss-Laquerre .............................................................. 40 4.4.2. Kuadratur Gauss-Chebysev ............................................................. 42 4.4.3. Kuadratur Gauss-Hermite................................................................ 42
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
5. Sistem Persamaan Linier ......................................................................................... 43 5.1. Eliminasi Gauss ................................................................................................ 43 5.2. Eliminasi Gauss–Jordan................................................................................... 45 5.3. Eliminasi Gauss–Jordan dengan ‘pivot’ maksimum................................... 46 5.3.1. Rekonstruksi pembentukan “scrambled inverse” ........................ 48 5.4. Metoda Iterasi ................................................................................................... 49 5.4.1. Metoda Jacobi ..................................................................................... 49 5.4.2. Metoda Gauss-Seidel......................................................................... 51 6. Matrik .......................................................................................................................... 52 6.1. Notasi dan Konsep-konsep Pendahuluan .................................................... 52 6.2. Determinan dan invers .................................................................................... 55 6.2.1. Menghitung determinan dengan eleminasi segitiga atas ............ 55 6.3. Matrik dan Vektor Eigen ................................................................................. 57 6.3.1. Metode ‘power’ untuk mendapatkan nilai & vektor eigen terbesar. ............................................................................................... 58 7. Persamaan Differensial Biasa ................................................................................. 61 7.1. Metoda Euler..................................................................................................... 63 7.2. Metoda ‘Multi–Step’......................................................................................... 64 7.2.1. Metoda Trapesium ............................................................................ 65 7.3. Metoda Runge-Kutta (RK) .............................................................................. 66 7.3.1. Metoda RK derajat dua ..................................................................... 66 7.3.2. Metoda RK berderajat tiga ............................................................... 67 7.3.3. Metoda RK berderajat empat ........................................................... 68 7.3.3.1. 7.3.3.2. 7.3.3.3.
Metoda Pertama ...................................................................................... 68 Metoda Kedua.......................................................................................... 68 Metoda Ketiga.......................................................................................... 68
7.4. Metoda ‘Predictor-Corrector’ ......................................................................... 69 7.4.1. Algoritma ‘Predictor-Corrector ....................................................... 70
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. iv Jack la Motta
DAFTAR ISI
Buku kuliah
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
DAFTAR PUSTAKA ....................................................................................................... 73
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. v Jack la Motta
DAFTAR GAMBAR
halaman Gambar 1 Teorema Nilai Antara ...................................................................................... 7 Gambar 2 Teorema Nilai Tengah ..................................................................................... 8 Gambar 3 Nilai Tengah Integral ....................................................................................... 9 Gambar 4 Interpretasi Deret Taylor secara geometris................................................. 10 Gambar 5 Metoda Bagi Paruh untuk mencari akar ..................................................... 17 Gambar 6 Metoda Newton untuk mencari akar .......................................................... 18 Gambar 7 Metoda Sekan untuk mencari akar .............................................................. 19 Gambar 8 Konsep integrasi trapesium .......................................................................... 28 Gambar 9 Konsep integrasi Simpson............................................................................. 31 Gambar 10 Fungsi y = w(x) untuk metoda Simpson.................................................... 32 Gambar 11 Cara pertama pemindahan kolom dengan elemen pivot....................... 48 Gambar 12 Cara kedua pemindahan kolom dengan elemen pivot........................... 48 Gambar 13 Cara kedua pemindahan kolom dengan elemen pivot........................... 49 D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Gambar 14 Gaya-gaya yang bekerja pada struktur ..................................................... 60 Gambar 15 Penyelesaian dengan Metoda Euler........................................................... 62
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. vi Jack la Motta
1.
Bab ERROR: ASAL & RAMBATANNYA
1.1. Pendahuluan Teorema 1.1.: ‘Nilai Antara’ (lihat Gambar 1) Jika f(x) suatu fungsi menerus pada x ∈ [a, b] dan m = Infimum f ( x) serta a ≤ x ≤b
M = Supremum f ( x) , maka untuk setiap bilangan ζ pada interval tertutup [m, M] a ≤ x ≤b
paling tidak ada satu titik ξ ∈ [a,b] sehingga f(ξ) = ζ. Khususnya ada dua titik u &
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
ū ∈ [a, b] dimana m = f(u) dan M = f(ū)
y
y = f(x)
ζ
M = f(ū)
m = f(u) a
u
x ξ
ū
b
Gambar 1 Teorema Nilai Antara
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 7 Jack la Motta
Error: Asal & Rambatannya
Buku kuliah
Teorema 1.2.: ‘Nilai Tengah’ (lihat Gambar 2) Jika f(x) menerus pada interval [a,b] serta turunan pertamanya ada dalam interval x ∈ (a,b). Maka paling tidak ada satu titik ξ ∈ (a,b) dimana: f ' (ξ ) =
f (b) − f (a ) b−a
Teorema 1.3: ‘Nilai Tengah Integral’ (lihat Gambar 3) Jika w(x) tidak negatif dan dapat dihitung integralnya pada interval [a,b] dan f(x) menerus pada [a,b], maka
b
b
a
a
∫ w( x) f ( x)dx = f (ξ )∫ w( x)dx untuk satu titik ξ ∈ [a,b]
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
tangennya = f(ξ)
f(ξ)
f(b)
ξ
b
f(a) a
x
Gambar 2 Teorema Nilai Tengah Teorema 1.4.: ‘Deret Taylor’ (lihat Gambar 4) Jika f(x) mempunyai n+1 turunan dan turunannya selalu menerus pada [a,b], dan jika x, x0 ∈ [a,b], maka: f ( x) = Pn ( x) + Rn +1 ( x)
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 8 Jack la Motta
Error: Asal & Rambatannya
Buku kuliah
( x − xo ) n ( n ) x − xo Pn ( x) = f ( xo ) + f ' ( x o ) + ... + f ( xo ) 1! n!
dengan
x
1 Rn +1 ( x) = ∫ ( x − t ) n f ( n +1) (t )dt n! x0 ( x − xo ) = (n + 1)!
n +1
f ( n +1) (ξ )
untuk ξ diantara xo dan x.
y = w(x)
w(xi)
f(xi)
garis referensi
y
a
xi
x
b
Gambar 3 Nilai Tengah Integral
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
x
Bukti:
∫ f ' (t )dt =
f ( x) − f ( xo )
xo
Jadi: x
f ( x) = f ( xo ) + ∫ f ' (t )dt xo
= ... + tf ' (t )]
t=x t = xo
ingat : ∫ udv = uv − ∫ vdu
x
− ∫ tdf ' (t ) xo
x
= ... + ( x − xo ) f ' ( xo ) + x( f ' ( x) − f ' ( xo ) − ∫ tf " (t )dt xo
x
x
xo
xo
= ... + ( x − xo ) f ' ( xo ) + x ∫ f " (t )dt − ∫ tf " (t )dt Akhirnya diperoleh:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 9 Jack la Motta
Error: Asal & Rambatannya
Buku kuliah
x
f ( x) = f ( x0 ) + ( x − xo ) f ' ( xo ) + ∫ ( x − t ) f " (t )dt ... dst. xo
x
= ... +
1 f " (t )d ( x − t ) 2 ... dst. ∫ 2 xo
f ( x) = f ( xo ) +
x − xo f (ξ ) 1!
α y
C p2(x
R3(x) R2(x)
p1(x
R1(x)
y= p0(x
A B f(x0
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
A’ x0
f(x
ξ
C’ x
x
Gambar 4 Interpretasi Deret Taylor secara geometris Secara geometris artinya: CC ' = AA'+ A' C ' tan α 14243 kesalahan pemotongan = AA'+ A' C '×
BC AB
AA' + BC = { { p0 ( x) R1 ( x)=( x− x0 ) f '(ξ )
Jadi kesalahan pemotongan
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 10 Jack la Motta
Error: Asal & Rambatannya
Buku kuliah
n +1
( x − xo ) Rn +1 ( x) = (n + 1)!
f ( n +1) (ξ )
= konstanta × ( x − x o ) = konstanta × Δx
n +1
n +1
Rn+1 (x) disebut sebagai kesalahan pemotongan order n+1 atau O(Δxn+1).
1.2. Bilangan Dalam Komputer Komputer menyajikan bilangan dalam 2 mode yaitu (1) Integer, (2) Floating point. Basis bilangan yang digunakan dalam komputer jarang sekali yang decimal (basis 10). Hampir semua komputer memakai basis 2 (binari) atau variannya seperti basis 8 (octal) dan basis 16 (hexadecimal). Contoh: a. Pada basis 2, semua bilangan terdiri dari 2 angka yaitu 0 dan 1. Jadi (11011.01)2 = 1.24 + 1.23 + 0.22 + 1.21 + 1.20 + 0.2-1 + 1.2-2 = 27. 25
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
b. Pada basis 16, semua bilangan terdiri dari/dinyatakan dengan angka 0, 1, …, 9, A,B, … ,F Jadi (56C.F)16 = 5.1612+ 6.161 + 12.160 + 15.16-1 = 1338.9375 Jika basis bilangan suatu komputer adalah β, maka suatu bilangan non-zero x disimpan didalam bentuk. x = σ (.a1a2a3 … at) β . βe dengan σ = -1 atau + 1, 0 ≤a1≤ β-1, e = integer, dan a (• a1 a2 a3 L at ) = a11 + a 22 + L + tt
β
β
β
dengan σ disebut tanda, e disebut eksponen Æ L ≤ e ≤ U, (•a1a2a3 … at) disebut mantissa, dan • disebut radix. Akurasi dari sajian ‘floating-point’ suatu komputer. Unit pembulatan, δ , suatu komputer adalah suatu bilangan positip terkecil yang mempunyai sifat bahwa 1+δ>1 Nilai nol, ε , suatu komputer adalah suatu bilangan positip terkecil dimana 1+ε>1 Secara praktis ε dan δ dapat dihitung sbb:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 11 Jack la Motta
Error: Asal & Rambatannya
10
Buku kuliah
ε = 1.0 ε = ε/2.0
If (1.0 + ε .GT. 1.0) GOTO 10 δ = ε * 2.0
1.2.1. Underflow and Overflow Jika suatu bilangan tidak mampu direpresentasikan oleh komputer karena e < L atau e > U, maka akan terjadi under/overflow. Jadi setiap bilangan harus berada dalam interval xL ≤ |x| ≤ xU L-1 -t L-1 dengan xL = β dan xU = (1 – β )β Dalam FORTRAN: ♦ Jika suatu hasil hitungan, |x|≥ xU, maka akan terjadi ‘overflow error’ dan program akan berhenti. ♦ Jika suatu hasil hitungan, |x|≥ xL, maka akan terjadi ‘underflow error’ biasanya x nilainya menjadi nol dan hitungan terus berlanjut.
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
1.3. Definisi dan Asal ‘error’ Dalam penyelesaian suatu masalah, dikehendaki jawaban yang sejati, yang disimbolkan sebagai xT, tetapi biasanya jawaban pendekatanlah yang didapat (ini disimbolkan sebagai xA). Error ( x A ) ≡ xT − x A Untuk banyak keperluan, bukan ‘error’ mutlak yang dikehendaki melainkan ‘error’ relatif dari xA yang dibutuhkan: x − xA Rel( x A ) ≡ T , xT ≠ 0 xT 19 Contoh: xT = e = 2.7182818... , x A = = 2.7142857... 7 Jadi: Error ( x A ) = 0.003996... , Rel( x A ) = 0.00147...
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 12 Jack la Motta
Error: Asal & Rambatannya
Buku kuliah
1.3.1. Angka signifikan (significant digits) Nilai xA dikatakan mempunyai m angka signifikan terhadap xT, jika kesalahan (xT-xA) mempunyai nilai ≤ 5 pada (m+1) angka dihitung ke kanan dari angka nonzero didalam xT. Contoh: 1 1 2 3 4 1 = 0. 3 3333... , x A = 0.333 , xT − x A = 0. 0 0 0 3 ⇑ ⇑ 3 karena pada angka ke 4 kesalahannya < 5, maka x A dikatakan mempunyai 3 angka signifikan, sehingga x A = 0.333 .
1) xT =
1
2) xT = 2 3.496 x A = 23.494 ⇑
1 2
3 4 5
xT − x A = 0 0 . 0 0 2 ⇑
karena pada angka ke 5 kesalahannya < 5, maka x A dikatakan mempunyai 4 angka signifikan, sehingga x A = 23.49 . 1
3) xT = 0.0 2138 x A = 0.02144 ⇑
1 2 3
xT − x A = 0.0 0 0 0 6 ⇑
karena pada angka ke 3 kesalahannya < 5, maka x A dikatakan mempunyai 2 angka signifikan, sehingga x A = 0.021 .
1.3.2. Asal dari ‘error’
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
1. Simplifikasi dan asumsi yang digunakan untuk merubah peristiwa alam ke dalam formula matematik. 2. Kesalahan/keteledoran : kesalahan aritmatik dan programming. 3. Ketidakpastian dalam data. 4. Kesalahan mesin. 5. Kesalahan matematis dalam kesalahan pemotongan
1.4. Rambatan ‘Error’ ♦ Ditentukan semua operasi aritmatik digantikan dengan tanda ω. Jadi ω: + - × / dan ŵ: + - × / → versi komputer ♦ Misalkan xA dan yA adalah bilangan yang akan digunakan dan kesalahannya terhadap xT dan xT adalah ε = xT − x A dan η = yT − y A
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 13 Jack la Motta
Error: Asal & Rambatannya
Buku kuliah
♦ Jika dilakukan hitungan xA ŵ yA, maka kesalahannya adalah xT ω yT − x A wˆ y A = ( xT ω yT − x A ω y A ) + ( x A ω y A − x A wˆ y A ) 144424443 144424443 I
II
I = kesalahan karena rambatan (‘propagated error’) II = kesalahan karena ‘rounding’ ataupun ‘choping’
1.4.1. ‘Propagated Error’ pada Perkalian xT yT − x A y A = xT yT − ( χ T − ε )( yT − η ) = xT η + yT ε − εη x y − xA yA η ε εη Rel( x A y A ) ≡ T T = + − xT y T y T xT χ T y T = Rel( x A ) + Rel( y A ) − Rel( x A ).Rel( y A ) Jika Rel( x A ) , Rel( y A ) <<1, maka Rel( x A y A ) = Rel( x A ) + Rel( y A )
1.4.2. ‘Propagated Error’ pada Pembagian xT x T − ε xT x A − − yT y T − η yT y A xA = Rel( ) ≡ xT xT yA yT yT
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
= 1− Rel(
x y − η . xT − x T y T + ε . y T ( xT − ε ) y T = T T x T y T − η . xT ( y T − η ) xT
xA Rel( x A ) − Rel( y A ) )= yA 1 − Rel( y A )
Jika Rel( y A ) <<1, maka Rel(
xA ) = Rel( x A ) − Rel( y A ) yA
Tampak bahwa pada perkalian dan pembagian ‘kesalahan relatif’ tidak membesar secara cepat.
1.4.3. ‘Propagated Error’ pada Penjumlahan dan Pengurangan ( xT ± y T ) − ( x A ± y A ) = ( xT − x A ) ± ( y T − y A ) = ε ± η Error ( x A ± y A ) = Error ( x A ) ± Error ( y A )
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 14 Jack la Motta
Error: Asal & Rambatannya
Buku kuliah
Penjabaran ini tampak logis tetapi tidak baik karena tidak dinyatakan dalam ‘kesalahan relatif’ 22 Contoh: xT = π , x A = 3.1416, yT = , y A = 3,1429 7 xT − x A = −7.35 × 10 −6 Rel( x A ) = −2.34 × 10 −6 yT − y A = −4.29 × 10 −5 Rel( y A ) = −1.36 × 10 −5 ( xT − yT ) − ( x A − y A ) = −0.0012645 − (−0.0013)
= 3.55 × 10 −5 Rel( x A − y A ) = −0.028 Jadi meskipun kesalahan pada ( x A − y A ) adalah kecil, tetapi ‘kesalahan relatif’nya
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
cukup besar melebihi Rel(xA) ataupun Rel(yA).
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 15 Jack la Motta
2.
Bab PERSAMAAN NON-LINIER
Di dalam matematika aplikasi pencarian akar persamaan f(x)=0 sering dijumpai. Biasanya jawaban analitis dari persamaan diatas tidak ada, sehingga harus dicari jawaban numeriknya yang biasa dilaksanakan dengan metode iterasi.
2.1. Metode Bagi Paruh (Bisection) Jika terdapat suatu f(x) yang menerus ∈ [a,b] dan f(a)f(b) < 0, maka menurut Teorema 1.1 paling tidak f(x) mempunyai satu akar f(x) mempunyai satu akar ∈ [a,b].
♦ Algoritma Bisect(f, a, b, akar, ε) 1. Hitung c := (a+b)/2 2. Jika b – c ≤ ε , maka akar:= c, dan ‘ exit’ D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
3. Jika {tanda f(b)•tanda f(c)} ≤ 0, maka a := c , jika tidak b := c 4. Kembali ke langkah nomor 1.
♦ Definisi Suatu deret hasil suatu iterasi {xn|n≥0} dikatakan menuju ke titik dengan derajat p ≥1 , jika
α − x n +1 ≤ c α − x n
p
n≥0
untuk beberapa nilai c>0. Jika p=1, deretnya disebut menuju ke titik secara linier. Pada kasus ini diperlukan nilai c<1; c disebut laju linier dari xn menuju . Ada beberapa metode yang membutuhkan definisi yang agak berbeda dengan diatas yaitu α − x n +1 ≤ c n α − x0 n ≥ 0
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 16 Jack la Motta
Persamaan Non-Linier
Buku kuliah
♦ Tingkat kelajuan metode bagi paruh dinyatakan dalam 1 α − c n ≤ ( ) n (b − a) 2 y
c = (a+b)/2
c = (a+b)/2
nilai awal
abaru
f(a)<0 y=f(x)
c = (a+b)/2
metode ini diulang-ulang sampai abs (c-b) < ε
abaru
α
f(b)>0
bbaru b nilai
x
akar sesungguhnya yang akan dicari
Gambar 5 Metoda Bagi Paruh untuk mencari akar
2.2. Metode Newton
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Deret Taylor:
1 f ( x) = f ( x n ) + ( x − x n ) f ' ( x n ) + ( x − x n ) 2 f " ( x n ) + ... 2 atau menurut Teorema 1.4 1 f ( x) = f ( x n ) + ( x − x n ) f ' ( x n ) + ( x − x n ) 2 f " (ξ ) 2 dengan ξ diantara xn dan x. Jika akar dari f(x), salah satunya adalah α, maka f (x = α ) = 0 1 f ( x) = f ( x n ) + (α − x n ) f ' ( x n ) + (α − x n ) 2 f " (ξ ) = 0 2 jadi
α = xn −
f ( xn ) 1 f " (ξ ) − (α − x n ) 2 f ' ( xn ) f ' ( xn ) 2
maka α dapat didekati dengan Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 17 Jack la Motta
Persamaan Non-Linier
Buku kuliah
x n +1 = x n −
dengan ‘errornya’
α − x n +1 = −
f ( xn ) f ' ( xn )
n≥0
f " (ξ ) (α − x n ) 2 2 f ' ( xn )
n≥0
Untuk nilai n→∞, ξ→α dan xn→α, jadi f " (α ) α − x n +1 = − (α − x n ) 2 2 f ' (α )
= konstanta × (α − x n ) 2 Sehingga metode Newton dikatakan mempunyai derajat kelajuan = 2 y y=f(x) grs singgung
akar sesungguhnya
α
nilai awal x1
x0
x
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Gambar 6 Metoda Newton untuk mencari akar
♦ Algoritma Newton (f, df, x0, ε, akar, itmax, ierr) 1. Keterangan : df adalah f’(x), itmax adalah iterasi maximum, ierr adalah ‘error flag’ 2. noiter:=1 3. penyebut:=df(x0) 4. jika penyebut = 0 maka ierr:=2, dan ‘exit’ 5. x1:= x0 - f(x0)/penyebut 6. jika |x1 – x0|≤ε, maka ierr:= 0, akar:= x1, dan ‘exit’ 7. jika noiter = itmax maka ierr:= 1, dan ’exit’ 8. noiter:= noiter +1, x0:=x1, dan ulangi langkah 3.
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 18 Jack la Motta
Persamaan Non-Linier
Buku kuliah
2.3. Metode Sekan Dengan menggunakan sifat segitiga sebangun diperoleh
BD CD = BA CE f ( x1 ) − f ( x0 ) f ( x1 ) − 0 = x1 − x0 x1 − x 2 x1 − x 0 atau Jadi x 2 = x1 − f ( x1 ) f ( x1 ) − f ( x0 ) x n − x n −1 x n +1 = x n − f ( x n ) f ( x n ) − f ( x n −1 ) y
n ≥1
y=f(x) D akar sesungguhnya
nilai awal x0
x2
x3 E
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
A
C x α x1 nilai awal B
Gambar 7 Metoda Sekan untuk mencari akar Metode sekan dapat dijabarkan dari metode Newton dimana f ' ( xn ) =
f ( x n ) − f ( x n −1 ) x n − x n −1
Derajat Konvergensi:
• untuk metode Newton p = 2 • untuk metode sekan p = ½(1+√5) = 1,618 • untuk metode ‘bisection’ p = 1
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 19 Jack la Motta
Persamaan Non-Linier
Buku kuliah
2.4. Akar dari Persamaan Polinomial p(x) = a0 + x(a1 + x(a2+…+x(an-1 + anx)…)
(A)
p(x) = a0 + a1x + a2x2+…+ anxn
(B)
atau
Pada Pers.(A) terdapat n perkalian & pertambahan, sedangkan dalam Pers.(B) terdapat: (2n–1) perkalian & pertambahan. Oleh karena itu dalam pemrograman komputer lebih disukai bentuk dalam Pers.(A), karena lebih efisien. Pers. (A) jika ditulis dalam FORTRAN menjadi
10
p = a(n) do 10 i = n,1,-1 p = p*x + a(i–1)
Untuk menghitung akar dari persamaan p(x) = 0 akan digunakan Metoda Newton. Untuk keperluan itu polinomial p(x) akan dimodifikasi sebagai berikut Disyaratkan: bn = an
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
bk = ak + zbk+1 , k = n–1, n–2,… , 0 Dari syarat ini p(x) dapat ditulis sebagai p(x) = b0 + (x– z)q(x) dengan q(x) = b1 + b2x +…+ bnx n-1 sehingga p’(x) = (x– z)q’(x) + q(x) → p’(z) = q(z)
♦ Algoritma: Polynew (a, n, x, ε, itmax, akar, b, ier) 1. Keterangan: a adalah vektor coef. dengan dimensi n, itmax adalah iterasi maksimum, b adalah vektor coef. dari polinomial yang baru, ier adalah indikator adanya error. 2. noiter: = 1 3. x: = x0 , bn := c := an 4. Untuk k = n-1, …, 1, bk := ak+ zbk+1, c := bk + zc 5. b0:= a0 + zb1 6. Jika c = 0, ier := 2, dan ‘exit 7. x1:= x0 – b0/c 8. Jika ⏐x1 – x0⏐≤ε , ier :=0 ,akar:= x1, dan ‘exit’ 9. Jika noiter:= itmax, ier:= 1, dan ‘exit’ Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 20 Jack la Motta
Persamaan Non-Linier
Buku kuliah
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
10. noiter:= noiter + 1, x0:= x1 , ulangi langkah ketiga.
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 21 Jack la Motta
3.
Bab TEORI INTERPOLASI
Jika kita mempunyai satu set data: (x0, y0) , (x1, y1), …, (xn, yn) maka dalam bab ini akan di jelaskan bagaimana harus mencarti polinomial yang melalui, data di atas. Jika polinomial ini ditulis sebagai: p(x) = a0 + a1x +… + anxn maka jika data diatas disubstitusikan akan didapat (n+1) persamaan dengan (n+1) variabel tidak diketahuinya yaitu: a 0 + a1 x 0 + L + a n x 0n = y 0
M a0
+ a1 x n
+ L + an x
n n
=
M yn
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Persamaan diatas jika diselesaikan akan menghasilkan a0, …, an sehingga polinomial p(x) dapat dicari .
3.1. Metoda Beda Terbagi Newton Notasi yang digunakan: f [xo ,x1 , ... ,x n ] =
f [x1 , ...,x n ] − f [x1 , ...,x n ] x n -x0
Contoh Order 0:
ƒ[x0] =ƒ[xn] f [x1 ] − f [x0 ] Order 1: f [x0 , x1 ] = x1 − x0 f [x1 , x 2 ] − f [x0 , x1 ] Order 2: f [x0 , x1 , x 2 ] = x 2 − x0 f [x1 , x 2 , x3 ] − f [x 0 , x1 ] Order 3: f [x0 , x1 , x 2 , x3 ] = x 2 − x0
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 22 Jack la Motta
Teori Interpolasi
Buku kuliah
Rumus beda terbagi Newton: pn(x) = ƒ[x0]+ (x3-x0)ƒ[x0,x1]+(x3-x0)(x-x1)ƒ[x0,x1]+… +(x–x0)…(x–xn-1)ƒ[ x0,x1,…,xn ] Contoh: Kita buat tabel beda terbagi berdasarkan polinomial ƒ(x) = x3 – 2x2 + 7x – 5 i 0 1 2 3
xi 0.0 1.0 3.0 4.0
Keterangan:
1 − (− 5) =6 1− 0 55 − 25 C= = 30 4−3 30 − 12 E= =6 4 −1 A=
f[xi] -5.0 1.0 25.0 55.0
f[xi, xi+1] 6 12 30
f2[ ] 2 6
f3[ ] 1
25 − 1 = 12 3 −1 B− A D= =2 3−0 E−D F= =1 4−0 B=
Contoh hitungan pn(x=0.5) = ?
• p1(x) = -5 + (x-0)6 = 6x – 5 ∴ p1 (0.5) = -2
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
• p2(x) = -5 + (x-0)6 + (x-0)(x–1)2 = 2x2 + 4x – 5 ∴ p2 (0.5) = -2,5
• p3(x) = -5 + (x-0)6 + (x-0)(x–1)2 +(x-0)(x-1)(x-3)1 = x3 – 2x2 + 7x – 5 ∴ p2 (0.5) = -15/8 = - 1.875
♦ Algoritma metoda beda terbagi Newton Divdif (d, x , n) 1. Keterangan: d dan x adalah vektor f(xi) dan xi = 0,1,…,n. Pada saat ‘exit’ di akan terisi oleh nilai f [x0,…,xi] 2. Kerjakan s/d langkah 4 untuk i = 1, 2 ,… ,n 3. Kerjakan sampai dengan langkah 4 untuk j = n, n-1, i 4. dj := (dj -dj-1) /(x-xj-1) 5. ’exit’
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 23 Jack la Motta
Teori Interpolasi
Buku kuliah
Interp(d, x, t, p) 1. Keterangan: Pada awalnya d dan x adalah vektor dari ƒ[x0,…,xi] dan xi, i = 0, 1 , …, n. Pada saat ‘exit’ p akan berisi pn(t). 2. p := dn 3. Kerjakan s/d langkah 4 untuk i = n-1 , n--2, …, 0 4. p := di + (t– xi)p 5. ‘ exit’
3.2. Interpolasi dengan tabel beda hingga 3.2.1. Beda Maju Notasi: Δƒ(xi) = ƒ(xi+1) - ƒ(xi) dengan xi = x0 + ih, i = 0, 1, 2, 3, … Untuk r ≥ 0, Δr+1ƒ(z) = Δrƒ(z+h) - Δr ƒ(z) Δrƒ(z) disebut ‘ beda maju order r ,‘ Δ disebut ‘operator beda maju ‘
Δ0ƒ(x) Δƒ(x)
= ƒ(x) = Δ0ƒ(z+h) - Δ0ƒ(z) = ƒ(x+h) - ƒ(x) 2 Δ ƒ(x) = Δƒ(x+h) - Δƒ(x) Contoh hitungan : Kita gunakan polinomial x3 – 2x2 + 7x – 5 dengan h = 1,0
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Contoh:
i 0 1 2 3 4
xi 0.0 1.0 2.0 3.0 4.0
f(xi) -5.0 1.0 9.0 25.0 55.0
Δƒ 6 8 16 30
Δ2ƒ 2 8 14
Δ3ƒ 6 6
Δ4ƒ 0
Korelasi antara ‘beda maju’ dengan ‘ beda terbagi’ f [x0 , x1 ] =
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
f [x1 ] − f [x0 ] f ( x 0 + h) − f ( x 0 ) Δf ( x 0 ) = = x1 − x 0 h h
hal. 24 Jack la Motta
Teori Interpolasi
Buku kuliah
f ( x 2 ) − f ( x1 ) f ( x1 ) − f ( x0 ) − x 2 − x1 x1 − x0 Δ2 f ( x0 ) f [x0 , x1 , x 2 ] = = x 2 − x0 2h 2 Secara umum:
f [x0 , x1 ,..., x n ] =
Δn f ( x 0 ) n! h n
Akan dijabarkan rumus interpolasi ‘beda maju’ dari rumus interpolasi ‘beda x − x0 terbagi’ Newton. Didefinisikan α = yang menunjukkan letak titik x h terhadap x0. Jadi misalnya α = 1.6, maka x terletak pada jarak 6/10 dari x1 ke arah x2 . Diinginkan rumus untuk: (x – x0) (x – x1) … (x – xk) dinyatakan dalam α x – xj = x0 + αh – (x0 + jh) = (α-j)h Jadi (x – x0) (x – x1) … (x – xk) = α(α -1)… (α -k)hk+1 sehingga 2 n Δf 0 2 Δ f0 n Δ f0 + α (α − 1)h + ... + α (α − 1)...(α − n + 1)h p n ( x) = f 0 + αh h n! h n 2! h 2
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Jika didefinisikan koefisien binomial sbb: ⎛ α ⎞ α (α − 1)...(α − k + 1) ⎜⎜ ⎟⎟ = , k>0 dan k! ⎝k⎠
⎛α ⎞ ⎜⎜ ⎟⎟ = 1 ⎝0⎠ maka didapat rumus interpolasi ‘beda maju’ sbb: n ⎛α ⎞ x − x0 p n ( x ) = ∑ ⎜⎜ ⎟⎟Δ j f ( x 0 ) dengan α = h j=0 ⎝ j ⎠ Contoh hitungan: p(x=1.5) = ? x − x 0 1.5 − 1.0 α= = = 1.5 h 1.0 1) p1(x) = ƒ(x0) + αΔƒ(x0) = -5 + 1.5 (6) = 4 2) p2(x) = ƒ(x0) + αΔƒ(x0) + α(α-1)∆2f(x0)/2! = -5 + 1.5 (6) + 1.5 (0.5)2/2! = 4.75
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 25 Jack la Motta
Teori Interpolasi
Buku kuliah
3.2.2. Beda Mundur ∇ƒ(z) = ƒ(z) - ƒ(z-h) ∇r+1ƒ(z) = ∇rƒ(z) - ∇rƒ(z-h)
Notasi:
r≥1
Rumus interpolasinya n x − x ⎛−1−α ⎞ ⎛ j −1−α ⎞ j ⎟⎟ = 1 ⎟⎟∇ f ( x0 ) dengan α = 0 , ⎜⎜ p n ( x) = ∑ ⎜⎜ j h j =0 ⎝ ⎝ 0 ⎠ ⎠ i -4 -3 -2 -1 0
xi 0.0 1.0 2.0 3.0 4.0
f(xi) -5.0 1.0 9.0 25.0 55.0
∇ƒ 6 8 16 30
∇ 2ƒ 2 8 14
∇ 3ƒ 6 6
∇ 4ƒ 0
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Contoh hitungan: p(x=3.5) = ? x − x − 3.5 + 4.0 α= 0 = = 0.5 h 1.0 1) p1(x) = ƒ(x0) + (-α)∇ƒ(x0) = 55 + (-0.5) 30 = 40 2) p2(x) = p1(x) + (-α)(-α+1)∇2f(x0)/2! = 40 + (-0.5)(0.5)14/2! = 38.25 3) p3(x) = p2(x) + (-α)(-α+1) (-α+2)∇3f(x0)/3! = 38.25 + (-0.5)(0.5)(1.5)6/3! = 37.875
3.3. Lagrange Polinomial Lagrange dibentuk dengan fomulasi berikut: n
p n ( x ) = ∑ Li ( x ) f (xi ) i =0 n
Li ( x ) = ∑ j =0 j ≠i
x − xj xi − x j
i = 0,1,..., n
Contoh:
x − x1 x − x0 f ( x0 ) + f ( x1 ) x0 − x1 x − x01 (x − x1 )(x − x2 ) f (x ) + (x − x0 )(x − x2 ) f (x ) + (x − x0 )(x − x1 ) f (x ) p2 ( x ) = (x0 − x1 )(x0 − x2 ) 0 (x1 − x0 )(x1 − x2 ) 1 (x2 − x0 )(x2 − x12 ) 2 pi ( x ) =
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 26 Jack la Motta
Teori Interpolasi
Buku kuliah
Contoh: hitung p2 ( x ) yang melalui titik-titik (0,15); (1,1); (3,25)
( ( x − x1 )( x − x2 ) x − 1)( x − 3) x 2 − 4 x + 3 = = L0 ( x ) = (x0 − x1 )(x0 − x2 ) (0 − 1)(0 − 3) 3 2 (x − x0 )(x − x2 ) = (x − 0)(x − 3) = x − 3x L1 ( x ) = (x1 − x0 )(x1 − x2 ) (1 − 0)(1 − 3) −2 2 (x − x01 )(x − x1 ) = (x − 0)(x − 1) = x − x L2 ( x ) = (x2 − x0 )(x2 − x1 ) (3 − 0)(3 − 1) 6 Jadi p 2 (x ) = L0 ( x ) × (− 5) + L1 ( x ) × (1) + L2 ( x ) × (25) = 2 x 2 + 4 x − 5
3.4.
Beberapa fakta penting dari’beda terbagi’
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
1.
f (m ) (ξ ) untuk ξ ∈ X {x0 , x1 ,..., x n } m! artinya interval terkecil dimana x 0 , x1 ,..., x m tercakup! Contoh: f (0 ) (ξ ) f [x0 ] = = f ( x0 ) 0! f (x1 ) − f ( x0 ) f [x0 , x1 ] = = f ' (ξ ) ξ ∈ [x0 , x1 ] x1 − x0 1 f [x0 , x1 , x2 ] = f ' (ξ ) ξ ∈ X {x0 , x1 , x 2 } 2 f [x0 , x1 ,....., xm ] =
dimana
X {x0 ,..., x m }
2. Jika f ( x ) adalah polynomial derajad m, maka ⎧polinomial derajad m − n − 1 n < m − 1 ⎪ f [x 0 ,..., x n , x ] = ⎨a m n = m −1 ⎪0 n > m −1 ⎩ dengan f ( x ) = a m x m + a m −1 x m −1 + ... + a1 x + a0
3. Kesalahan dalam interpolasi
(x − x0 )(x − x1 )...(x − xn ) (n+1) (ξ x ) f (n + 1)! ξ x ∈ Η{x0 ,...., x n , x}
f (x ) − p n (x ) = dengan 4.
d f [x 0 ,..., x n , x ] = f [x 0 , x,..., x n , x, x ] dx
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 27 Jack la Motta
Bab
4. INTEGRASI NUMERIS
4.1. Rumus trapesium dan Simpson Pada bab ini akan dibicarakan cara menghitung integral secara numeris dari b
∫
I(f) =
f(x)dx
a
dimana [a,b] berhingga y
y=f(x) (b,f(b))
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
y =p1(x) (a,f(a)) x a
b Gambar 8 Konsep integrasi trapesium
Rumus trapesium pada dasarnya adalah mendekati f(x) dengan garis lurus yang melalui (a,f(a)) dan (b,f(b)) b−a [ f (a) − f (b)] I1 ( f ) = 2
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 28 Jack la Motta
Integrasi Numeris
Buku kuliah
Error: x −a ⎧x − b ⎫ f (x) − p1(x) = f (x) − ⎨ f (a) + f (b)⎬ b −a ⎩a −b ⎭ = (x − a)(x −b) f [a,b, x]
ingat definisi ‘beda terbagi f[a,b,x]. Jadi ‘error’: b
1 E1( f ) = ∫ f (x)dx− (b − a)[f (a) + f (b)] 2 a b
= ∫(x − a)(x −b) f [a,b, x]dx a
dengan harga tengah integral, didapat: b
E1( f ) = f [a, b,ξ]∫ (x −a)(x −b)dx a ≤ξ ≤ b a
⎛1 ⎞⎛ 1 ⎞ = ⎜ f "(η)⎟ ⎜− (b −a)3 ⎟ ⎝2 ⎠⎝ 6 ⎠ 3 (b −a) =− f "(η) 12
η∈[a, b]
Jika interval [a,b] dibagi menjadi n pias sehingga untuk n≥1, h = (b-a)/n dan xj = a + jh, j = 0,1,…,n, didapat: b
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
I(f) =
∫
f(x)dx =
xj
∑ ∫
f(x)dx
j =1 x j −1
a
=
n
n
∑
j =1
⎧h ⎫ h3 − − f ( x ) f ( x ) f " (η j ) ⎬ ⎨ j −1 j 12 ⎩2 ⎭
[
]
dengan xj-1≤nj≤xj. Sehingga integralnya dapat didekati dengan 1 ⎤ ⎡1 I n ( f ) = h ⎢ f 0 + f1 + ... + f n =1 + f n ⎥ 2 ⎦ ⎣2 Kesalahan In(f) terhadap I(f) adalah En ( f ) = I ( f ) − En ( f ) n
= ∑− j =1
=−
n ≥1
h3 f " (η j ) 12
⎤ h3n ⎡ 1 n ⎢ ∑ f " (η j )⎥ 12 ⎣ n j =1 ⎦
Perlu diingat bahwa
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 29 Jack la Motta
Integrasi Numeris
Buku kuliah
Min f " ( x) ≤
a ≤ x ≤b
1 n f " ( x) ∑ f " (η j ) ≤ Max a ≤ x ≤b n j =1
karena f”(x) menerus pada a ≤ x ≤ b, maka h3n En ( f ) = − f " (η ) η ∈ [a, b] 12 h 2 (b − a) f " (η ) =− 12 (b − a ) 3 =− f " (η 12n 2 ~ Estimasi kesalahan asimtotis E n ( f )
Limit n →∞
En ( f ) = Limit n →∞ h2
⎡ 1 n ⎤ ⎢− ∑ f " (η j )h⎥ ⎣ 12 j =1 ⎦
=−
⎡ n ⎤ 1 Limit ⎢∑ f " (η j )h ⎥ 12 n→∞ ⎣ j =1 ⎦
=−
1 f " ( x)dx 12 ∫a
b
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
h2 ~ maka En ( f ) ≡ − { f ' (b) − f ' (a )} 12 Definisi: ~ Jika En(f) adalah kesalahan eksak, sedangkan En ( f ) adalah estimasi darinya, maka ~ En ( f ) disebut estimasi kesalahan asimtotis dari En(f) jika: ~ ~ En ( f ) En ( f ) − En ( f ) = 1 atau Limit =0 Limit n →∞ E ( f ) n →∞ En ( f ) n
4.1.1. Rumus trapesium terkoreksi ~ Dengan menggunakan En ( f ) , rumus trapesium dapat ditingkatkan menjadi: ~ CTn ( f ) = I n ( f ) + E n ( f )
1 ⎤ h2 ⎡1 = h ⎢ f 0 + f1 + ... + f n =1 + f n ⎥ − { f ' (b) − f ' (a )} 2 ⎦ 12 ⎣2
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 30 Jack la Motta
Integrasi Numeris
Buku kuliah
4.1.2. Rumus Simpson y
(c,f(c))
y=f(x)
(b,f(b)) y =p2(x)
(a,f(a)) x a
c=(a+b)/2
b
Gambar 9 Konsep integrasi Simpson Dalam metode Simpson fungsi f(x) didekati dengan p2(x) yang melalui 3 titik (a,f(a)), (c,f(c)) dan (b,f(b)) dimana c = (a+b)/2. b ⎡ ( x − c )( x − b) ⎤ ( x − a )( x − b) ( x − a )( x − c) I2( f ) = ∫ ⎢ f (a ) + f (c ) + f (b)⎥ dx (a − c )(a − b) (c − a )(c − b) (b − a )(b − c) ⎦ a⎣
=
h [ f (a) + 4 f (c) + f (b)] dengan h = b − a 3 2
Kesalahannya: E2 ( f ) = I ( f ) − I 2 ( f ) D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
b
= ∫ ( x − a )( x − b)( x − c) f [a, b, c, x]dx a
Harga tengah integral tidak dapat digunakan karena (x-a)(x-c)(x-b) berganti tanda pada x = c. x
Didefinisikan: w( x) = ∫ (t − a)(t − c)(t − b)dt a
Beberapa fakta mengenai w(x): ( a − b) 4 , w( x) > 0 untuk a < x < b 64 w' ( x) = ( x − a )( x − c)( x − b)
w(a ) = w(b) = 0, w(c) =
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 31 Jack la Motta
Integrasi Numeris
Buku kuliah
(a − b )4
y
64
y=w(x)
a
c = ( a + b) / 2
b
x
Gambar 10 Fungsi y = w(x) untuk metoda Simpson Jadi E2(f) dapat ditulis sebagai: b
E 2 ( f ) = ∫ w' ( x) f [ a, b, c, x]dx a
b
= w( x) f [a, b, c, x]]a − ∫ w( x) b
a
d f [ a, b, c, x] dx
b
E 2 ( f ) = − ∫ w( x) f [a, b, c, x, x]dx a
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
b
= − f [a, b, c, ξ , ξ ]∫ w( x)dx
ξ ∈ [ a, b]
a
=−
f ( 4) (η ) ⎡ 4 5 ⎤ h 24 ⎢⎣15 ⎥⎦
h 5 ( 4) f (η ) η ∈ [a, b] 90 Jika interval [a,b] dibagi menjadi n pias, n ≥ 2, h = (b-a)/n, xj=a+jh untuk =−
j = 1,2,3,…,n, sehingga n / 2 x2 j
I( f ) = ∑
∫ f ( x)dx
n = genap
j =1 x2 j − 2
⎧h ⎫ h 5 ( 4) f (η j )⎬ = ∑ ⎨ f 2 j − 2 + 4 f 2 j −1 + f 2 j − 90 j =1 ⎩ 3 ⎭ dengan x2j-2 ≤ nj ≤ x2j n/2
[
]
Rumus Simpson:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 32 Jack la Motta
Integrasi Numeris
Buku kuliah
h [ f 0 + 4 f1 + 2 f 2 + 4 f 3 + 2 f 4 + ... + 2 f n−2 + 4 f n−1 + f n ] 3 Kesalahan estimasi: h 5 (n / 2) 2 n / 2 ( 4) En ( f ) = I ( f ) − I n ( f ) = − ∑ f (η j ) 90 n j =1 In ( f ) =
h 4 (b − a) ( 4) =− f (η ) η ∈ [a, b] 180 (b − a ) 5 ( 4 ) =− f (η ) 180n 4 h4 ~ Estimasi kesalahan asimtotis: E n ( f ) = − f ( 3) (b) − f ( 3) (a) 180
[
]
4.2. Rumus Newton–Cotes Rumus trapesium dan Simpson sebetulnya merupakan dua buah rumus pertama dari rumus Newton-Cotes. Untuk n≥1, h = (b-a)/n, xj=a+jh untuk j = 0,1,2,3,…,n. Didefinisikan In(f) dengan mengganti f(x) dengan polinomial pn(x) pada titik-titik x0, x1,…, xn: b
b
a
a
I ( f ) = ∫ f ( x)dx =& ∫ p n ( x)dx
Dengan interpolasi Lagrange untuk pn(x), maka b
n
n
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
I n ( f ) = ∫ ∑ l j ,n ( x) f ( x j )dx = ∑ w j ,n ( x) f ( x j ) a j =0
j =0
b
dengan w j ,n = ∫ l j ,n ( x)dx
j = 0,1,..., n
a
Untuk nilai n = 1 dan 2 telah disajikan sebagai rumus trapesium dan Simpson. Sekarang untuk n = 3, contoh untuk menghitung w0 adalah: b ( x − x1 )( x − x 2 )( x − x3 ) w0 = ∫ dx ( x0 − x1 )( x0 − x 2 )( x 0 − x3 ) a Jika x = x0+μh, 0 ≤ μ ≤ 3, maka:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 33 Jack la Motta
Integrasi Numeris
Buku kuliah
x3
w0 = ∫ − x0
1 ( x − x1 )( x − x 2 )( x − x3 )dx 6h 3 3
=−
1 ( μ − 1)h ( μ − 2)h ( μ − 3)h hdμ 6h 3 ∫0 3
3h h = − ∫ ( μ − 1)( μ − 2)( μ − 3)dμ = 60 8
Jika w1, w2, w3 dihitung dengan cara di atas, akhirnya akan didapat untuk n = 3 3h I 3 ( f ) = [ f ( x0 ) + 3 f ( x1 ) + 3 f ( x 2 ) + f ( x3 )] 8 Kesalahan pada In(f) dinyatakan sebagai berikut: a) Untuk n genap: E n ( f ) = C n h n +3 f ( n + 2 ) (η )
η ∈ [ a, b]
n
dengan C n =
1 μ 2 ( μ − 1)...( μ − n)dμ (n + 2)! ∫0
b) Untuk n gasal: E n ( f ) = C n h n + 2 f ( n +1) (η )
η ∈ [ a, b]
n
dengan C n =
1 μ 2 ( μ − 1)...( μ − n)dμ ∫ (n + 1)! 0
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
4.2.1. Rumus Newton-Cotes Tertutup n = 1, rumus trapesium b h h3 ∫a f ( x)dx = 2 [ f (a) + f (b)] − 12 f " (ξ ) n = 2, rumus Simpson b 5 h⎡ a+b ⎤ h ( 4) f ( x ) dx = f ( a ) + 4 f ( ) + f ( b ) − ∫a ⎥ 90 f (ξ ) 3 ⎢⎣ 2 ⎦ n=3 b
∫ a
3h 3h 5 ( 4 ) f ( x)dx = [ f (a) + 3 f (a + h) + 3 f (b − h) + f (b)] − f (ξ ) 8 80
n = 4, rumus Boole b 7 2h ⎡ a+b ⎤ 8h ( 6 ) ∫a f ( x)dx = 45 ⎢⎣7 f (a) + 32 f (a + h) + 12 f ( 2 ) + 32 f (b − h) + 7 f (b)⎥⎦ − 945 f (ξ )
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 34 Jack la Motta
Integrasi Numeris
Buku kuliah
~ Definisi: Integrasi numerik I ( f ) yang mendekatri I(f) disebut mempunyai ~ derajat ketepatan m jika: (a) I ( f ) = I ( f ) untuk semua polinomial f(x) derajat ~ ≤ m , (b) I ( f ) ≠ I ( f ) untuk beberapa polinomial f(x) derajat m + 1 Contoh: Pada rumus Newton-Cotes untuk n = 1, 3 dikatakan mempunyai derajat ketepatan m = 1, 3. Sedangkan untuk n = 2, 4 mempunyai derajat ketepatan m = n + 1 = 3, 5. Tampak bahwa rumus Newton–Cotes dengan n genap menghasilkan derajat ketepatan ekstra dibandingkan dengan n gasal.
4.2.2. Rumus Newton–Cotes terbuka Ada rumus Newton–Cotes yang tidak menggunakan salah satu atau kedua titik di ujung interval. Contoh yang paling sederhana adalah rumus titik tengah: b a + b (b − a ) 3 f ( x ) dx = ( b − a ) f ( )+ f " (η ) η ∈ [a, b] ∫a 2 24 Rumus kompositnya: b
∫ f(x)dx = Ι (f) + Ε (f) n
n
a
Ι n(f) = h[f(x1 ) + L + f(x n )] h 2(b − a) f"(η) η ∈ [a, b] 24 (b − a ) 1 dengan h = , x j = a + ( j − )h sebagai titik tengah dari titik-titik n 2 (a + ( j − 1)h, a + jh ) untuk j = 1,2,K, n .
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
E n(f) =
Rumus Newton–Cotes yang sedemikian ini disebut dengan rumus terbuka, sedangkan rumus yang terdahulu disebut tertutup. n = 2: x2 h3 f ( x ) dx = 2 hf ( x ) + f " (ξ ) 1 ∫x 3 0 n = 3: x3
∫
x0
3h 3h 3 f ( x)dx = [ f ( x1 ) + f ( x 2 )] + f " (ξ ) 2 4
n = 4: x4
∫
f ( x)dx =
x0
4h 14h 5 [2 f ( x1 ) − f ( x 2 ) + 2 f ( x 3 )] + f 3 45
( 4)
(ξ )
n = 5: Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 35 Jack la Motta
Integrasi Numeris
x5
∫
f ( x)dx =
x0
dimana h =
Buku kuliah
5h 95h 5 ( 4) [11 f ( x1 ) + f ( x 2 ) + f ( x 2 ) + f ( x3 ) + 11 f ( x 4 )] + f (ξ ) 144 24
x n − x0 , x0 = batas bawah, xn = batas atas. n
4.3. Kuadratur Gaussian Pada metoda integrasi sebelumnya, rumus integrasinya berdasarkan polinomial derajat rendah yang merupakan pendekatan f (x) dengan jumlah pias semakin besar. Kuadratur Gaussian, rumus integrasinya menggunakan polinomial yang derajatnya makin tinggi. b
n
a
j =1
I ( f ) = ∫ w( x) f ( x)dx = ∑ w j ,n f ( x j ,n ) = I n ( f )
Sebagai ilustrasi: 1
∫
−1
n
f ( x)dx = ∑ w j f ( x j ) j =1
dengan w(x) ≡ 1 Faktor pemberat {w j }dan titik nodal {x j }dipilih sedemikian sehingga kesalahan 1
n
∫
En ( f ) = D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
−1
f ( x)dx − ∑ w j f ( x j ) j =1
sama dengan nol untuk suatu polinomial f(x) dengan derajat setinggi mungkin. E n (a0 + a1 x + L + a m x m ) = a 0 E n (1) + a1 E n ( x) + L + a m E n ( x) Jadi E n ( f ) = 0 untuk setiap polinomial derajat ≤ m, bila dan hanya bila
En ( x i ) = 0
i = 0,1, K, m .
♦ Kasus 1. n = 1 . Karena hanya 2 parameter, w1 dan x1 sehingga diperlukan 2 persamaan: E n (1) = 0 E n ( x) = 0 1
∫ 1dx − w
1
1
= 0 → w1 = 2
−1
∫ xdx − w x
1 1
= 0 → x1 = 0
−1
1
sehingga
∫ f ( x)dx =& 2 f (0)
−1
♦ Kasus 2. n = 2 . Ada 4 parameter w1, w2, x1, x2, sehingga dibutuhkan 4 persamaan:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 36 Jack la Motta
Integrasi Numeris
Buku kuliah
1
E n ( x i ) = ∫ x i dx − ( w1 x1 + w2 x 2 ) = 0 untuk i = 0,1,2,3 i
i
−1
atau w1 + w2 = 2
w1 x1 + w2 x 2 = 0
2 3 3 w1 x1 + w2 x 2 = 0 3 menghasilkan rumus: 1 1 1 ∫−1 f ( x)dx = f (− 3 3 ) + f ( 3 3 ) w1 x1 + w2 x 2 = 2
2
mempunyai derajat ketelitian 3. Bandingkan dengan rumus Simpson yang menggunakan tiga titik.
♦ Kasus 3. untuk n, Terdapat 2n parameter {w1} dan {x1} sehingga terdapat 2n persamaan: E n ( x i ) = 0 i = 0,1, K,2n − 1 atau i = 1,3,K,2n − 1 ⎧⎪ 0 i 2 w x = ⎨ ∑ j j i = 0,2,K ,2n − 2 j =1 ⎪⎩ i + 1 n
persamaan diatas merupakan sistem persamaan non-linier penyelesaiannya tidak selalu jelas. Oleh karena itu digunakan cara lain.
yang
4.3.1. Kuadratur Gauss-Legendre Untuk w( x) ≡ 1 , rumus Gauss pada interval [-1,1] adalah: D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
1
∫
−1
n
f ( x)dx = ∑ w j f ( x j ) j =1
dengan titik xj adalah akar dari polinomial Legendre derajat n dalam interval [-1,1]. Faktor pemberatnya adalah −2 wi = i = 1,2, K , n ' (n + 1) Pn ( xi ) Pn +1 ( xi ) 2 2 n +1 (n!) 4 f ( 2 n ) (η ) En ( f ) = (2n + 1)[(2n)!]2 (2n)!
dan
Tabel 1 Gauss–Legendre titik dan pembobot n 2 3 4 Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
x1 ±0.5773502692 ±0.7745966692 0.0 ±0.8611363116
w1 1.0 0.5555555556 0.8888888889 0.3478546451 hal. 37 Jack la Motta
Integrasi Numeris
Buku kuliah
n 5
6
7
8
x1 ±0.3399810436 ±0.9061798459 ±0.5384693101 0.0 ±0.9324695142 ±0.6612093865 ±0.2386191861 ±0.9491079123 ±0.7415311856 ±0.4058451514 0.0 ±0.9602898565 ±0.7966664774 ±0.5255324099 ±0.1834346425
w1 0.6521451549 0.2369268851 0.4786286705 0.5688888889 0.1712344924 0.3607615730 0.4679139346 0.1294849662 0.2797053915 0.3818300505 0.4179591837 0.1012285363 0.2223810345 0.3137066459 0.3626837834
Contoh:
1
♦
n = 1:
∫
1 1 − 3 3 3 3 644744 8 6447448 f ( x)dx = f (0.5773502692) + f (− 0.5773502692)
−1
n = 3: 1
∫ f ( x)dx =
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
−1
0.5555555556 f (0.7745966692) + 0.5555555556 f (−0.7745966692) + 0.8888888889 f (0)
♦ Untuk integral pada interval umum [a,b], maka digunakan transformasi sebagai berikut: b 1 b−a ⎛ a + b + (b − a ) x ⎞ ⎟dx ∫a f (t )dt = 2 −∫1 f ⎜⎝ 2 ⎠ = 3
♦
∫ (x 1
3
∫ 1
3
⎛ a + b + (b − a) x j b−a n w j f ⎜⎜ ∑ 2 j =1 2 ⎝
⎞ ⎟⎟ ⎠
2 + x 2 + x + 1)dx = 34 , dihitung dengan kuadratur Gaussian menghasilkan: 3
⎡ ⎛ ⎛ 1 ⎞ ⎞⎤ 1 ⎛ ⎞ 3 ⎟ ⎟⎥ ⎜ 1 + 3 + (3 − 1)⎜ − 1 + 3 + (3 − 1) 3⎟ ⎢ ⎜ 3 −1 ⎢ 3 ⎠ ⎟⎥ ⎝ ⎜ 3 ⎜ ⎟ 1.0 × f f ( x)dx = + 1.0 × f ⎟⎥ ⎜ 2 ⎢ 2 2 ⎜ ⎟ ⎜ ⎟ ⎟⎥ ⎜ ⎢ ⎝ ⎠ ⎠⎦ ⎝ ⎣ 1 1 = f (2 + 3 ) + f (2 − 3 ) = 34.66666667 3 3
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 38 Jack la Motta
Integrasi Numeris
Buku kuliah
4.4. Polinomial Orthogonal Kuadratur Gauss-Legendre menggunakan polinomial orthogonal Legendre. Ada banyak famili polinomial yang orthogonal. Secara umum suatu famili polinomial gk(x) disebut orthogonal terhadap fungsi pemberat w(x) , jika : b
∫ w( x) g
( x) g m ( x)dx = 0
n
n≠m
a
b
∫ w( x)[g
( x)] dx = c(n) ≠ 0 2
n
a
Contoh : set {sin (kx)} dan {cos (kx)}
♦ Polinomial Legendre. Pn(x) → orthogonal pada interval [-1,1] terhadap w(x) = 1 1
∫ P ( x) P n
m
( x)dx = 0
n≠m
−1 1
∫ [P ( x)] dx = c(n) ≠ 0 2
n
−1
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Beberapa Pn(x): P0(x) = 1 P1(x) = x P2(x) = ½(3x2-1) P3(x) = ½(5x3-3x) ) P4(x) = ⅛(35x4-30x2+3) Rumus rekursiv: 2n − 1 n −1 Pn ( x) = xPn −1 ( x) − Pn − 2 ( x) n n
♦ Polinomial Laquerre. Ln(x) → orthogonal pada interval [0,∞] terhadap w(x) = e-x ∞
∫ L ( x) L n
m
( x)dx = 0
n≠m
0
∞
∫ [L
( x)] dx = c(n) ≠ 0 2
n
0
Beberapa Ln(x): L0(x) = 1 L1(x) = -x+1 L2(x) = x2-4x+2 L3(x) = -x3+9x2-18x+6 Rumus rekursiv:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 39 Jack la Motta
Integrasi Numeris
Buku kuliah
Ln ( x) = (2n − x − 1) Ln −1 ( x) − (n − 1) 2 Ln − 2 ( x) ♦ Polinomial Chebysev. Tn(x) → orthogonal pada interval [-1,1] terhadap w(x) = 1/√(1-x2) 1 Tn ( x)Tm ( x) ∫−1 1 − x 2 dx = 0 n ≠ m 1 [Tn ( x)]2 dx = c(n) ≠ 0 ∫ 2 −1 1 − x Beberapa Tn(x): T0(x) = 1 T1(x) = x T2(x) = 2x2-1 T3(x) = 4x3-3x Rumus rekursiv: Tn ( x) = 2 xTn −1 ( x) − Tn − 2 ( x)
♦ Polinomial Hermite. Hn(x) → orthogonal pada interval [-∞,∞] terhadap w( x) = e − x
2
∞
∫e
−∞ ∞
− x2
H n ( x) H m ( x)dx = 0
∫ e [H − x2
n≠m
( x)] dx = c(n) ≠ 0 2
n
−∞
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Beberapa Hn(x): H0(x) = 1 H1(x) = 2x H2(x) = 4x2-2 H3(x) = 8x3-12x Rumus rekursiv: H n ( x) = 2 xH n −1 ( x) − 2(n − 1) H n − 2 ( x)
4.4.1. Kuadratur Gauss-Laquerre ∞
n
0
j =1
−x ∫ e f ( x)dx = ∑ w j f ( x j )
Kuadratur ini dapat digunakan untuk menghitung:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 40 Jack la Motta
Integrasi Numeris
Buku kuliah
∞
∞
a
0
−t −a −x ∫ e f (t )dt = e ∫ e f ( x + a)dx n
= e −a ∑ w j f ( x j + a) j =1
dengan wj = faktor pemberat, xj = akar dari polinomial Laquerre
Tabel 2. wj dan xj dari kuadratur Gauss-Laguerre n 2 3
4
5
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
6
10
14
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
wj 0.85355 33905 0.14644 66094 0.71109 30099 0.27851 77335 0.01038 92565 0.60315 41043 0.35741 86924 0.03888 79085 0.00053 92947 0.263560319718 1.413403059107 3.596425771041 7.085810005859 12.640800844276 0.222846604179 1.188932101673 2.992736326059 5.775143569105 9.837467418383 15.982873980602 0.13779347054 0.729454549503 1.80834290174 3.401433697855 5.552496140064 8.330152746764 11.8437858379 16.279257831378 21.996585811981 29.920697012274 0.093307812017 0.492691740302 1.215595412071) 2.269949526204)
xj 0.58578 64376 3.41421 35623 0.41577 45567 2.29428 03602 6.28994 50829 0.32254 76896 1.74576 11011 4.53662 02969 9.39507 09123 0.521755610583 0.398666811083 0.0759424496817 0.00361175867992 2.33699723858E-05 0.45896467395 0.417000830772 0.113373382074 0.0103991974531 0.000261017202815 8.9854790643E-07 0.308441115765 0.401119929155 0.218068287612 0.0620874560987 0.00950151697518 0.000753008388588 0.000028259233496 4.24931398496E-07 1.83956482398E-09 9.91182721961E-13 0.21823488594 0.342210177923 0.263027577942 0.126425818106
hal. 41 Jack la Motta
Integrasi Numeris
Buku kuliah
n
wj 3.667622721751 5.425336627414 7.565916226613 10.120228568019 13.130282482176 16.65440770833 20.776478899449 25.623894226729 31.407519169754 38.530683306486 48.026085572686
xj 0.040206864921 0.00856387780361 0.00121243614721 0.000111674392344 6.45992676202E-06 2.2263169071E-07 4.22743038498E-09 3.92189726704E-11 1.45651526407E-13 1.48302705111E-16 1.60059490621E-20
4.4.2. Kuadratur Gauss-Chebysev 1
∫
−1
dengan
π
wj =
n
1− x2
b−a 2 −∫1 1
f (t )dt =
a
=
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
f(x)dx = ∑ w j f ( x j ) j =1
2 j −1 π ) , j = 1, 2, ... , n 2n ⎡ 1 ⎛ a + b + (b − a ) x ⎞⎤ 1− x2 f ⎜ ⎟⎥ dx ⎢ 2 ⎝ ⎠⎦ 1− x2 ⎣
dan x j = cos(
b
∫
n
1
(b − a )π 2n
n
∑ j =1
⎛ a + b + (b − a) x j 1 − x 2j f ⎜⎜ 2 ⎝
⎞ ⎟⎟ ⎠
4.4.3. Kuadratur Gauss-Hermite ∞
n
−x ∫ e f (x)dx = ∑ w j f ( x j ) 2
j =1
−∞
Tabel 3. wj dan xj dari kuadratur Gauss - Laguerre n 2 3 4 5
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
wj 0.88622 69255 0.29540 89752 1.18163 59006 0.08131 28354 0.80491 40900 0.01995 32421 0.39361 93232 0.94530 87205
xj ± 0.70710 67811 ± 1.22474 48714 0.0 ± 1.65068 01239 ± 0.52464 76233 ± 2.02018 28705 ± 0.95857 24646 0.0
hal. 42 Jack la Motta
Bab
5. SISTEM PERSAMAAN LINIER
5.1. Eliminasi Gauss Eliminasi Gauss digunakan mencari akar sistem persamaan linier. f1 ( x1 , x 2 ,..., x n ) = 0
f 2 ( x1 , x 2 ,..., x n ) = 0
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
M f n ( x1 , x 2 ,..., x n ) = 0 Contoh: Ditinjau sistem persamaan: 2x1 - 7x2 + 4x3 = 9 x1 + 9x2 - 6x3 = 1 -3x1 + 8x2 + 5x3 = 6 yang akarnya adalah x1 = 4, x2 = 1, dan x3 = 2 Persamaan diatas dalam bentuk matrik dapat ditulis sebagai berikut: [B]{x}={u} 4⎤ ⎧ x1 ⎫ ⎧9⎫ ⎡ 2 −7 ⎪ ⎪ ⎪ ⎪ ⎢ 1 9 − 6⎥⎥ ⎨ x 2 ⎬ = ⎨1⎬ ⎢ ⎢⎣− 3 8 5⎥⎦ ⎪⎩ x3 ⎪⎭ ⎪⎩6⎪⎭ Untuk menjelaskan eliminasi Gauss, maka dibentuk suatu matrik sebagai berikut: ⎡ 2 −7 4 9 1 0 0⎤ ⎢ [B u I ] = ⎢ 1 9 − 6 1 0 1 0⎥⎥ ⎢⎣− 3 0 0 1⎥⎦ 8 5 6 Kita kalikan baris 1 dengan 1/2, tambahkan (-1 x baris 1 yang baru) kepada baris 2, dan tambahkan (3x baris 1 yang baru) kepada baris 3.
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 43 Jack la Motta
Sistem Persamaan Linier
⎡1 − 7 / 2 2 ⎢ ⎢0 25 / 2 − 8 ⎢⎣0 − 5 / 2 11
Buku kuliah
9/2 −7/2 39 / 2
1 / 2 0 0⎤ ⎥ − 1 / 2 1 0⎥ 3 / 2 0 1⎥⎦
Operasi diatas sama dengan pembentukan/pengubahan sistem persamaan asli menjadi 7 9 x1 − x 2 + 2 x3 = 2 2 25 7 x 2 − 8 x3 = − 2 2 5 39 − x 2 + 11x3 = 2 2 Perhatikan bahwa operasi di atas jika ditulis dalam bentuk matrik adalah ⎡ 1 ⎤ ⎢ 2 0 0⎥ ⎡ 2 − 7 1 0 0⎤ 9 4 ⎢ 1 ⎥ ⎢ ⎥ ⎢− 0 1 0⎥ 1 9 −6 1 0⎥ ⎢ 1 ⎢ 2 ⎥ ⎢ 0 0 1⎥⎦ 6 8 5 ⎢ 1 0 1⎥ ⎣− 3 ⎢⎣ 2 ⎥⎦
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Selanjutnya dilakukan operasi sebagai berikut: kalikan baris 2 dengan 2/25 dan tambahkan (5/2 x baris 2 yang baru) kepada baris 3 ⎡1 − 7 / 2 2 9/2 1/ 2 0 0⎤ ⎢ ⎥ 1 − 16 / 25 − 7 / 25 − 1 / 25 2 / 25 0⎥ ⎢0 ⎢⎣0 0 47 / 5 94 / 25 7/5 1 / 5 1⎥⎦ Operasi terakhir mengubah sistem persamaan menjadi: 7 9 x1 − x 2 + 2 x3 = 2 2 16 7 x2 − x3 = − 25 25 47 94 x3 = 25 5 Kalikan baris 3 dengan 5/47. Tambahkan ke baris 2: (16/25 x baris 3 yang baru). Tambahkan ke baris 1: (-2 x baris 3 yang baru) ⎡1 − 7 / 2 0 1/ 2 19 / 24 − 2 / 47 − 10 / 47 ⎤ ⎢ ⎥ 1 0 1 13 / 235 22 / 235 16 / 235⎥ ⎢0 ⎢⎣0 0 1 2 7 / 47 1 / 47 5 / 47⎥⎦ Akhirnya tambahkan ke baris 1: (7/2 x baris 2)
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 44 Jack la Motta
Sistem Persamaan Linier
Buku kuliah
⎡1 0 0 ⎢ ⎢0 1 0 ⎢⎣0 0 1
4 1 2
6 / 235⎤ ⎥ 13 / 235 22 / 235 16 / 235⎥ 7 / 47 1 / 47 5 / 47⎥⎦
93 / 235 67 / 235
Jadi sistem persamaan menjadi x1 = 4, x2 = 1, x3 = 2 dan inverse matrik [B] adalah ⎡93 / 235 67 / 235 6 / 235 ⎤ ⎢13 / 235 22 / 235 16 / 235⎥ ⎥ ⎢ ⎢⎣ 7 / 47 1 / 47 5 / 47 ⎥⎦
5 ⎞ ⎛1 2 Dari pengamatan: det B = ⎜ × × ⎟ ⎝ 2 25 47 ⎠ Jadi kalau di ‘resume’
[B
−1
= 235
u I] ⇓
⎡ 2 −7 4 ⎢ 9 −6 ⎢ 1 ⎢⎣− 3 8 5
1 0 0⎤ ⎥ 0 1 0⎥ 0 0 1⎥⎦
9 1 6 ⇓
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
⎡1 0 0 ⎢ ⎢0 1 0 ⎢⎣0 0 1
4 1 2
93 / 235 67 / 235 6 / 235⎤ ⎥ 13 / 235 22 / 235 16 / 235⎥ 7 / 47 1 / 47 5 / 47⎥⎦
[I
x B −1
]
5.2. Eliminasi Gauss–Jordan Pada eliminasi Gauss di atas secara garis besar terdiri dari beberapa langkah: a. operasi normalisasi: elemen diagonal diubah menjadi bernilai 1 b. operasi reduksi: elemen non-diagonal diubah menjadi bernilai 0 Pada eleminasi Gauss – Jordan operasi a & b dikerjakan bersamaan. Contoh: 2x1 - 7x2 + 4x3 = 9 x1 + 9x2 - 6x3 = 1 -3x1 + 8x2 + 5x3 = 6 ⎡ 2 −7 4 9 1 0 0⎤ ⎢ [B u I ] = ⎢ 1 9 − 6 1 0 1 0⎥⎥ ⎢⎣− 3 8 5 6 0 0 1⎥⎦
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 45 Jack la Motta
Sistem Persamaan Linier
Buku kuliah
Normalisasi baris 1 dengan membaginya dengan elemen ‘pivot’ = 2, kemudian: a. baris 2 - baris 1 yang baru b. baris 3 + 3x baris 1 yang baru ⎡1 − 7 / 2 2 9/2 ⎢ −7/2 ⎢0 25 / 2 − 8 ⎢⎣0 − 5 / 2 11 39 / 2
1 / 2 0 0⎤ ⎥ − 1 / 2 1 0⎥ 3 / 2 0 1⎥⎦
Normalisasi baris 2 dengan membaginya dengan elemen ‘pivot’ = 25/2, kemudian: a. kurangi (-7/2 x baris 2 yang baru) dari baris 1 b. kurangi (-5/2 x baris 2 yang baru) dari baris 3 ⎡ 1 0 − 6 / 25 88 / 25 9 / 25 7 / 5 0⎤ ⎢ ⎥ − 7 / 25 − 1 / 25 2 / 25 0⎥ ⎢0 1 − 16 / 25 ⎢⎣0 0 47 / 5 94 / 25 7/5 1 / 5 1⎥⎦ Normalisasi baris 3 dengan membaginya dengan elemen ‘pivot’ = 47/5, kemudian: a. kurangi (-6/25 x baris 3 yang baru) dari baris 1 b. kurangi (-16/25 x baris 3 yang baru) dari baris 2 ⎡1 0 0 4 93 / 235 67 / 235 6 / 235⎤ ⎢ ⎥ 1 13 / 235 22 / 235 16 / 235⎥ = I x B −1 ⎢0 1 0 ⎢⎣0 0 1 2 7 / 47 1 / 47 5 / 47⎥⎦
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
[
]
25 47 ⎞ ⎛ det B = ⎜ 2 × × ⎟ = 235 2 5 ⎠ ⎝
5.3. Eliminasi Gauss–Jordan dengan ‘pivot’ maksimum Jika matrik [B] mempunyai salah satu elemen yang mempunyai nilai kecil sekali dibandingkan elemen yang lain, maka cara’pivoting’ yang sebelumnya dapat memberikan hasil yang tidak akurat. Oleh karena itu dipilih elemen ‘pivot’ yang mempunyai nilai terbesar. Contoh: -3x1 + 8x2 + 5x3 = 6 2x1 - 7x2 + 4x3 = 9 x1 + 9x2 - 6x3 = 1
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 46 Jack la Motta
Sistem Persamaan Linier
Buku kuliah
⎡− 3 8 5 ⎢ [B u I ] = ⎢ 2 − 7 4 ⎢⎣ 1 9 −6 Dipilih elemen b32 = 9 sebagai ‘pivot’ ⎡−3 8 5 6 ⎢ 4 9 ⎢ 2 −7 ⎢⎣1 / 9 1 − 6/9 1/ 9 selanjutnya direduksi sebagai berikut: ⎡− 35 / 9 0 31 / 3 46 / 9 ⎢ 88 / 9 ⎢ 25 / 9 0 − 2 / 3 ⎢⎣ 1 / 9 1 − 2 / 3 1/ 9
1 0 0⎤ ⎥ 0 1 0⎥ 0 0 1⎥⎦
6 9 1
0 ⎤ ⎥ 0 1 0 ⎥ 0 0 1 / 9⎥⎦
1 0
1 0 − 8 / 9⎤ ⎥ 0 1 7/9 ⎥ 0 0 1 / 9 ⎥⎦
(I)
Operasi 2:
⎡− 35 / 93 0 1 ⎢ ⎢ 235 / 93 0 0 ⎢⎣ − 13 / 93 1 0
46 / 93 940 / 93 41 / 93
−8/9⎤ ⎥ 2 / 31 1 67 / 93⎥ 2 / 31 0 5 / 93 ⎥⎦ 3 / 31 0
(II)
Operasi 3
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
⎡0 0 ⎢ ⎢1 0 ⎢⎣0 1 [ ?1
1 0
2 4
0
1
x
?2
7 / 47 1 / 47 ⎤ ⎥ 93 / 235 67 / 235⎥ 16 / 235 13 / 235 22 / 235⎥⎦ 5 / 47 6 / 235
]
(III) (IIIa)
Dari hasil terakhir (III) terlihat bahwa akar persamaan {x} dapat dislesaikan, tetapi bagaimana matrik ?1 & ?2. Sebetulnya [?1] elemennya adalah elemen [B]-1, hanya letaknya tidak betul sehingga perlu diatur untuk mendapatkan inverse [B] yang sesungguhnya. Ada butir yang sangat penting dari hasil diatas:
• Akar dari [B]{x}={u} dapat dicari tanpa menghitung [B]-1. • Hitungan inverse suatu matrik lebih baik dihindari karena mahal beayanya. Untuk menghemat memori komputer, maka pada cara terakhir (eliminasi Gauss–Jordan dengan ‘pivot’ maksimum) hasil dari inverse dimasukan kedalam matrik [B].
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 47 Jack la Motta
Sistem Persamaan Linier
Buku kuliah
5.3.1. Rekonstruksi pembentukan “scrambled inverse” ♦ Pembentukan I Yang Semula I ⎡− 35 / 9 0 31 / 3 ⎢ ⎢ 25 / 9 0 − 2 / 3 ⎢⎣ 1 / 9 1 − 2 / 3
46 / 9 88 / 9 1/ 9
1 0 − 8 / 9⎤ ⎥ 0 1 7/9 ⎥ 0 0 1 / 9 ⎥⎦
dipindah kolom 3 dari [I] dipindah
b32
kolom 2 dari [B]
Gambar 11 Cara pertama pemindahan kolom dengan elemen pivot Yang Baru I ⎡− 35 / 9 − 8 / 9 31 / 3 46 / 9⎤ ⎢ 25 / 9 7 / 9 − 2 / 3 88 / 9⎥⎥ ⎢ ⎢⎣ 1 / 9 1 / 9 − 2 / 3 1 / 9⎥⎦
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
♦ Pembentukan II Yang Semula II
dipindah
⎡− 35 / 93 0 1 ⎢ ⎢ 235 / 93 0 0 ⎢⎣ − 13 / 93 1 0
46 / 93 940 / 93 41 / 93
−8/9⎤ ⎥ 2 / 31 1 67 / 93⎥ 2 / 31 0 5 / 93 ⎥⎦ 3 / 31 0
kolom 1 dari [I]
b13
dipindah kolom 3 dari [B]
Gambar 12 Cara kedua pemindahan kolom dengan elemen pivot Yang Baru II ⎡− 35 / 93 − 8 / 9 3 / 31 46 / 93⎤ ⎢ 235 / 93 67 / 93 2 / 31 940 / 93⎥ ⎥ ⎢ ⎢⎣ − 13 / 93 5 / 93 2 / 31 41 / 93⎥⎦
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 48 Jack la Motta
Sistem Persamaan Linier
Buku kuliah
♦ Pembentukan III Yang Semula III
dipindah
⎡0 0 1 ⎢ ⎢1 0 0 ⎢⎣0 1 0
2 4 1
1 / 47 ⎤ ⎥ 6 / 235 93 / 235 67 / 235⎥ 16 / 235 13 / 235 22 / 235⎥⎦ 5 / 47
7 / 47
kolom 2 dari [I] dipindah
b21
kolom 1 dari [B]
Gambar 13 Cara kedua pemindahan kolom dengan elemen pivot Yang Baru III 1 / 47 5 / 47 2⎤ ⎡ 7 / 47 ⎢93 / 235 67 / 235 6 / 235 4⎥ ⎥ ⎢ ⎢⎣13 / 235 22 / 235 16 / 235 1 ⎥⎦
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
“scrambled inverse”
{x} → akar dari sistem persamaan
perlu diatur dahulu untuk menghasilkan inverse yang sesungguhnya
5.4. Metoda Iterasi
5.4.1. Metoda Jacobi Kita bahas sistem persamaan: [B]{x} = {u} atau b11 x11 + b12 x 2 + L + b1n x n = u1 b21 x1 + b22 x 2 + L + b2 n x n = u 2 (A)
M M bn1 x1 + bn 2 x 2 + L + bnn x n = u n Metoda Jacobi membentuk persamaan untuk mendekati persamaan di atas:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 49 Jack la Motta
Sistem Persamaan Linier
Buku kuliah
u1 − b12 x 2 − b13 x3 − ... − b1n x n b11 u − b21 x1 − b23 x3 − ... − b2 n x n x2 = 2 b22 u n − bn1 x1 − bn 2 x 2 − ... − bn ,n −1 x n −1 xn = bnn
x1 =
(B)
n
u i − ∑ bij x j xi =
j =1 j ≠i
, i = 1,2,L, n bii jika terjadi bii = 0 atau nilainya kecil, maka harus diadakan pengaturan sehingga bii ≠ 0 atau
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Metoda ini dimulai dengan “tebakan” nilai awal {x0} kemudian dimasukkan kedalam persamaan (B) untuk menghitung {x} baru. Contoh: 4 x1 + 2 x 2 + x3 = 11⎫ ⎧1⎫ ⎪ ⎪ ⎪ − x1 + 2 x 2 = 3 ⎬ ⇒ {x} = ⎨2⎬ ⎪3 ⎪ 2 x1 + x 2 + 4 x3 = 16⎪⎭ ⎩ ⎭ Persamaan di atas ditulis lagi: 11 1 1 x1 = − x 2 − x3 4 2 4 3 1 x 2 = + x1 2 2 1 1 x3 = 4 − x1 − x 2 2 4 t Vektor awal {x0} = [1, 1, 1] 11 1 1 x11 = − × 1 − × 1 = 2 4 2 4 3 1 x 21 = + × 1 = 2 2 2 1 1 13 x31 = 4 − × 1 − × 1 = 2 4 4 13 t Jadi , {x1 } = [2, 2, ] 4 {x 2 } = [0.9375 2.5 2.5]t
(C)
{x3 } = [0.875 1.96875 2.90625]
t
{x 4 } = [1.03906 1.9375 3.0703]
t
{x5 } = [1.01367 2.0195 2.9961]
t
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 50 Jack la Motta
Sistem Persamaan Linier
Buku kuliah
M M t {x14 } = [1.0000 2.0000 3.0000] Dalam metoda ini hitungan elemen vektor yang baru menggunakan elemen vektor yang lama.
5.4.2. Metoda Gauss-Seidel Dibandingkan dengan metoda Jacobi, metoda Gauss-Seidel menghitung elemen vektor baru dengan menggunakan elemen yang baru saja dihitung. Contoh: digunakan sistem persamaan yang digunakan sebelumnya, jadi persamaan (C) dapat digunakan. Vektor awal {x0} = [1, 1, 1]t 11 1 1 x11 = − × 1 − × 1 = 2 4 2 4 3 1 5 x 21 = + × 2 = 2 2 2 1 1 5 19 x31 = 4 − × 2 − × = 2 4 2 8 t
5 19 ⎤ ⎡ Jadi, {x1 } = ⎢2, , 2 8 ⎥⎦ ⎣ t {x 2 } = [0.9063 1.9531 3.0586] {x3 } = [1.0088 2.0044 2.9945]
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
t
{x 4 } = [0.9992 1.9996 3.0005]
t
{x5 } = [1.0000 2.0000 3.0000]
t
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 51 Jack la Motta
Bab
6. MATRIK
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
6.1. Notasi dan Konsep-konsep Pendahuluan ⎡ a11 ⎢ A = (aij ) = ⎢ M ⎢ ⎢⎣a m1 ⎡b11 ⎢ B = (bij ) = ⎢ M ⎢ ⎢⎣bn1 ⎡ c11 ⎢ C = (cij ) = ⎢ M ⎢ ⎢⎣c n1 ⎡ d11 ⎢ D = (d ij ) = ⎢ ⎢d p1 ⎣
a12
am2 b12
bn 2 c12
cn 2
d12 d p2
L a1n ⎤ ⎥ M ⎥ m×n ⎥ L a mn ⎥⎦ L b1 p ⎤ ⎥ n× p M ⎥ ⎥ L bnp ⎥⎦ L c1 p ⎤ ⎥ n× p M ⎥ ⎥ L c np ⎥⎦ L d1q ⎤ ⎥ n×q M ⎥ L d pq ⎥⎦
Dengan notasi diatas, maka hal-hal dibawah ini berlaku: 1. E = B + C adalah matriks n x p dengan eij = bij + cij B+C =C + B n
2. F = AB adalah matriks m x p dengan f ij = ∑ aik bkj k =1
A( B + C ) = AB + AC , ( B + C ) D = BD + CD Jika AB = BA maka A & B disebut ‘commute’ Contoh:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 52 Jack la Motta
Matrik
Buku kuliah
⎡− 1 3⎤ ⎡ 1 2⎤ ⎢ ⎥ ⎢ ⎥ ⎡1 2 − 1⎤ A=⎢ ⎥ dan B = ⎢ 2 1⎥ dan C = (cij ) = ⎢− 3 1⎥ ⎣3 4 2 ⎦ ⎢ 0 3⎥ ⎢ 2 4⎥ ⎣ ⎦ ⎣ ⎦ maka ⎡ 0 5⎤ ⎡15 20 10⎤ ⎢ ⎢ ⎥ ⎥ B + C = ⎢− 1 2⎥ (B + C) A = ⎢ 5 6 5 ⎥ ⎢ 2 7⎥ ⎢23 32 12⎥ ⎣ ⎣ ⎦ ⎦ ⎡− 4 2 ⎤ A( B + C ) = ⎢ ⎥ ⎣ 0 37⎦ 3. ⎡1 3 1 ⎤ ⎡ 2 1 3⎤ ⎥ ⎢ A = ⎢2 0 − 1⎥ dan B = ⎢⎢− 3 1 2⎥⎥ ⎢⎣4 1 2 ⎥⎦ ⎢⎣ 3 1 2⎥⎦
maka 4 0⎤ ⎡− 4 5 5 0⎥⎥ AB = ⎢⎢ 1 1 ⎢⎣ 11 7 14 7 ⎥⎦ 7⎤ ⎡ 10 6 14⎤ B 2 = ⎢⎢− 3 0 − 3⎥⎥ 0⎥⎥ ⎢⎣ 9 6 15⎥⎦ 6⎥⎦ ⎡− 19 − 6 − 10⎤ ( A − B)( A + B) = ⎢⎢ − 5 13 7 ⎥⎥ dan A 2 ⎢⎣ 3 4 4⎥⎦
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
⎡ 11 A 2 = AA = ⎢⎢− 2 ⎢⎣ 14 ⎡16 9 BA = ⎢⎢ 7 − 7 ⎢⎣13 11
11⎤ 4 ⎥⎥ 18⎥⎦
⎡ 1 − 2 − 14⎤ − B 2 = ⎢⎢ 1 5 3⎥⎥ ⎢⎣5 8 − 8⎥⎦
tampak ternyata bahwa (A-B)(A+B) ≠ (A2-B2), ternyata (A-B)(A+B) = A2-BA+AB-B2) 4. At = A tranpos → G = At , gij = aji (B + C)t = Bt + Ct, (ABD)t = DtBtAt Jika A = At maka A adalah simetris ⎡1 3 1 ⎤ ⎡1 2 4⎤ ⎥ ⎢ t Contoh: A = ⎢2 0 − 1⎥ A = ⎢⎢3 0 1 ⎥⎥ ⎢⎣4 1 2 ⎥⎦ ⎢⎣1 − 1 2⎥⎦ H adalah matrik diagonal jhj, H adalah n x n dan hij = 0 untuk i ≠ j. Jika hii ditulis sebagai hi, matrik diagonal H sering ditulis [ h1, h2, …, hn] Perhatikan: AH adalah matrik mxn → (hjaij)
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 53 Jack la Motta
Matrik
Buku kuliah
HB adalah matrik nxp → (hjbij) H adalah matrik skalar jhj elemennya hii mempunyai nilai sama yaitu h. hA = Ah = AH, hB = Bh = HB Jika h =1, maka H disebut matrik identitas dan biasanya ditulis sebagai I IA = AI = A ⎡h11 0 L 0 ⎤ ⎡1 0 L 0⎤ ⎥ ⎢0 h ⎢0 1 L 0 ⎥ L 0⎥ 22 ⎥ Contoh: H = ⎢ dan I = ⎢ ⎢ M ⎢M M ⎥ M⎥ ⎥ ⎢ ⎢ ⎥ 0 L hnn ⎦ ⎣0 0 L 1 ⎦ ⎣0 5. Jika AK = I, maka K adalah unik dan disebut invers dari A yang biasa ditulis sebagai A-1. Matrik A disebut non-singuler jhj A-1 ada 6. Determinan Jika terdapat suatu matrik, A, nxn, maka dapat dihitung determinan A yang diberi notasi det [A], det A, atau |A|. Beberapa sifat determinan: det (At) = det (A), det (A-1) = [det(A)]-1 det (AB) = det (A) det (B) Untuk menghitung determinan dibutuhkan beberapa definisi: 1. Minor dari aij adalah determinan dari suatu (n-1) x (n-1) matrik yang dibentuk dari matrik A, nxn, dimana baris dan kolom yang berisi aij dibuang
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
2. Kofaktor dari aij adalah suatu bilangan hasil perkalian antara (-1)i+j dikalikan dengan minor dari aij , dan diberi notasi Aij Determinan A dapat dihitung sebagai: n
n
j =1
i =1
det (A) = ∑ a rj Arj = ∑ a is Ais
Contoh:
atau
a11 kofaktor a11 1 4 ⎡ ⎤ det ⎢ ⎥ = 1 × 5 + 4(− 2) = −3 ⎣2 5⎦ kofaktor a12 a12 = 5 x 1 + 4 (-2) = -3 a22
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
kofaktor a22
hal. 54 Jack la Motta
Matrik
Buku kuliah
kofaktor a13 kofaktor a11 kofaktor a12 644744 8 647 4 48 4 644744 8 3 − 1⎤ ⎡5 −1 2 2 1+1 − 1 1+ 2 1 1+ 3 1 + 3 × (− 1) + (− 1) × (− 1) det ⎢⎢ 1 − 1 2⎥⎥ = 5 × (− 1) −3 4 2 4 2 −3 ⎢⎣2 − 3 4⎥⎦ = 5 × (+1) × (−4 + 6) + 3 × (−1) × (4 − 4) − 1 × (+1) × (−3 + 2) = 11 7. Vektor adalah matrik dengan kolom tunggal Contoh: ⎧ u1 ⎫ ⎧ v1 ⎫ ⎪u ⎪ ⎪v ⎪ ⎪ 2⎪ ⎪ ⎪ u = ⎨ ⎬ v = ⎨ 2⎬ ⎪M⎪ ⎪M⎪ ⎪⎩u n ⎪⎭ ⎪⎩vn ⎪⎭ n
' inner product ' (dotproduct ) : (u, v) = ∑ u i vi = u t v = v t u i =1
Beberapa sifat: (u, v+w) = (u,v)+(u,w) (u+v, w) = (u,w)+(v,w) (αu,v) = α(u,w)+(v,w) (αu,v) = α(u,v) α skalar (u, u) ≥ 0, (u,u) = 0 jhj {u} = 0
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
6.2. Determinan dan invers Secara numeris determinan dan invers dapat diselesaikan dengan metode Gauss beserta variannya yang sudah dijelaskan pada Bab 5. Jika matrik A adalah nxn, maka pernyataan-pernyataan dibawah ini adalah equivalen: 1. [A]{x} = {B} mempunyai penyelesaian yang unik 2. [A]{x} = 0 berarti {x} = 0 3. [A]-1 ada 4. det (A) ≠ 0
6.2.1. Menghitung determinan dengan eleminasi segitiga atas Dengan metode kofaktor dihitung determinan suatu matrix A, nxn, dengan ekspansi terhadap kolom pertama.
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 55 Jack la Motta
Matrik
Buku kuliah
n
det A = ∑ ai1 Ai1 = a11 A11 + a 21 A21 + ... + a n1 An1 i =1
Jika Adalah matrik segitiga atas, maka ai1 = 0, i = 2,3,…,n ∴ det A = a11A11 = a11a22 … ann Jika A matrik tidak segitiga, maka ⎡u11 u12 L u1n ⎤ ⎡ a11 a12 L a1n ⎤ ⎢0 a ⎥ ⎢a a 22 L a 2 n ⎥ L u 2 n ⎥⎥ 21 22 ⎢ ⎢ =K det A = det ⎢ M ⎢ M M ⎥ M ⎥ ⎥ ⎢ ⎥ ⎢ 0 L u nn ⎦ ⎣0 ⎣a n1 a n 2 L a nn ⎦ Contoh:
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
1 4 2 −1 3 −5 6 1 − 3B1 + B2 A= 1 0 4 5 − B1 + B3 1 9 7 6 B1 + B4 −6 1 4 2 −1 0 − 17 0 4 B2 / (− 17 ) A= 0 0 2 6 0 1 21 1 1 4 2 −1 0 1 0 − 4 / 17 A = −17 4 B2 + B3 0 −4 2 6 0 1 21 1 − 25B2 + B4 1 0 A = −17 0 0
4 2 −1 1 0 − 4 / 17 1 0 2 86 / 17 B3 2 0 21 117 / 17
1 0 A = (2)(−17) 0 0 1 4 0 1 A = −34 0 0 0 0
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
4 2 −1 1 0 − 4 / 17 0 1 43 / 17 0 21 117 / 17 − 21B3 + B4 2 −1 0 − 4 / 17 = (− 34)(1)(1)(1)(− 786 / 17 ) = 1572 1 43 / 17 0 − 786 / 17
hal. 56 Jack la Motta
Matrik
Buku kuliah
6.3. Matrik dan Vektor Eigen Pada setiap matrik A, nxn, terdapat satu set vektor yang disebut vektor eigen dan satu set skalar yang disebut nilai eigen. Vektor u disebut eigen dari matrik A jhj u vektor tidal nol dan λ adalah suatu skalar (yang mungkin nol nilainya), sehingga [A] {u} = λ{u} (I) Skalar λ disebut nilai eigen dari matrik A Pers. (I) dapat ditulis sebagai: [A-λI]{u} = {0} (II) Pers. (II) mempunyai penyelesaian dengan {u} tidak nol, jika Φ(λ) ≡ det (A-λI) = 0 Φ(λ) disebut fungsi karakteristik dari matrik A. a12 a1n ⎤ L ⎡a11 − λ ⎢ a a 22 − λ L a 2 n ⎥⎥ A − λI = ⎢ 21 ⎢ M M L M ⎥ ⎥ ⎢ an2 L a nn − λ ⎦ ⎣ a n1
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Contoh: Menentukan vektor & nilai eigen matrik ⎡ 2 − 1 0⎤ A = ⎢⎢− 1 3 − 1⎥⎥ ⎢⎣ 0 − 1 2⎥⎦ 2−λ −1 0 3−λ −1 −1 −1 φ (λ ) = − 1 3 − λ − 1 = (2 − λ ) + −1 2 − λ 0 2−λ 0 −1 2 − λ = (2 − λ ){(3 − λ )(2 − λ ) − 1} − (2 − λ + 0) = λ3 − 7λ2 + 14λ − 8
Akar dari Φ(λ)=0 adalah nilai eigen λ1=1, λ2=2, λ3=4. Vektor eigen untuk λ1=1 0⎤ ⎧0⎫ ⎡2 − 1 − 1 ⎪ ⎪ ⎥ ⎢ −1 3 −1 − 1⎥{u} = ⎨0⎬ ⎢ ⎪0⎪ ⎢⎣ 0 − 1 2 − 1⎥⎦ ⎩ ⎭ Dengan eliminasi Gauss-Jordan ⎡ 1 − 1 0 0 ⎤ B2 + B1 ⎡ 1 −1 0 0 ⎤ ⎥ ⎢ ⎥ ⎢ 1 −1 0 ⎥ ⎢− 1 2 − 1 0 ⎥ B1 + B2 ⇒ ⎢0 ⎢⎣0 − 1 1 0 ⎥⎦ B2 + B3 ⎢⎣ 0 − 1 1 0 ⎥⎦
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 57 Jack la Motta
Matrik
⎡1 0 −1 0 ⎤ ⎥ ⎢ ⎢0 1 − 1 0 ⎥ ⎢⎣0 0 0 0 ⎥⎦
Buku kuliah
ini berarti terdapat 2 persamaan linier untuk 3 bilangan tak u1 − u 3 = 0⎫ diketahui: ⎬ u1 = u 2 = u 3 u 2 − u 3 = 0⎭
⎧1⎫ ⎪⎪ Jadi vektor eigen untuk λ1=1 adalah: {u}1 = u 3 ⎨1⎬ ⎪1⎪ ⎩⎭
Vektor eigen untuk λ2=2 ⎧0⎫ tampak bahwa baris 1 & 3 identik, sehingga sehingga ⎡ 0 − 1 0⎤ ⎪ ⎪ ⎥ ⎢− 1 1 − 1⎥{u} = ⎨0⎬ hanya terdapat 2 persamaan: ⎢ ⎪ ⎪ ⎢⎣ 0 − 1 0⎥⎦ − u 2 = 0⎫ ⎩0⎭ ⎬ u1 = −u 3 , u 2 = 0 u1 + u 2 − u 3 = 0⎭ Vektor eigen untuk λ3=4 0⎤ ⎧0⎫ ⎡− 2 − 1 ⎢ − 1 − 1 − 1⎥{u} = ⎪0⎪ Gauss ⎨ ⎬ ⇒ ⎥ ⎢ ⎪0⎪ Jordan ⎢⎣ 0 − 1 − 2⎥⎦ ⎩ ⎭
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
⎧0⎫ ⎡ 1 0 − 1⎤ ⎢0 1 2⎥{u} = ⎪0⎪ ⎨ ⎬ ⎥ ⎢ ⎪0⎪ ⎢⎣0 0 0⎥⎦ ⎩ ⎭ ⎧1 ⎫ u1 − u 3 = 0⎫ ⎪ ⎪ ⎬ u1 = u 3 , u 2 = −2u 3 ⇒ {u}3 = u 3 ⎨− 2⎬ u 2 + 2u 3 = 0⎭ ⎪1 ⎪ ⎩ ⎭
6.3.1. Metode ‘power’ untuk mendapatkan nilai & vektor eigen terbesar. Langkah-langkah: 1. Rumuskan fenomena fisik ke bentuk [A]{u}=λ{u} 2. Prakirakan vektor awal {u}0≠{0} 3. Hitung vektor baru {u}1=[A]{u}0 4. Faktorkan koefisien terbesar {u}1’=λ1{u}1 5. Kembali ke langkah 3 dengan {u}0={u}1’ sampai {u}1≈{u}0 Contoh: Tentukan nilai & vektor eigen dari matrik ⎡30 6 5 ⎤ A = ⎢⎢ 6 30 9 ⎥⎥ ⎢⎣ 5 9 30⎥⎦ Vektor awal: {u}t0 = {0, 0, 1}
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 58 Jack la Motta
Matrik
Buku kuliah
⎧0⎫ ⎧5⎫ ⎧0.166⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ [ A]⎨0⎬ = ⎨ 9 ⎬ = 30⎨0.300⎬ ⎪1⎪ ⎪ ⎪ ⎪1.000 ⎪ ⎩ ⎭ 0 ⎩30⎭1 ⎩ ⎭1 ⎧0.166⎫ ⎧11.780 ⎫ ⎧0.351⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ [ A]⎨0.300⎬ = ⎨18.996 ⎬ = 33.530⎨0.566⎬ ⎪1.000 ⎪ ⎪33.530⎪ ⎪1.000 ⎪ ⎩ ⎭1 ⎩ ⎭2 ⎩ ⎭2 ⎧0.351⎫ ⎧18.926 ⎫ ⎧0.514⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ [ A]⎨0.566⎬ = ⎨29.886⎬ = 36.849⎨0.811⎬ ⎪1.000 ⎪ ⎪ ⎪ ⎪1.000 ⎪ ⎩ ⎭ 2 ⎩36.849 ⎭ 3 ⎩ ⎭3 ⎧0.514⎫ ⎧25.286⎫ ⎧0.634⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ [ A]⎨ 0.811⎬ = ⎨36.414 ⎬ = 39.869⎨0.913⎬ ⎪1.000 ⎪ ⎪ ⎪ ⎪1.000 ⎪ ⎩ ⎭ 3 ⎩39.869 ⎭ 4 ⎩ ⎭4
dst … akhirnya didapat
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
⎧0.800⎫ ⎪ ⎪ λ1 = 43,49 dan {u}1 = ⎨1.000 ⎬ ⎪0.965⎪ ⎩ ⎭
Contoh aplikasi: Suatu elemen berbentuk piramida dalam suatu benda yang dikenai gaya-gaya luar. Gaya normal dan geser yang sejajar sumbu-sumbu koordinat telah diketahui, maka diinginkan bidang patahan yang mungkin terjadi dan besarnya gaya normal yang bekerja pada bidang itu. Kesetimbangan gaya-gaya dalam Gambar 14 sebagai berikut: 1 1 1 ∑ Fx = 0 ⇒ 2 dydzσ x + 2 dxdzτ yx + 2 dxdyτ zx + Rx = 0 ldAσ x + mdAτ yx + ndAτ zx + ldAσ = 0 lσ x + mτ yx + nτ zx = −σl
1 1 1 dydxτ xy + dxdzσ y + dxdyτ zy + R y = 0 2 2 2 lτ xy + mσ y + nτ zy = −σm 1 1 1 ∑ Fz = 0 ⇒ 2 dydzτ xz + 2 dxdzτ yz + 2 dxdyσ z + Rz = 0 lτ xz + mτ yz + nσ z = −σn
∑F
y
=0⇒
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 59 Jack la Motta
Matrik
Buku kuliah
z τxz
Rn
σx
τyz
dz
σy τyx
Rz= Rnn Rn=σdA
τxy
x
dx
Rx= Rnℓ
dy τzx
τzy y
σz
Ry= Rnm
luas bidang miring = dA ½dydz = dA cos(N,x) = ℓdA ½dxdz = dA cos(N,y) = mdA ½dxdy = dA cos(N,z) = ndA Gambar 14 Gaya-gaya yang bekerja pada struktur
Sekarang dalam bentuk matrik: ⎡σ x τ yx τ yz ⎤ ⎧ l ⎫ ⎧l ⎫ ⎢ ⎥⎪ ⎪ ⎪ ⎪ ⎢τ xy σ y τ zy ⎥ ⎨m⎬ = −σ ⎨m⎬ ⎪n⎪ ⎢τ xz τ yz σ z ⎥ ⎪⎩ n ⎪⎭ ⎩ ⎭ ⎣ ⎦ atau [A]{u} = λ{u} ← problem nilai dan vektor eigen, dengan:
• [Α] matrik yang elemennya terdiri dari gaya geser & normal sejajar sumbu D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
koordinat ← diketahui.
• {u} vektor yang elemennya terdiri atas cosinus sudut bidang patahan dengan sumbu kordinat ← dicari.
• λ skalar yang merupakan/menyatakan gaya normal/tegangan normal yang bekerja pada bidang patahan ← dicari.
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 60 Jack la Motta
7.
Bab PERSAMAAN DIFFERENSIAL BIASA
Dalam bidang teknik sering dijumpai persamaan suatu fenomena alam yang dinyatakan dalam persamaan diferensial biasa (PDB) Contoh:
♦ Problem nilai awal:
y’ = f(x,y) dengan y(x0) = Y0
♦ Problem nilai batas:
y” = g(x,y,y’) dimana a<x
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
dengan A & B adalah matrik 2x2 dan γ1 & γ2 konstanta yg telah diketahui. Taylor series Taylor mengatakan bahwa suatu fungsi dengan sifat tertentu dapat dinyatakan sebagai h h2 h n (n ) f ( x0 ) + ... f ( x ) = f ( x0 ) + f ′( x0 ) + f ′′( x0 ) + ... + 1! 2! n! atau h h2 h n (n ) ′ ′ ′ y 0 + ... h + x − x0 y(x ) = y0 + y0 + y 0 + ... + 1! 2! n! Deret ini akan digunakan dalam bab ini 1 Contoh: y ′ − y = 0, y (0) = 1 2 dy 1 dy 1 = y → = dx dx 2 y 2 Secara analitis 1 1 ∫ y dy = 2 ∫ dx x 1 Jadi 1n y = x = c → y = e 2 + c 2 0 y (0) = 1 → 1 = e 2 + c → c = 0 Jika Maka secara analitis y = e
x
2
Dengan deret Taylor dapat diselesaikan sbb:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 61 Jack la Motta
Persamaan Differensial Biasa
1 y 2 1 y ′′ = y ′ 2 Pers. Asli 1 y ′′′ = y ′′ 2 1 y ′′′' = y ′′′ 2 1 1 Jadi y (x ) = 1 + x + x 2 2 8
Buku kuliah
y (0 ) = 1
y′ =
x0 = 0
y ′(0) =
1 2 1 y ′′(0 ) = 4 1 y ′′′(0) = 8 1 3 1 4 + x + x + ... 48 384
y
y = ex/2
y = 1+
1 1 1 3 x + x2 + x 2 8 384
y = 1+
1 1 x + x2 2 8
y = 1+
1 x 2 x
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Gambar 15 Penyelesaian dengan Metoda Euler Cara numeris untuk menyelesaikan problem nilai awal adalah diferensial hingga. Pada metoda diferensi hingga penyelesaian pendekatan didapat pada titik-titik hitung x0 < x1 < x 2 < ... < x n < ... dan nilai pendekatan pada setiap x n diperoleh dengan menggunakan nilai-nilai yg didapat sebelumnya. Ditinjau PDB:
y ′ = f ( x, y ), y ( x0 ) = Y0
Penyelesaian sesungguhnya ditulis Y(x), sehingga pers. diatas menjadi: Y’(x) = f(x, Y(x)), Y(x0) = Y0
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 62 Jack la Motta
Persamaan Differensial Biasa
Buku kuliah
Penyelesaian pendekatannya ditulis y(x) dan nilai y(x0), y(x1), …, y(xn), … atau ditulis sebagai y0, y1, …, yn, … sebuah pias h 〉 0 digunakan untuk mendefinisikan titik-titik hitung x j = x0 + jh
j = 0,1, ...
Jika akan diadakan perbandingan penyelesaian pendekatan untuk beberapa nilai h, maka yh(x) digunakan untuk menyatakan y(x) dengan pias h.
7.1. Metoda Euler Dengan deret Taylor hitung Y(xn+1) dengan menggunakan Y(xn) h2 Y ( x n +1 ) = Y ( x n ) + h Y ′(x n ) + Y ′′(ξ n ) 2 dengan x n ≤ ξ ≤ x n +1 Rumus Euler menjadi: dengan
y n +1 = y n + hf ( x n , y n ) y 0 = Y0
n = 0,1,2, ...
Kesalahan diskritisasi adalah
Y ( x n +1 ) − y n +1 =
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Contoh:
h2 Y ′′(ξ n ) 2
PDB y ′ = y, y (0 ) = 1 → Y ( x ) = e x
Rumus Euler uth h = 0.2 y n +1 = y n + 0.2 y n = 1.2 y n y1 = 1.2 y 0 = 1.2 × 1 = 1.2
y 2 = 1.2 y1 = 1.2 × 1.2 = 1.44 Jika ditabelkan: h 0.2
x 0.4 0.8 1.2 1.6 2.0
yh(x) 1.44000 2.07360 2.98598 4.29982 6.19174
0.1
0.4 0.8 1.2
1.46410 2.14356 3.13843
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
Yh(x) Yh(x)- yh(x) 1.49182 0.05182 2.22554 0.15194 3.32012 0.33414 4.95303 0.65321 7.38906 1.19732 1.49182 2.22554 3.32012
0.02772 0.08198 0.18169
hal. 63 Jack la Motta
Persamaan Differensial Biasa
h
Buku kuliah
x 1.6 2.0
yh(x) 4.59497 6.72750
Yh(x) Yh(x)- yh(x) 4.95303 0.35806 7.38906 0.66156
Perhatikan bahwa kesalahan menurun ½ dari nilai pertama karena h dikecilkan ½ kali
7.2. Metoda ‘Multi–Step’ Secara umum rumus langkah majemuk dapat ditulis sebagai:
y n +1 = ∑ a j y n − j + h ∑ b j f (x n − j , y n − j ) p
p
j =0
j = −1
n≥ p
Koefisien a 0 , ..., a p , b−1 , b0 , ..., b p adalah suatu konstanta, dan p ≥ 0 Jika a p ≠ 0 dan b p ≠ 0 , metoda ini disebut metoda langkah (p+1), karena (p+1) nilai pendekatan sebelumnya digunakan untuk menghitung y n +1 . Nilai y1 , ..., y p harus dihitung dengan cara lain. Metoda Euler adalah metoda langkah tunggal karena p = 0, dan a 0 = 1 , b−1 = 0 , b0 = 1 . Jika b−1 = 0 maka Y( xn +1 ) hanya terdapat pada ruas kiri, sehingga rumusnya disebut rumus eksplisit. Jika b−1 ≠ 0 , maka Y( xn +1 ) terdapat diruas kanan maupun kiri, sehingga disebut D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
rumus implisit. Koefisien a j dan b j dapat dihitung dari p
p
j =0
j =0
∑ a j = 1 − ∑ ja j + p
p
j =0
j = −1
p
∑b j = −1
j
i i −1 ∑ (− j ) a j + i ∑ (− j ) b j
=1 =1
i = 2, ..., m
Rumus terakhir menjamin bahwa Y ( x ) dapat diderivikasikan
(m + 1)
kali. Jika a 0 = 0 , a1 = 1 , b−1 = 0 , b0 = 2 , maka didapat rumus untuk metoda titik tengah
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
y n +1 = y n −1 + 2hf ( x n , y n )
n≥1
hal. 64 Jack la Motta
Persamaan Differensial Biasa
Buku kuliah
Merupakan metoda langkah ganda yang explisit. 1 Jika a 0 = 1 , b−1 = b0 = , maka didapat rumus trapesium yang implisit dan 2 merupakan metoda langkah tunggal: 1 y n +1 = y n + h[ f (x n , y n ) + f ( x n +1 , y n +1 )] , n ≥ 0 2
7.2.1. Metoda Trapesium Metoda ini dapat pula dijabarkan dari: Y’(t) = f(t, Y(t)) Diintegrasikan dari [x n , x n +1 ] x n +1
x n +1
xn
xn
∫ Y ′(t )dt = ∫ f (t , Y (t ))dt
1 h3 h[ f ( x n , Y ( x n )) + f ( x n +1 , Y ( x n +1 ))] − Y ′′′(ξ n ) x n ≤ ξ n ≤ x n +1 2 12 sehingga pendekatannya menjadi 1 y n +1 = y n + h[ f ( x n , Y ( x n )) + f ( x n +1 , Y (x n +1 ))] 2 Karena rumusnya implisit, maka yn+1 dapat dihitung dengan iterasi, jadi secara Y ( x n +1 ) − Y ( x n ) =
umum:
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
y n( +j +11) = y n +
[
]
1 h f ( x n , y n ) + f (x n +1 , y n( +j )1 ) 2
j=0, 1, …
Langkah-langkah hitungan: 1. xn,yn telah diketahui/dihitung pada langkah sebelumnya 2. y n(0+)1 diprakirakan dgn rumus eksplisit, misalkan
y n(0+)1 = y n + hf ( x n , y n )
3. y n(0+)1 dimasukkan kedalam ruas kanan sehingga y n(1+)1 dapat dihitung 4. langkah 2 diulang s/d ketelitian yang dikehendaki Walaupun secara umum dapat diselesaikan dengan iterasi, tetapi mungkin dapat diselesaikan dengan cara lain atau bahkan tanpa iterasi tergantung dari ƒ(x,y). Contoh:
y’= y
Karena
ƒ(x,y) = y, maka y n +1 = y n +
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
y(0) = 1
1 h( y n + y n +1 ) 2
hal. 65 Jack la Motta
Persamaan Differensial Biasa
Buku kuliah
h 2y = h n 1− 2 1+
y n +1
Untuk
1.1 y n = 1.2222 y n 0.9 y1 = 1.2222 y 0 = 1.2222 × 1
h = 0.2 → y n +1 =
y 2 = 1.2222 y1 = 1.2222 × 1.2222 = 1.49383 Jika ditabelkan: x 0.4 0.8 1.2 1.6 2.0
yh(x) 1.49383 2.23152 3.33350 4.97968 7.43878
Yh(x) 1.49182 2.22554 3.32012 4.95303 7.38906
Yh(x)- yh(x) -0.00200 -0.00598 -0.01339 -0.02665 -0.04972
7.3. Metoda Runge-Kutta (RK)
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
Metoda Runge-Kutta merupakan metoda langkah tunggal yang lebih teliti dibandingkan metoda Euler. Semua metoda RK dapat ditulis sebagai: y i +1 = y1 + hΦ( xi , yi , h ) dengan Φ( xi , y i , h ) disebut fungsi penambah
7.3.1. Metoda RK derajat dua dengan
Φ = ak1 + bk 2 k1 = f ( xi ,yi ) k 2 = f ( xi + ph,qhf ( xi ,y i ) + y i )
= f ( xi + ph,qhk1 + y i )
dimana a,b,p,q akan ditentukan kemudian. Ditinjau deret Taylor untuk 2 variabel (x,y):
f ( x + r , y + s ) = f ( x, y ) + rf x ( x, y ) + sf x ( x, y ) + rsf xy ( x, y ) +
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
[
1 2 r f xx ( x, y ) + 2
1 2 3 s f yy (x, y ) + O ( r + s ) 2
]
hal. 66 Jack la Motta
Persamaan Differensial Biasa
Buku kuliah
∴ k 2 = f ( xi ,y i ) + phf x (xi ,y i ) + k1 qhf y ( xi ,y i ) + O(h 2 )
Jadi
y i +1 = y i + hΦ ( xi ,y i ,h )
[
]
= y i + h[af (xi ,y i ) + bf ( xi ,y i )] + h 2 pbf x ( xi ,y i ) + bqf ( xi ,y i ) f y ( xi ,y i )
(A)
Dengan deret Taylor: y i +1
= y i + hy ′(xi ,y i ) +
h2 y ′x ( xi ,y i ) + ... 2!
[
]
h2 y i′′ + y ′y y ′ 2 h2 = y i + hf ( xi ,y i ) + f x (xi ,y i ) + f (xi ,y i ) f y ( xi ,y i ) 2 A = B → a + b = 1, bp = 12 , bq = 12 = y i + hy i′ +
[
(B)
]
Tidak dapat diselesaikan karena hanya ada 3 persamaan dengan 4 bilangan anu yaitu a,b,p,q. Biasanya nilai b adalah ½ atau 1.
♦ Untuk b = ½, a = ½, p =q = 1 Euler utk y
y i +1
i +1 6474 8 h = y i + [ f ( xi , y i ) + f {xi + h, y i + hf ( xi , y i )}] 424 3 1444424444 3 2 1 slope di xi
(C)
slope di xi +1 dihitung dgn metoda Euler
Jadi metoda RK dapat dipandang sebagai metoda ‘predictor-corrector’
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
1. Langkah predictor: a. prakirakan yi +1 dengan metoda Euler b. slope (y’) dititik (xi +1 , yi +1 ) adalah f ( xi +1 , y i +1 ) 2. Langkah corrector: a. hitung slope (y’) dititik (xi, yi) yaitu f(xi, yi) b. hitung slope rerata = (slope 1.b + slope 2.a)/2 c. hitung yi+1 = yi + h x hasil 2.b
♦ Untuk b = 1, a = 0, p = q = ½ y i +1 = y i + hf (xi + 12 h, y i + 12 hf ( xi + y i ))
7.3.2. Metoda RK berderajat tiga y i +1 = y i + 16 h[k1 + 4k 2 + k 3 ] k1 = f ( x i , y i )
k 2 = f (xi + h 2 , y i 12 hk1 )
k 3 = f ( xi + h, y i + 2hk 2 − hk1 )
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 67 Jack la Motta
Persamaan Differensial Biasa
Buku kuliah
7.3.3. Metoda RK berderajat empat 7.3.3.1.
Metoda Pertama
(k1 + 2k 2 + 2k 3 + k 4 ) f ( xi , y i ) f (xi + 12 h, y i + 12 hk1 ) f (xi + 12 h, y i + 12 hk 2 ) f (xi + 12 h, y i + 12 hk 3 )
y i +1 = y i + k1 = k2 = k3 = k4 =
7.3.3.2.
h 6
Metoda Kedua
(k1 + 3k 2 + 3k 3 + k 4 ) f ( xi , y i )
y i +1 = y i + k1 =
h 8
k 2 = f (xi + h 3 , y i + 13 hk1 )
k 3 = f (xi + 23 h, y i − 13 hk1 + hk 2 )
k 4 = f ( xi + h, y i + hk1 − hk 2 + hk 3 )
7.3.3.3.
Metoda Ketiga Metoda inilah yang paling banyak digunakan y i +1 = y i + 16 h k1 + 2 − 2 k 2 + 2 − 2 k 3 + k 4
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
[ (
)
k1 = f ( x i , y i )
(
)
)
( )hk )
k 2 = f (xi + 12 h, y i + hk1 )
( = f (x
(
k 3 = f xi + 12 h, y i + 12 − 1 + 2 hk1 + 1 1 − k4
Contoh:
i
+ h, y i −
1 2
(
hk 2 + 1 +
1 2
y (0 ) = 1 ¼ eksak Y ( x ) = e
y′ = y 1 2
x
1 2
]
)hk ) 2
3
2
Y (1) = e 1 2 = 1.648721271 Diselesaikan dengan RK4: f ( x, y ) = 12 Y
k1 = f (x 0 ,y 0 ) = f (0,1) = 12 × 1 =
1 2
k 2 = f (x 0 + 12 h,y 0 + 12 hk1 ) = f (0 + 12 .1,1 + 12 .1. 12 ) = f ( 12 , 54 ) = 12 . 54 =
5 8
( ( ) ( )hk ) [1 + (− 1 + 2 )(1)( ) + (1 − )(1)( )] = 0.64331
k 3 = f x0 + 12 h,y 0 + 12 − 1 + 2 hk1 + 1 − =
1 2
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
1 2
1 2
1 2
5
1 2
2
8
hal. 68 Jack la Motta
Persamaan Differensial Biasa
( [1 −
Buku kuliah
k 4 = f x 0 + h,y 0 − =
1 2
1 2
(1)( 8 ) + (1 +
1 2
( )hk ) )(1)(0.64331)] = 0.828125
hk 2 + 1 + 1 2
5
1 2
3
[ + (2 + 2 )(0.64431) + 0.828125] = 1.6484375
y (1) = y (0 ) + 16
1 2
7.4. Metoda ‘Predictor-Corrector’ Metoda langkah majemuk berdasarkan rumus integrasi. Secara umum metoda ini mengintegrasi PDB pada interval [xi-k, xi+1] sebagai berikut: y ′ = f ( x, y ) y i +1
xi + 1
yi − k
xi − k
∫ dy =
∫ f (x, y )dx
y i +1 = y i − k +
xi +1
∫ f (x, y )dx
xi − k
r
f(x, y) didekati polinomial derajat r, Φ ( x ) = ∑ a j x j j =0
Integrasi terbuka dan beda terbagi mundur menghasilkan:
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
♦ untuk k = 0, r = 3 y i +1 = y i +
h 24
(55 f i − 59 f i −1 + 37 f i −2 − 9 f i −3 )
♦ untuk k = 1, r = 1 y i +1 = y i −1 + 2hf i
O(h5) O(h3)
♦ untuk k = 3, r = 3 y i +1 = y i −3 + 43 h(2 f i − f i −1 + 2 f1+ 2 )
(I)
♦ untuk k = 5, r = 5 y i +1 = y i −5 + 103 h(11 f i − 14 f i −1 + 26 f i − 2 − 14 f i −3 + 11 f i − 4 )
O(h5)
O(h7)
Integrasi tertutup dan beda terbagi mundur menjadi:
♦ untuk k = 0, r = 3 y i +1 = y i +
h 24
♦ untuk k = 5, r = 5 y i +1 = y i −1 +
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
(9 f i +1 + 19 f i − 5 f i −1 +
h 3
( f i +1 + 4 f i +
f i −1 )
f i−2 )
O(h5) (II)
O(h5)
hal. 69 Jack la Motta
Persamaan Differensial Biasa
♦ untuk k = 3, r = 5 y i +1 = y i −3 +
Buku kuliah
2h 45
(7 f i +1 + 32 f1 + 12 f i −2 + 7 f i −3 )
O(h7)
Kesulitan metoda langkah majemuk adalah pada saat permulaan y i − k belum terhitung sehingga harus harus dihitung dengan cara lain, misalnya metoda Euler. Dari hasil di atas tampak bahwa integrasi terbuka memberikan rumus eksplisit; sehingga hitungan tidak menggunakan iterasi. Integrasi tertutup menghasilkan rumus implisit, sehingga membutuhkan iterasi. Walaupun menggunakan iterasi, integrasi tertutup lebih disukai karena ketelitiannya lebih tinggi. Contoh: 14 5 ( 4 ) h f (ξ ) 15 1 5 (4) (II) kesalahan: − h f (ξ 90
Pada
xi −3 ≤ ξ ≤ xi +1
(I) kesalahan:
)
xi −3 ≤ ξ ≤ xi +1
♦ Rumus Adam-Bashforth (eksplisit) 1. Yn +1 = Yn + hYn′ + 12 h 2Y ′′(ξ n )
[3Y ′ − Yn′−1 ] + 125 h 3Y (3) (ξ n ) Yn +1 = Yn + 12h [23Yn′ − 16Yn′−1 + 5Yn′− 2 ] + 83 h 4Y (4 ) (ξ n ) 251 5 (5 ) Yn +1 = Yn + 24h [55Yn′ − 59Yn′−1 + 37Yn′− 2 − 9Yn′−3 ] + 720 h Y (ξ n )
2. Yn +1 = Yn + 3. D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
4.
h 2
♦ Rumus Adam – Moulton (implisit) 1. Yn +1 = Yn + hYn′+1 − 12 h 2Y ′′(ζ n )
[Yn′+1 + Yn′ ] − 121 h 3Y (3) (ζ n ) Yn +1 = Yn + 12h [5Yn′+1 + 8Yn′ − Yn′−1 ] − 121 h 4Y (4 ) (ζ n ) 19 Yn +1 = Yn + 24h [9Yn′+1 + 19Yn′ − 5Yn′−1 + Yn′− 2 ] − 720 Y (5 ) (ζ n )
2. Yn +1 = Yn + 3. 4.
h 2
7.4.1. Algoritma ‘Predictor-Corrector Rumus A-M membutuhkan penyelesaian iterasi, sedangkan rumus A-B tidak, tetapi A-M ketelitiannya lebih tinggi. Algoritma ‘predictor-corrector’ berusaha menggabungkan keuntungan kedua rumus diatas, sebagai berikut:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 70 Jack la Motta
Persamaan Differensial Biasa
Buku kuliah
1. Gunakan rumus A-B untuk memperkirakan Yn +1 (predictor; rumus A-B) 2. Hitung Yn +1 memakai rumus A-M tanpa iterasi dengan memakai Yn +1 nilai dari a (corrector; rumusA-M) Contoh: 1. A-B-M derajat 4: a. Predictor A-B
y n +1 = y n + h [55 f n − 59 f 24
+37 f
−9 f
n −1
n−2
n −3
]
b. Corrector A-B
y n +1 = y n + h [9 f n +1 − 19 f n−5 f 24
+f
n −1
n−2
]
2. Rumus Milne derajat 4: a. Predictor
y n +1 = y n −3 + 4h (2 f n − f n −1 + 2 f n − 2 ) 3 b. Corrector
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
y n +1 = y n −1 + h ( f n +1 + 4 f n + f n −1 ) 3 Contoh penggunaan: dy 1 − y=0 dx 2 Penyelesaiaan:
y(0)=1 Æ hitung y(1) = ?
y ′ = 12 y → f ( x, y ) = 12 y
Catatan: solusi eksak Y = e
x
2
Digunakan h = 0.25 Euler y n +1 = y n + hf n i n-3 n-2 n-1 n n+1 A – B –4:
xi 0.00 0.25 0.50 0.75 1.00
yi 1.0000 1.1250 1.2656 1.4238 1.6018
f(xi,yi) 0.5000 0.5625 0.6328 0.7119 0.8009
[55 f n − 59 f n−1 + 37 f n−2 − 9 f n−3 ] .24 = y n + 024 [55(0.7119) − 59(0.6328) + 37(0.5625) − 9(0.50)]
y n +1 = y n +
h 24
= 1.4238 + 0.18887 = 1.6127
(predictor)
A-B-M-4:
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 71 Jack la Motta
Persamaan Differensial Biasa
Buku kuliah
f n +1 = 12 y n +1 = 12 × 1.6127 = 0.8063
[9 f n+1 + 19 f n − 5 f n−1 + f n−2 ] .25 [9(0.8063) + 19(0.7119) − 5(0.6328) + 0.5625] = 1.4238 + 024
y n +1 = y n +
h 24
= 1.4238 + 0.1894 = 1.6132
A-M-4: Iterasi 2:
f n +1 = 12 × 1.6132 = 0.8066
.25 y n +1 = 1.4238 + 024 [9(0.8066) + 19(0.7119) − 5(0.6328) + 0.5625]
f n +1 Iterasi 3: Iterasi 9:
(predictor - corrector)
= 0.613216 = 0.8066
y n +1 = 1.613217
M f n +1 = 0.80661 y n +1 = 1.613217
Tampak bahwa algoritma ‘predictor-corrector sudah mencukupi dibandingkan dengan iterasi. Tabel hasil
D:\My Documents\Publikasi\Metoda Numerik\Metoda Numerik.doc (2028 Kb)
x 0,00 0,25 0,50 0,75 1,00
ex/2 1,0000 1,13315 1,28403 1,45499 1,64872
Euler 1,0000 1,1250 1,2656 1,4238 1,6018 -2,846%
A-B-4 1,6127 -2,185%
A-B-M-4 1,6132 -2,154%
Untuk latihan: hitung y(1) = ? dengan h = 0,125, bandingkan dengan hasil h = 0,25
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 72 Jack la Motta
DAFTAR PUSTAKA
Carnahan, Brice, H.A. Luther, James O. Wilkes, Applied Numerical Methods, John Wiley & Sons, New York, 1969. Spiegel, R. Murray, Theory and Problems of Statistics, Schaum’s Outline Series, McGraw-Hill International Book Company, Singapore, 1981. Al-Khafaji, Amir Wahdi, John R.Tooley, Numerical Methods in Engineering Practice, Holt, Rinehart and Winston, Inc., New York, 1986. Anonim, fx–7000G Owner’s Manual, CASIO® Atkinson, Kendall E., An Introduction to Numerical Analysis, John Wiley & Sons, New York, 1989. James, M.L., G.M. Smith, J.C. Wolford, Applied Numerical Methods for Digital Computation with Fortran and CSMP, 2nd Edition, Harper International Edition, New York, 1977.
Metoda Numerik Ir. Djoko Luknanto, M.Sc., Ph.D.
hal. 73 Jack la Motta