Penyelesaian numerik Persamaan differensial Poisson 2D, dengan syarat batas Dirichlet, Neumann, Robin dilengkapi Analisis Real dan Program Matlab

SOLUSI NUMERIK PERSAMAAN DIFFERENSIAL PARSIAL ELIPTIK: PERSAMAAN POISSON DAN LAPLACE 2D

Tesis Untuk memenuhi sebagian persyaratan Mencapai derajat sarjana S2 Program Studi Matematika Jurusan Ilmu-ilmu Mipa

Diajukan oleh : La Ode Muhammad Umar Reky Rahmat Rikardo 20553 / I-4 / 1623 / 03

Kepada PROGRAM STUDI MATEMATIKA SEKOLAH PASCASARJANA UNIVERSITAS GADJAH MADA YOGYAKARTA 2005

ii

ii

PERNYATAAN

Dengan ini saya menyatakan bahwa dalam tesis ini tidak terdapat karya yang pernah diajukan untuk memperoleh gelar kesarjanaan disuatu perguruan tinggi dan sepanjang pengetahuan saya tidak terdapat karya atau pendapat yang pernah ditulis atau diterbitkan oleh orang lain, kecuali yang secara tertulis diacu dalam naskah ini dan disebutkan dalam daftar pustaka

Yogyakarta, 31 Desember 2005

La Ode Muhammad Umar Reky R.R

iii

Pelajarilah ilmu dan ajarlah manusia dan rendahkan diri kepada guru-gurumu serta berlaku lemah lembutlah terhadap murid-muridmu. ( HR. Thabrani )

Untuk buah hati & penyejuk mataku , La Ode Muhammad Fayyad Faid Fidai Wa Ode Naurah Najlaa’ Kamila Mutmainnah

PRAKATA

Assalamualaikum Wr. Wb.

Puji syukur kehadirat Allah SWT, atas segala berkat,rahmat dan hidayah-Nya, penulis dapat menyelesaikan tugas akhir ini dengan segenap tenaga dan pikiran. Dalam penyelesaian tugas akhir ini, penulis banyak menerima bimbingan dan bantuan dari berbagai pihak, maka dalam kesempatan ini penulis mengucapkan terima kasih kepada :

1.

Ibu Dr.Lina Aryati, MS., selaku Dosen Pembimbing yang telah meluangkan waktu yang sangat berharga dan membimbing penulis tanpa pamrih,

semoga Allah

mencatatnya menjadi amal ibadah dan diberikan balasan yang berlipat ganda. 2.

Bapak Dr. Widodo, MS., selaku Ketua Jurusan Matematika, UGM.

3.

Ibu Dr. Ch. Rini Indrati, M.Si., Bapak Dr. Ari Suparwanto, M.Si, Bapak Dr. Widodo, MS.; penulis mengucapkan terimakasih yang sebanyak-banyaknya atas arahan, petunjuk untuk kesempurnaan tulisan ini.

4.

Ibu Prof. DR. Sri Wahyuni, MS., selaku Ketua Pengelola Program Pascasarjana Jurusan matematika, FMIPA, UGM dan pengasuh matakuliah Aljabar.

5.

Bapak Prof. Dr. Soeparna Darmawijaya, selaku pengasuh matakuliah Analisis Real.

6.

Ibu Emiliana Sunaryani Y., selaku Tata Usaha Jurusan Matematika dan Pelaksana Akademik Program S2 Matematika FMIPA-UGM.

7.

Seluruh staf dosen Jurusan Matematika, FMIPA, UGM, yang telah memberikan banyak bekal ilmu, semoga Allah

mencatanya menjadi amal ibadah dan

diberikannya balasan yang berlipat ganda.

v

8.

Bapak Abdul Rahman, S.Si, M.Si, (Alumni Fisika UGM angk. 2002, dosen Jurusan Fisika, Fmpa, UNCEN, Papua); Penulis mengucapkan terimakasih yang tidak terkira atas bantuanya dalam perintisan program komputer pada tulisan ini.

9.

Bapak Supiyanto,S.Si,M.Kom.,(Ilkom UGM 2003, dosen Jur. Matematika, FMIPA, UNCEN, Papua), Bapak Ruliansyah, S.kom,M.kom., (Ilkom UGM 2003, Dosen STIMIK, Palembang), Bapak Hidayat, S.kom,M.kom., (Ilkom UGM 2002, Dosen Jur. Matematika, IKIP, Gorontalo), Pak Chan (Ilkom UGM 2002, staf PEMDA Sumatera Barat), Bapak Rais,S.Si, M.Si. (Matematika UGM 2003, dosen Universitas Sulawesi Tenggara), Pak La Suha (Matematika UGM 2003, dosen UNPATI, Ambon), Yahya Khairun, Spd, M.Si, (Matematika UGM 2002, Dosen matematika UNHAIR, Ternate), Pak Karman Lanani (Matematika UGM 2002, dosen matematika UNHAIR, Ternate), Pak Nursalam (Matematika UGM 2003, dosen matematika Universitas Lambung Mangkurat), dan rekan-rekan di program S2 ilmu komputer angkatan 2003, terimakasih atas bantuan dan motifasinya selama bersama di Universitas Gadjah Mada.

10.

Sdr. Muhammad Hidayat

(S1, Jurusan kimia, FMIPA, UGM, angkatan 2003);

Penulis mengucapkan banyak terimakasih

atas fasilitas dan software untuk

visualisasi dalam tulisan ini. Dengan segala keyakinan dan kerendahan hati penulis menyadari bahwa penulisan ini masih jauh dari kesempurnaan. Oleh karena itu kritik dan saran yang membangun senantiasa sangat penulis harapkan dari para pembaca. Akhirnya, harapan penulis semoga tulisan ini dapat bermanfaat bagi pihak-pihak yang membutuhkan. Wassalamualaikum Wr. Wb

Penulis

vi

DAFTAR ISI

halaman i

HALAMAN JUDUL LEMBAR PENGESAHAN

ii

HALAMAN PERNYATAAN

iii

HALAMANPERSEMBAHAN

iv

PRAKATA

v

DAFTAR ISI

vii

DAFTAR GAMBAR

xi

ARTI LAMBANG DAN SINGKATAN

xiv

INTISARI

xiii

ABSTRACT

xiv

BAB I.

PENDAHULUAN

1

I.1. Latar Belakang Masalah

1

1.2 Manfaat Penelitian

2

1.3 Tujuan Penelitian

2

1.4 Tinjaun Pustaka

3

1.5 Metode Penelitian

4

1.6 Sistematika Penulisan

4

BAB II.

LANDASAN TEORI

6

2.1 Persamaan Differensial Poisson dan Laplace 2D

6

2.2 Bentuk Umum Syarat Batas

7

2.3 Deret Taylor

7

vii

2.4 Deret Taylor Fungsi Dua Peubah untuk Mendapatkan Rumus Beda Hingga Turunan Parsial

8

2.5 Barisan Bilangan Real

12

2.6 Norm

19

2.6.1 Norm Vektor

19

2.6.2 Matrik Transpose Conjugate

20

2.6.3 Norm Matriks

21

BAB III. MENYELESAIKAN PD POISSON DAN LAPLACE 2D DENGAN METODE BEDA HINGGA DAN METODE ITERASI SOR

22

3.1 Hubungan PD Poisson dan Laplace 2D, Metode Beda Hingga dan Barisan Bilangan Real

22

3.2 Kekonvergenan Metode Iterasi

24

3.3 Skema Iterasi Metode SOR

29

3.4 Parameter Relaksasi Metode SOR

30

3.5 Menyelesaikan Persamaan Differensial Poisson dan Laplace 2D Menggunakan Metode Beda Hingga dan Metode Iterasi SOR 3.6 Mengatur Nilai Awal Iterasi 3.6.1

Mengatur Nilai Awal pada Batas domain

3.6.2

Skema Untuk Mengatur Nilai Awal Pada Grid Bagian Dalam

3.6.3

Domain

34 42 42

43

Urutan Langkah Mengatur Nilai Awal Iterasi

3.7 Bentuk Matrik PD Poisson dan Laplace 2D

viii

45 46

3.7.1

PD Poisson dan Laplace 2D dengan syarat batas domain seluruhnya bertipe Dirichlet

3.7.2

48

PD Poisson dan Laplace 2D dengan syarat batas tipe Dirichlet terdiri dari dan bertipe Neumman

3.7.3

53

PD Poisson dan Laplace 2D dengan syarat batas domain seluruhnya bertipe Neumman

3.7.4

55

PD Poisson dan Laplace 2D dengan syarat batas

terdapat

tipe Robin

57

3.8 Sifat skema beda hingga PD Poisson 2D dan Laplce 2D

59

3.9 Urutan Penyelesaian

61

BAB IV. PENGGUNAAN METODE DAN APLIKASINYA

62

4.1 PD Poisson 2D dengan Syarat Batas Seluruhnya Bertipe Dirichlet 4.1.1 PD Poisson 2D dengan Nilai Awal Iterasi u ( 0)

=

0

63 64

4.1.2 PD Poisson 2D dengan Nilai Awal Iterasi u ( 0 ) dimodifikasi

67

4.2 PD Laplace 2D dengan Syarat Batas Seluruhnya bertipe Dirichlet 4.2.1 PD Lplace 2D dengan Nilai Awal Iterasi u ( 0 )

=

0

69 70

4.2.2 PD Laplace 2D dengan Nilai Awal Iterasi u ( 0 ) dimodifikasi

72

4. 3 PD Poisson 2D dengan Domain Mempunyai Batas bertipe Dirichlet dan Batas Lainnya Bertipe Neumann

ix

74

4.4 Contoh Kasus, Metode SOR Dapat Menghasilkan Solusi yang Tidak Tunggal

75

4.5 Persamaan Differensial Laplace dengan Batas Domain Terdapat Tipe Dirichlet, Neumman dan Robin

77

4.6 Konduksi Panas Steady State

82

4.7 Transformasi Koordinat

87

BAB V PENUTUP

91

5.1 Kesimpulan

91

5.2 Masalah lanjutan

91

LAMPIRAN

92

DAFTAR PUSTAKA

x

DAFTAR GAMBAR

halaman Gambar 3.1. (a) Ilustrasi kekonvergenan secara toritis dan secara numeris

(b) Ilustrasi barisan Cauchy untuk menghentikan iterasi

22

Gambar 3.2. Domain dan partisinya

34

Gambar 3.3. Vektor normal satuan

36

Gambar 3.4. Domain dengan syarat batas Robin atau Neumann

37

Gambar 3.5. Mengatur tebakan awal pada batas bertipe Robin

42

Gambar 3.6. Skema untuk menentukan nilai awal iterasi

44

Gambar 3.7. Sketsa mencari tebakan awal pada level grid ke-1 menggunakan

skema I atau skema II

45

Gambar 3.8. Sketsa mencari tebakan awal pada level grid ke-2 menggunakan

(a) skema I

(b) skema II

46

Gambar 3.9. Domain dan partisinya

47

Gambar 3.10. Ilustrasi letak nilai-nilai u pada domain

51


53


57

Gambar 3.13. Ilustrasi letak nilai-nilai u dan koefisiennya

59

Gambar 3.14. Skema penyelesaian masalah

61

Gambar 4.1. Grafik contour u ( x, y ) = x 2

+ y

2

63

Gambar 4.2. Grafik solusi numerik

64

Gambar 4.3. Grafik error global

66

xi

Gambar 4.4. Pengaruh parameter SOR (

) terhadap banyak iterasi untuk

mencapai kekonvergenan ( v 0 ) menggunakan toleransi ε = 10 −6 Gambar 4.5. Grafik solusi numerik

66 59

Gambar 4.6. Grafik error global dengan parameter SOR

= 1.2

60

Gambar 4.7. Grafik error global dengan parameter SOR

=1

60

Gambar 4.8. Grafik u ( x, y ) = x

2

− 3 xy − y

2

+ x − 2 y

61


62

Gambar 4.10. Grafik error global

63

Gambar 4.11. Pengaruh parameter SOR (

) terhadap banyak iterasi untuk

mencapai kekonvergenan (n) menggunakan toleransi ε = 10 −6 Gambar 4.12. Grafik solusi numerik

64 62

Gambar 4.13. Grafik error global dengan

= 1.2 ,

grid 33x33

65

Gambar 4.14. Grafik error global dengan

= 1,

grid 2049x2049

65

Gambar 4.15. Grafik error solusi numerik

67

Gambar 4.16. Grafik solusi numerik dengan u ( 0 )

=

Gambar 4.17. Grafik solusi numerik dengan u ( 0 )

= 10,

0, ω = 1.5 ω = 1. 333

68 69


70

Gambar 4.19. Sketsa solusi numerik

71

Gambar 4.20. Grafik error relatif

71

Gambar 4.21. Error teredam

72

Gambar 4.22. Error relatif tidak teredam

73

Gambar 4.23. Suhu pada lempengan logam

75

Gambar 4.24. Sketsa suhu pada lempengan logam

76

xii

Gambar 4.26. Aproksimasi suhu dan bentuk lelehan

78

Gambar 4.27. Error relatif suhu pada lempengan logam

79

xiii

ARTI LAMBANG DAN SINGKATAN

PD

:

Persamaan Differensial

2D

:

Dua dimensi

Ω

:

Domain

i

:

Indeks pada arah horizontal

j

:

Indeks pada arah vertikal

I

:

Indeks maksimum pada arah horizontal

I

:

Matrik identitas

J

:

Indeks maksimum pada arah vertikal

h

:

Ukuran kisi pada arah horizontal

k

:

Ukuran kisi pada arah vertikal

v

:

Nomor iterasi

0

:

Vektor nol

0

:

Matriks nol

e (v )

:

Error global pada iterasi v

ε

:

Nilai toleransi kekonvergenan

v0

:

Banyaknya iterasi untuk mencapai kekonvergenan

ω

:

Parameter metode iterasi SOR

D

:

Matrik diagonal

L

:

Matrik segitiga bawah ( strictly lower triangular )

U

:

Matrik segitiga atas ( strictly upper triangular )

M

:

Matrik iterasi (iterative matrix)

xiv

M ω

:

Matrik iterasi SOR

λ

:

Nilai eigen

M

:

Spectral radius matrik iterasi M

ω

:

Spectral radius matrik iterasi SOR M ω

|| . ||

:

Norm

⇒

:

Jika ... maka ...

⇔

:

Jika dan hanya jika

:

Akhir bukti

▄ , ●

xv

INTISARI

Pada tulisan ini, metode beda hingga dan metode iterasi Succesive Over Relaxation (SOR) dengan memanfatkan syarat batas Dirichlet sebagai nilai awal iterasi digunakan untuk menyelesaikan persamaan differensial Poisson dan Laplace 2D pada domain segiempat. Metode ini dapat meningkatkan akurasi numerik dan effisien. Lebih lanjut, diberikan beberapa contoh dan aplikasinya.

Kata Kunci : Persamaan Poisson dan Laplace, metode beda hingga, SOR, syarat batas Dirichlet,

xvi

ABSTRACT

In this paper, the finite difference methods and sucesive over relaxation (SOR) using Dirichlet boundary condition for initial guess are used to determine the solutions of Poisson and Laplace equation on a rectangular domain. These methods preserve the accuracy and provide the efficiency. Moreover, some applications are presented. Key Words

: Poisson and Laplace equations, finite difference methods, SOR, Dirichlet boundary condition.

xvii

BAB I PENDAHULUAN

I.1. Latar Belakang Masalah

Dalam banyak hal pada bidang terapan, sering ditemukan permasalahan yang melibatkan persamaan differensial Poisson dan Laplace. Bentuk ini biasanya dikaitkan dengan suatu keadaan steady state yaitu keadaan yang tidak bergantung pada waktu. Persamaan differensial Poisson dan Laplace mempunyai aplikasi yang luas pada bidang rekayasa. Persamaan differensial ini dapat dijumpai pada masalah elektrostatik, dinamika fluida, konduksi panas, masalah air tanah dan lain-lain. Umumnya persamaan differensial yang dihadapi pada bidang aplikasi tidak selalu mudah diselesaikan secara analitis. Demikian halnya dengan persamaan differensial Poisson dan Laplace, kesulitannya bergantung pada syarat batas yang diberikan. Karena kendala utama ini, maka penyelesaian menggunakan metode numerik menjadi alternatif penting. Salah satu metode numerik yang sering digunakan adalah metode beda hingga (Finite Difference Methods). Dibandingkan dengan metode numerik lainya, metode

ini memiliki konsep yang lebih sederhana dan relatif mudah diaplikasikan dalam program komputer. Pendekatan persamaan differensial Poisson dan Laplace 2D menggunakan metode beda hingga, umumnya menghasilkan sistem persamaan linier yang 1

berukuran besar tetapi mempunyai struktur teratur dan elemennya kebanyakan nol. Penyelesaian selanjutnya sering dipilih metode iterasi. Hal - hal yang akan dikaji dalam tulisan ini antara lain -

Bagaimana menyelesaikan persamaan differensial Poisson dan Laplce 2D dengan syarat batas bertipe Dirichlet, Neumann, dan Robin menggunakan metode beda hingga dan metode iterasi SOR.

-

Hubungan antara error dan kekonvergenan, serta hal yang menjamin metode iterasi dapat konvergen.

-

Strategi

untuk

meningkatkan

akurasi

numerik

dan

mempercepat

konvergensi dengan memanfaatkan nilai-nilai pada syarat batas Dirichlet sebagai nilai tebakan pada awal iterasi. -

Aplikasi persamaan differensial Poisson dan Laplace 2D.

Untuk menyederhanakan masalah maka bentuk domain persamaan differensial Poisson dan Laplace yang ditinjau pada tulisan ini adalah segiempat.

1.2 Manfaat Penelitian

Hasil penelitian ini diharapkan dapat memberikan sumbangan pemikiran untuk penyelesaian masalah pada kasus-kasus yang berkaitan dengan persamaan differensial Poisson dan Laplace dalam bidang teknik, fisika, kelautan, dan bidang rekayasa lainnya.

2

1.3 Tujuan Penelitian

Penelitian ini bertujuan untuk mendapatkan strategi penyelesaian persamaan differensial Poisson dan Laplace 2D agar dihasilkan akurasi numerik yang baik dan proses komputasi yang cepat.

1.4 Tinjaun Pustaka

Basaruddin (1994) mengemukakan bentuk umum syarat batas persamaan differensial parsial orde dua dan cara mendapatkan rumus pendekatan beda hingga. Teori tentang barisan oleh Royden (1989), sedangkan Sturler (2003) mengemukakan teorema-teorema yang berhubungan dengan kekonvergenan metode iterasi Skema iterasi SOR dikemukakan oleh Bailey (2003), sedangkan O’brien (1986) memberikan teorema tentang nilai parameter yang digunakan metode SOR. Marek (2005) mengemukakan definisi tentang matriks konvergen. Luenberger (1969) memberikan aksioma tentang norm. Nakamura (2003) dan Constantinides (1987) memberikan petunjuk tentang pembuatan program komputer untuk menyelesaikan persamaan Poisson dan Laplace. Abayomi (2003) memberikan contoh persamaan differensial Poisson dan Laplace dengan solusi analitisnya telah diketahui, untuk dijadikan pembanding terhadap solusi numerik. Dalam penerapannya persamaan differensial Poisson dan Laplace sering dijumpai pada model matematika di bidang teknik dan fisika. Salah satu contohnya diberikan oleh Constantinides (1987) melalui persamaan panas steady state. 3

1.5 Metode Penelitian

Penelitian ini dilakukan dengan metode studi pustaka yakni dengan mengkaji karya ilmiah, jurnal, bulletin dan buku- buku yang berkaitan dengan materi penelitian. Selanjutnya dilakukan uji coba metode dan visualisasi melalui program komputer menggunakan software MATLAB versi 5.3 dan versi 7.

1.6 Sistematika Penulisan

Bab I adalah pendahuluan memuat latar belakang masalah, manfaat penelitian, tinjauan pustaka, metode penelitian, dan sistematika penulisan. Bab II berisi teori tentang persamaan differensial Poisson dan Laplace 2D dalam klasifikasi persamaan differensial parsial orde dua, bentuk umum syarat batas, teori tentang metode beda hingga, barisan bilangan real, dan norm. Bab III membicarakan teori tentang metode iterasi untuk menyelesaikan sistem persamaan linier, teorema kekonvergenan yang menjadi jaminan berhasil atau tidaknya metode iterasi bila diterapkan, penurunan skema iterasi SOR dan teorema mengenai nilai

parameter relaksasinya. Dikemukakan pula cara menyelesaikan

masalah persamaan differensial Poisson dan Laplace 2D dengan syarat batas bertipe Dirichlet , Neumann, dan Robin menggunakan metode beda hingga dan proses iterasi

SOR. Selain itu, dikemukakan pula cara

memanfatkan nilai pada batas bertipe

Dirichlet untuk tebakan awal iterasi.

Bab IV berisi penggunaan metode dan aplikasinya. Sebagian contoh pada bab ini diketahui solusi eksaknya, kemudian dijadikan pembanding terhadap solusi 4

numerik yang diperoleh. Bab ini juga membicarakan kasus khusus yakni kasus solusi tidak tunggal dan kasus metode iterasi gagal mendekati penyelesaian eksaknya. Bab ini diakhiri dengan contoh aplikasi yaitu mengaproksimasi bentuk lelehan logam. Bab V adalah penutup yang berisi kesimpulan dan saran untuk penelitian lebih lanjut yang belum dikaji dalam tulisan ini.

5

BAB II LANDASAN TEORI

2.1 Persamaan Differensial Poisson dan Laplace 2D (Nakamura, 1991)

Bentuk umum persamaan differensial parsial order dua adalah ∂ 2 u ( x, y ) ∂ 2 u ( x, y ) ∂ 2 u ( x, y ) ∂ u ( x, y ) ∂ 2 u ( x, y ) ∂ 2 u ( x, y ) A + B + C + D + E + F = f ( x, y ) ∂ x∂ y ∂ x 2 ∂ y 2 ∂ x ∂ x 2 ∂ x 2

dengan A,B,C,D,E,F, f adalah fungsi x dan y. Klasifikasi dibedakan dalam tiga tipe yaitu 2 1. Parabolik, jika B - 4AC = 0 2 2. Eliptik, jika B - 4AC < 0 2 3. Hiperbolik, jika B - 4AC > 0

Persamaan differensial Poisson dan Laplace merupakan bentuk khusus persamaan differensial bertipe eliptik. Bentuk persamaan differensial Laplace dan Poisson 2D adalah PD Laplace :

∂ 2 u ( x, y ) ∂ 2 u ( x, y ) + = 0 ∂ x 2 ∂ y 2

PD Poisson :

∂ 2 u ( x, y ) ∂ 2 u ( x, y ) = f(x,y) + ∂ x 2 ∂ y 2

jadi PD Laplace adalah bentuk homogen dari PD Poisson. Fungsi f ( x, y ) pada PD Poisson disebut fungsi sumber (source function). Makna fungsi sumber bergantung pada masalah yang dimodelkan oleh PD Poisson.

6

2.2 Bentuk Umum Syarat Batas (Basaruddin, 1994)

Bentuk umum syarat batas adalah Pu+Q

∂u =G ∂n

(2.2.1)

dengan n adalah vektor normal satuan (vektor tegak lurus pada kurva batas), P dan Q adalah bilangan konstan yang tidak bernilai nol secara serentak. Syarat batas Dirichlet diperoleh jika Q=0, sedangkan jika P=0 diperoleh syarat batas Neumann. Syarat batas (2.2.1) adalah syarat batas Robin atau syarat batas campuran (mixed condition).

2.3 Rumus Taylor

Diberikan fungsi dua peubah

u( x,y). Jika fungsi u( x,y) mempunyai

∂u ( x, y ) ∂ 2 u ( x, y ) ∂ n u ( x, y ) ∂ n +1u ( x, y ) turunan sampai ke (n+1), yaitu , ,..., , dan ∂ x ∂ x 2 ∂ x n ∂ x n +1 ∂u ( x, y ) ∂ 2 u ( x, y ) ∂ n u ( x, y ) ∂ n +1u ( x, y ) , ,..., , ∂ x ∂ x 2 ∂ x n ∂ x n +1

maka

rumus

Taylor

fungsi

u( x+ h0 , y+ k 0 ) disekitar ( x, y ) adalah u ( x + h0 , y + k 0 ) = u ( x, y ) + A1 + A2 + A3 + ...+ An+ Rn

(2.3.1)

dengan A1 = h0

∂u ( x, y ) ∂u ( x, y ) + k 0 ∂ x ∂ y

2 ⎛ 2 ∂ 2 u ( x , y ) ∂ u ( x , y ) ∂ u ( x , y ) 2 ∂ u ( x , y ) A2 = ⎜ h 0 + 2 h 0 k 0 + h0 2 ⎜ ∂ ∂ x y ∂ ∂ y 2 x ⎝

⎞ ⎟⎟ / 2 ! ⎠

2 2 3 ⎛ 3 ∂ 3u ( x, y) ∂ u ( x, y ) ∂u ( x, y ) 2 2 ∂u ( x, y ) ∂ u ( x, y ) 3 ∂ u ( x, y ) A3 = ⎜ h0 + 3h0 k 0 + 3h0 k 0 + k 0 3 2 2 ⎜ ∂ ∂ ∂ ∂ ∂ ∂ y 3 x x y x y ⎝

⎞ ⎟⎟ / 3 ! ⎠

7

⎛ ⎛ n ⎞ n ∂ nu( x, y) ⎛ n ⎞ n −1 ∂ n −1u( x, y) ∂u( x, y) ⎛ n ⎞ n ∂ nu( x, y) ⎞ ⎟/ n! An = ⎜ ⎜⎜ ⎟⎟h0 ⎜ ⎟ .... + ⎜ ⎟h0 k 0 + + ⎜⎜ ⎟⎟k 0 n n −1 n ⎜ 0 ⎟ 1 n ∂ x ∂ x ∂ y ∂ y ⎝ ⎠ ⎝ ⎠ ⎝ ⎝ ⎠ ⎠ ⎛ ⎛ n + 1 ⎞ n +1 ∂ n +1u (c1 , c 2 ) ⎛ n + 1 ⎞ n ∂ n u (c1 , c 2 ) ∂u (c1 , c 2 ) ⎟⎟h0 ⎟⎟h0 k 0 Rn = ⎜⎜ ⎜⎜ + ⎜⎜ + .... n +1 n 0 1 y ∂ x x ∂ ∂ ⎝ ⎝ ⎠ ⎝ ⎠ n ⎛ n + 1 ⎞ ⎛ n + 1 ⎞ n +1 ∂ n +1u (c1 , c 2 ) ⎞ n ∂u (c1 , c 2 ) ∂ u ( c1 , c 2 ) ⎟⎟h0 k 0 ⎟⎟k 0 ⎟⎟ /(n + 1) ! + ⎜⎜ + ⎜⎜ ∂ x ∂ y n ∂ y n +1 ⎝ n + 1 ⎠ ⎠ ⎝ n ⎠ dengan c1 suatu titik diantara x dan x+ h0 , c2 suatu titik diantara y dan y+ k 0 .

2.4 Deret Taylor Fungsi Dua Peubah untuk Mendapatkan Rumus Beda

Hingga Turunan Parsial Berdasarkan rumus Taylor

(2.3.1), maka deret Taylor u( x+ h0 , y+ k 0 )

disekitar ( x, y ) , adalah

u ( x + h0 , y + k 0 )

= u ( x, y ) + h0 u x ( x, y ) + k 0 u y ( x, y ) h0 u xx ( x, y ) + 2h0 k 0 u xy ( x, y ) + k 0 u yy ( x, y ) 2

+

2! h0 u xxx ( x, y ) + 3h0 k 0 u xxy ( x, y ) + 3h0 k 0 u xyy ( x, y ) + k 0 u yyy ( x, y ) 3

+

2

2

2

3

3!

+ ... .

Untuk h0 = h, k 0 = k , deret Taylor u( x+ h0 , y+ k 0 ) disekitar ( x, y ) menjadi u ( x + h, y ) = u ( x, y ) + h u x ( x, y ) +

h 2 u xx ( x, y )

2!

+

h 3 u xxx ( x, y )

3!

+ ... .

(2.4.1)

Untuk h0 = - h, k 0 = 0, deret Taylor u( x+ h0 , y+ k 0 ) disekitar ( x, y ) menjadi u ( x − h, y ) = u ( x, y ) − h u x ( x, y ) +

h 2 u xx ( x, y )

2!

−

h 3u xxx ( x, y )

3!

+ ... .

(2.4.2)

8

Untuk h0 = 0, k 0 = k , deret Taylor u( x+ h0 , y+ k 0 ) disekitar ( x, y ) menjadi u ( x, y + k ) = u ( x, y ) + k u y ( x, y ) +

k 2 u yy ( x, y )

2!

+

k 3u yyy ( x, y )

3!

+ ... . (2.4.3)

Untuk h0 = 0, k 0 = - k , deret Taylor u( x+ h0 , y+ k 0 ) disekitar ( x, y ) menjadi u ( x, y − k ) = u ( x, y ) − k u y ( x, y ) +

k 2 u yy ( x, y )

2!

−

k 3u yyy ( x, y )

3!

+ ... . (2.4.4)

Pada tulisan ini, x=xi=(i-1)h, y=y j=( j-1)k dengan i=1,2,3,…I, dan j=1,2,3,…J, maka untuk mempermudah penulisan skema beda hingga pada

bagian selanjutnya, didefinisikan u ( x, y ) = u ( xi , y j ) = u i , j u ( x + h, y ) = u ( xi + h, y j ) = u i +1, j u ( x − h, y) = u ( xi − h, y j ) = u i −1, j u ( x, y + k ) = u ( xi , y j + k ) = u i , j +1 u ( x, y − k ) = u ( xi , y j − k ) = u i , j −1 .

Persamaan (2.4.1) dapat ditulis menjadi 2

u(x+h,y) = u(x,y) + h u x(x,y) + O(h )

( O(h2) adalah error pemotongan yang merupakan fungsi kuadrat dari ukuran kisi pada arah sumbu X ), selanjutnya diperoleh u x(x,y) =

u ( x + h, y ) − u ( x, y ) h

+ O(h)

sehingga diperoleh rumus pendekatan beda maju (forward difference) untuk u x ( x, y ), yaitu

u x ( x, y ) ≈

u ( x + h, y ) − u ( x, y ) h

9

atau (u x ) i , j ≈

u i +1, j − u i , j

.

h

(2.4.5)

Persamaan (2.4.2), dapat ditulis menjadi u ( x − h, y ) = u ( x, y ) − h u x ( x, y ) + O(h 2 )

atau u x ( x, y ) =

u ( x, y ) − u ( x − h, y ) h

+ O(h)

sehingga diperoleh rumus pendekatan beda mundur ( backward difference) untuk

(u x ) i , j , yaitu (u x ) i , j ≈

u i , j − u i −1, j

(2.4.6)

.

h

Persamaan (2.4.1) dikurangi dengan persamaan (2.4.2) diperoleh 3

u( x+h,y ) - u( x-h,y) = 2 h u x( x,y) + O(h )

atau u x( x,y) =

sehingga

diperoleh

rumus

u ( x + h, y ) − u ( x − h, y )

2h

pendekatan

2

+ O(h )

beda tengah

(central

difference)

untuk (u x ) i , j , yaitu (u x ) i , j ≈

u i +1, j − u i −1, j

2h

.

(2.4.7)

10

Persamaan (2.4.1) dijumlahkan dengan persamaan (2.4.2) diperoleh 2

u( x+h,y) + u( x-h,y)

4

= 2 u(x,y) + h u xx(x,y) + O(h )

atau u ( x − h, y ) − 2u ( x, y ) + u ( x + h, y )

u xx( x,y) =

h

2

2

+ O(h )

sehingga diperoleh rumus beda tengah (central difference) untuk u xx ( x, y), yaitu u ( x − h, y ) − 2u ( x, y ) + u ( x + h, y )

u xx( x,y) ≈

h

2

atau u i −1, j − 2u i , j + u i +1, j

(u xx ) i , j ≈

h2

(2.4.8)

.

Rumus-rumus pendekatan beda hingga lainnya dapat diperoleh melalui deret Taylor (2.4.3) dan deret Taylor (2.4.4). Rumus-rumus beda hingga tersebut adalah u i , j +1 − u i , j

:

(u )

≈

backward difference :

(u )

≈

central difference

:

(u )

≈

central difference

:

(u )

≈

forward difference

y i , j

y i , j

y i , j

yy i , j

(2.4.9)

k u i , j − u i , j −1

(2.4.10)

k u i , j +1 − u i , j −1

(2.4.11)

2k u i , j −1 − 2u i , j + u i , j +1 k 2

.

(2.4.12)

11

2.5

Barisan Bilangan Real (Royden, 1989)

Definisi 2.5.1 v Barisan {u ( ) } ⊂ R dikatakan konvergen (mempunyai limit), jika terdapat bilangan

u sehingga ∀ε > 0 , terdapat bilangan asli v 0 , dan jika v > v 0 berlaku

| u (v) − u | <

ε

(2.5.1)

disingkat lim u ( v ) = u , dan dikatakan {u ( v ) } konvergen ke u atau {u ( v ) } berlimit u v →∞

untuk n → ∞ .

Teorema 2.5.1

Jika barisan {u ( v ) } ⊂ R naik monoton dan terbatas ke atas maka {u ( v ) } konvergen ke batas atas terkecilnya.

Bukti.

(a) Barisan {u ( v ) } naik monoton, jadi u

(v )

≤ u ( v +1) , v=0,1,2,... .

(2.5.2)

(b) Barisan {u ( v ) } terbatas ke atas, jadi {u ( v ) } mempunyai batas atas terkecil atau ada

bilangan

M0 = sup{u ( v ) }, ( sup{u ( v ) } selalu ada untuk setiap barisan

{u ( v ) } yang terbatas ke atas, berbeda halnya dengan maks{u ( v ) } ). Bilangan M0 mempunyai sifat i)

u

(v )

≤ M0, v = 0,1,2,... .

(2.5.3)

12

(v ) ii) ∀ε > 0 , terdapat bilangan asli v 0 ( ∃ u 0 ∈ {u ( v ) } ), sehingga berlaku

M0 -

ε

< u ( v ) ≤ M0 (< M0 + ε ) 0

(2.5.4)

Dari ketaksamaan (2.5.2), (2.5.3), dan (2.5.4) disimpulkan, ∀v ≥ v0 berlaku M0 - ε < u (

v)

≤ M0 (< M0 + ε )

M0 - ε < u ( ) < M0 + ε ) v

v | u ( ) − M 0 | < ε .

Dengan kata lain terbukti lim u ( v ) = M 0 = sup{u (v ) } . ▄ v →∞

Teorema 2.5.2

Jika

barisan

{u ( v ) } ⊂ R

turun monoton dan terbatas ke bawah maka

{u ( v ) } konvergen ke batas bawah terbesarnya.

Bukti.

(a) Barisan {u ( v ) } turun monoton, jadi u ( v ) ≥ u ( v +1) , v=0,1,2,... .

(2.5.5)

(b) Barisan {u ( v ) } terbatas ke bawah, jadi {u ( v ) } mempunyai batas bawah terbesar atau ada bilangan B0 = inf{u ( v ) } , ( inf{u ( v ) } selalu ada untuk setiap barisan {u ( v ) } yang terbatas ke bawah, berbeda halnya dengan min{u ( v ) } ). Bilangan B0 bersifat i)

u(

v)

≥ B0 , v=0,1,2,... .

(2.5.6)

13

(v ) ii) ∀ε > 0 , terdapat bilangan asli v 0 ( ∃ u 0 ∈ {u ( v ) } ) sehingga berlaku

B0 + ε > u

( v0 )

≥ B0 ( > B0 - ε ).

(2.5.7)

Dari (2.5.5), (2.5.6), dan (2.5.7) disimpulkan bahwa ∀v ≥ v 0 berlaku v B0 + ε > u ( ) ≥ B0 ( > B 0 − ε )

B 0 + ε > u ( v ) > B 0 − ε v B 0 − ε < u ( ) < B 0 + ε

v | u ( ) − B 0 | < ε .

Dengan kata lain terbukti lim u ( v ) = B0 = inf {u v } . ▄ v →∞

Teorema 2.5.3 v Setiap barisan {u ( ) } ⊂ R yang terbatas mempunyai barisan bagian {u

( vk )

} yang

konvergen

Bukti.

Diambil sebarang barisan {u ( v ) } ⊂ R yang terbatas. Dari barisan {u ( v ) } dibentuk u

( v k )

= maks{u (1) , u ( 2 ) ,..., u ( k ) }, ∀k

u

( v k )

≤ u (v

diperoleh

( Jadi {u

( {u

vk )

vk )

k

+1)

,

∀k .

} ⊂ {u ( v ) } barisan naik monoton dan terbatas . Menurut Teorema 2.5.1, v

( ) } konvergen ke sup{u k } . ▄

14

Definisi 2.5.2 v Barisan {u ( ) } ⊂ R disebut barisan Cauchy, jika untuk setiap

ε

> 0 terdapat

bilangan asli v 0 sehingga, jika m, v ≥ v0 berlaku | u (v ) − u (m) | <

ε

.

(2.5.8)

Teorema 2.5.4

Setiap barisan Cauchy {u ( v ) } ⊂ R merupakan barisan terbatas.

Bukti.

Diambil sebarang barisan Cauchy {u ( v ) } ⊂ R . Jadi ∀ε > 0 , ( khususnya

ε

= 1 ), ada bilangan asli v0 sehingga, jika m, v ≥ v0

v m berlaku | u ( ) − u ( ) | < 1 , jadi khususnya

( | u (v) − u

v0 )

| < 1

atau

− 1 < u ( v ) < u ( v ) + 1 , untuk v ≥ v0 .

(2.5.9)

( v 1) (v ) M0 = maks {u (1) , u (1) ,..., u 0 − , u 0 + 1}

(2.5.10)

( v 1) (v ) B0 = min {u (1) , u (1) ,..., u 0 − , u 0 − 1}

(2.5.11)

u

( v0 )

0

Diambil

Dari (2.5.9), (2.5.10) dan (2.5.11) diperoleh B0 ≤ u (v ) ≤ M0, ∀v . Jadi {u ( v ) } terbatas.

▄

15

Teorema 2.5.5 v v Jika barisan {u ( ) } ⊂ R konvergen maka {u ( ) } barisan Cauchy.

Bukti.

Diambil sebarang barisan {u ( v ) } konvergen, jadi ada bilangan u sehingga

∀ε > 0 , terdapat bilangan asli v 0 sehingga, jika v ≥ v0 berlaku

| u (v ) − u | <

ε

/2.

Jadi, jika m, v ≥ v 0 berlaku | u (v ) − u ( m) | = | u (v) − u + u − u (v ) |

≤ | u ( v) − u | + | u − u (v ) | <

ε

/2 +

ε

/2 =

ε

.

Jadi {u ( v ) } barisan Cauchy. ▄

Teorema 2.5.6

Setiap barisan {u ( v ) } Cauchy di R merupakan barisan konvergen.

Bukti.

Diambil sebarang barisan Cauchy {u ( v ) } ⊂ R. Menurut Teorema 2.5.4, barisan {u ( v ) } terbatas. Oleh karena itu,

(v ) ada {u k } ⊂ {u ( v ) } yang konvergen, katakan

(v ) {u k } konvergen ke l.

16

Diambil bilangan

ε

> 0 sebarang.

(i) Karena {u (v ) } ⊂ R barisan Cauchy maka ada bilangan asli v1 sehingga , jika m, v ≥ v1 berlaku | u ( v ) − u ( m) | <

ε

/2 .

(2.5.12)

(v ) (ii) Karena {u k } konvergen ke l, maka ada bilangan asli v2 sehingga,

jika vk ≥ v2 berlaku |u

( vk )

−l | <

ε

/2.

(2.5.13)

Diambil bilangan asli v 0 = maks { v1, v2 } Jadi dari (2.5.12) dan (2.5.13) diperoleh bahwa, jika v ≥ v 0 berlaku (v ) (v ) | u (v) − l | = | u (v ) − u 0 + u 0 − l |

≤ | u (v ) − u (v ) | + | u (v ) − l | 0

<

ε

/2 +

ε

/2 =

0

ε

.

Jadi {u ( v ) } konvergen ke l, ( lim u (v ) = l ). ▄ v →∞

Teorema 2.5.7

Jika {u ( v ) } ⊂ R konvergen maka limitnya tunggal

Bukti.

Diketahui {u ( v ) } ⊂ R konvergen.

17

Andaikan ada bilangan l dan u sehingga lim u ( v ) = l v →∞

dan lim u ( v ) = u . v →∞

Jadi ∀ε > 0 , ( ε / 2 > 0 ) terdapat bilangan asli v1 dan v2 sehingga | u (v ) − l | <

ε

/ 2 , ∀v > v1

(2.5.14)

| u (v ) − u | <

ε

/ 2 , ∀v > v 2

(2.5.15)

dan

diambil bilangan asli v 0 = maks { v1, v2 }, berlaku 0 ≤ | l – u | = | l – u + u - u | (v)

(v)

≤ | l – u(v)| + | u(v)- u | < ε / 2 +

ε

/2 =

ε

.

Jadi ∀ε > 0 , berlaku 0 ≤ | l – u | < ε . Dengan kata lain l = u.

▄

Definisi 2.5.3

Ruang metrik ( X ,d) dikatakan lengkap ( complete) jika setiap barisan Cauchy di dalamnya konvergen. Berdasarkan Teorema 2.5.6 dan Definisi 2.5.3, disimpulkan bahwa ruang metrik ( R,d), dengan d( x, y)= | x – y | ∀ x, y ∈ R , adalah ruang metrik lengkap. Mengingat ( R,d) adalah ruang metrik lengkap, maka ketaksamaan (2.5.8) pada

bab

selanjutnya

digunakan

sebagai

kriteria

kekonvergenan

pada

ketaksamaan (3.5.18).

18

2. 6 Norm Definisi 2.6.1 (Luenberger,1969)

Ruang vektor X bernorm linier adalah ruang vektor yang padanya didefinisikan fungsi bernilai real yang memetakan setiap elemen x di X ke dalam bilangan real || x|| yang disebut norm dari vektor x. Norm memenuhi aksioma berikut 1) || x|| ≥ 0, ∀ x ∈ X , || x||=0, ⇔ x = 0 2) || x + y || ≤ || x|| + || y||, ∀ x, y ∈ X 3) || α x|| = | α | || x||, untuk setiap skalar α , ∀ x ∈ X

Teorema 2.6.1

Dalam ruang norm linier X , berlaku || x|| - || y|| ≤ || x – y ||, ∀ x, y ∈ X . Bukti.

|| x|| - || y|| = || x – y + y || - || y|| ≤ || x – y || + || y|| - || y|| = || x – y ||. ▄

2.6.1

Norm Vekor

Diberikan vektor n-dimensi

⎡ x1 ⎤ ⎢ x ⎥ ⎢ 2⎥ x = ⎢ x3 ⎥ . ⎢ ⎥ ⎢M ⎥ ⎢ x n ⎥ ⎣ ⎦

19

Secara umum norm vektor x ditulis || x|| atau sering pula ditulis dengan |x|, adalah bilangan nonnegatif yang didefinisikan sedemikian sehingga memenuhi aksioma pada Definisi 2.6.1. Norm || x || p , dengan p= 1,2,... didefinisikan sebagai

⎛ p ⎞ || x || p = ⎜ ∑ xi ⎟ ⎝ i ⎠

1 / p

.

Kasus khusus || x || ∞ didefinisikan sebagai || x || ∞ = maks xi . i

Norm vektor yang paling umum ( biasanya hanya ditulis norm saja) adalah l2-norm, yaitu

|| x || 2 =|| x ||= x12 + x 22 + ... + x n2

2.6.2

Matriks Transpose Conjugate

Diberikan matriks A

dengan komponen bilangan real atau kompleks.

Matrik A H adalah transpose conjugate matrik A. Nama lain dari matriks transpose conjugate adalah matriks adjoint , Hermitian conjugate, atau tranjugate.

Simbol yang biasanya digunakan A* , A H , dan A† ( pada teori kuantum). Contoh

⎡1 A = ⎢ ⎣2 + 3i

4 - 6i ⎤

⎥,

7 ⎦

⎡1 A H = ⎢ ⎣4 + 6i

2 − 3i ⎤

⎥.

7 ⎦

20

Definisi 2.6.2.1

Matrik E disebut sebagai unitary matrix, jika EE H = I , dengan I matrik identitas.

2.6.3

Norm Matriks

Diberikan matriks A bujursangkar dengan komponen bilangan real atau kompleks. Norm matriks A yaitu || A|| adalah bilangan nonnegatif yang diasosiasikan dengan matrik A sedemikian sehingga memenuhi aksioma pada Definisi 2.6.1. Adapun contoh-contoh norm matriks diberikan sebagai berikut. Definisi p-norm A p , dengan 1 ≤ p ≤ ∞ , yaitu A p = maks Ax

p

x dengan || x|| p =1

dengan || x || p menyatakan norm vektor. Norm A 1 , didefinisikan sebagai A 1 = maks j

Norm spectral, yaitu A

2

n

∑| a

i , j

|.

i =1

didefinisikan sebagai H

A 2 = ( maksimum nilai eigen A A )1 / 2 .

Norm Hilbert- Schmidt dari suatu matrik A didefinisikan sebagai A HS =

∑a

2 i , j

.

i , j

Norm Frobenius dari suatu matrik A didefinisikan sebagai

A

= F

m

n

∑∑ a i

2 i , j

= Trace( AA H ) .

j

21

BAB III MENYELESAIKAN PD POISSON DAN LAPLACE 2D DENGAN METODE BEDA HINGGA DAN METODE ITERASI SOR

3.1 Hubungan PD Poisson dan Laplace 2D, Metode Beda Hingga

dan

Barisan Bilangan Real. Barisan bilangan real muncul pada proses iterasi skema beda hingga PD poisson dan Laplce 2D. Secara numeris, untuk menghentikan proses iterasi dibutuhkan kriteria konvergensi yang merujuk pada barisan Cauchy. Iterasi dapat tidak berhenti bila barisan yang dihasilkan bukan barisan Cauchy menurut kriteria yang diberikan, dalam hal ini dinyatakan iterasi tidak konvergen. Proses iterasi dimulai dengan mengambil sebarang nilai awal iterasi u ( 0 ) . Karena pemilihan nilai awal iterasi u ( 0 ) sebarang, maka u ( 0 ) disebut pula tebakan awal iterasi.

Nilai awal iterasi u ( 0 ) sebarang

Skema iterasi menghasilkan barisan (1) (2) (3) (v-1) (v) u , u , u ,..., u , u ,...

Secara numeris : Secara teoritis : (v) Jaminan u → u ?

Jika proses iterasi memenuhi kriteria yang merujuk pada barisan Cauchy maka dianggap (v) u →u (a)

22

R

(2)

ui,j

ui,j

(4)

ui,j

(vo -1)

ui,j ui,j

(vo)

 (vo) (vo -1 ) |< ε  | ui,j – ui,j 

(3)

(1)

ui,j ui,j

(0)

0

1

2

3

4 ...

v0 -1

v0

V (Nomor Iterasi)

Keterangan : ε : Nilai toleransi kekonvergenan u : Solusi sebenarnya (v) u : Solusi iterasi ke-v yang digunakan untuk mendekati u v0 : Banyaknya iterasi untuk mencapai kriteria kekonvergenan (nomor iterasi terakhir )

u (v)

u1(,v1)  u1,1   ( v)    u1, 2  u  =   , u =  1, 2  , v = 1,2,3,... . Μ  Μ  u I(,vJ)  u I ,J     

(b) Gambar 3.1.

(a) Ilustrasi kekonvergenan secara teoritis dan secara numeris (b) Ilustrasi barisan Cauchy untuk menghentikan iterasi

Pembahasan kekonvergenan metode iterasi pada bagian berikut ini dimaksudkan untuk mencari syarat yang menjadi jaminan suatu proses iterasi konvergen.

23

3.2 Kekonvergenan Metode Iterasi Definisi 3.2.1, (Bailey, 2003) Misalkan nilai-nilai eigen dari matriks B adalah λ i , i=1,2,...m. Spectral radius matriks B adalah

µ B

maks | λ i | .

=

(3.2.1)

i =1, 2,...m

Definisi 3.2.2, (Marek, 2005) Diberikan B matriks bujursangkar. Matrik B disebut matriks konvergen, jika lim B v →∞

v

ada. Matrik B disebut matriks konvergen nol, jika lim B v = 0. v →∞

Hubungan kekonvergenan metode iterasi dengan nilai spectral radius suatu matriks yang disebut matriks iterasi akan dibicarakan berikut ini. Diberikan sistem persamaan linier Au

=

(3.2.2)

b

dengan A matriks berukuran mxm, b adalah matriks berukuran mx1, dan u adalah matriks berukuran mx1 (u adalah matriks yang akan dicari). Untuk mendapatkan skema iterasi, maka persamaan (3.2.2) ditulis menjadi Nu

=

(3.2.3)

( N − A)u + b

dengan N adalah matriks yang mempunyai invers. Matriks N menentukan jenis metode iterasi (metode Jacobi, Gauss-Seidel, SOR, dan lain-lain). Jika N –1 dikalikan dengan persamaan (3.2.3) diperoleh u

=

N 1 ( N − A)u + N 1b −

−

24

atau u = Mu + c ,

v = 1,2,3,Λ

(3.2.4)

dengan M=N –1(N-A) , c = N −1b . Dibentuk skema iterasi dari persamaan (3.2.4), yaitu u(

v)

= Mu

( v −1)

+c

,

v = 1,2,3,Λ

(3.2.5)

–1 dengan u ( 0 ) diberikan sebagai tebakan awal, M=N (N-A), c = N −1b .

Matriks M disebut matriks iterasi (iterative matrix). Skema iterasi (3.2.5) (0) (1) (2) ( 3) selanjutnya menghasilkan barisan vektor u , u , u , u ,Λ .

( ) Diketahui vektor u v adalah vektor pada iterasi v yang digunakan untuk

mendekati solusi sebenarnya u, jadi vektor errornya (error global) adalah e(

v)

=

u −u(

v)

(3.2.6)

.

Bila persamaan (3.2.4) dikurangi dengan persamaan (3.2.5) diperoleh u −u

(v)

= M (u − u

( v −1)

(3.2.7)

)

atau e(

sehingga, e (1)

= Me

(0)

v)

= Me

, e ( 2)

= Me

(v )

= M

( v −1)

(1)

;

v=1,2,3 ,

= M

(3.2.8)

...

2

e

( 0)

, ....

jadi e

v

e

(0)

(3.2.9)

.

Teorema 3.2.1, (Sturler, 2003)

lim M v x = 0, untuk sebarang x v →∞

⇔

lim M v = 0 v →∞

25

Bukti.

lim M v x = 0, untuk sebarang x ⇒

(a)

v →∞

lim M v = 0 v→∞

 x1   x  2 v Diketahui lim M x = 0, untuk sebarang x =   v →∞ Μ     x m  lim || M v xi || = 0 , untuk sebarang xi , i =1,2,...,m.

jadi

v→∞

Diambil ( E unitary matrix ) E = [e1 e2 Λ em ] . Karena lim || M v xi || = 0 , berlaku untuk sebarang xi , i =1,2,...,m, jadi berlaku v →∞

pula untuk xi = ei , i =1,2,...,m, yaitu berakibat

lim || M v ei || = 0, v →∞

i = 1,2, Λ , m,

lim || M v E || = 0 . v →∞

Jadi || M v || = || M V E E H || ≤ || M v E || || E H || = || M v E || . lim || M v || ≤ lim || M V E || = 0 . v →∞

Karena 0 ≤ || . || , maka

v →∞

lim || M v || = 0 , berakibat v →∞

lim M v = 0. v →∞

(b) lim M v = 0 ⇒ lim M v x = 0, untuk sebarang x v →∞

v →∞

v Diketahui lim M = 0 ⇒ lim || M v || = 0 . v →∞ v→∞

Karena

|| M v x || ≤ || M v || || x || , untuk sebarang x,

maka

lim || M v x || ≤ lim || M v || || x || = 0 . v →∞

v →∞

26

Karena 0

≤

lim || M v x ||

|| . || , maka

=

v →∞

0

berakibat lim M v x = 0, untuk sebarang x . • v →∞

Karena Teorema 3.2.1 berlaku untuk sebarang vektor x, jadi berlaku pula untuk x= e (0 ) , sehingga diperoleh

lim M v e ( 0 ) = 0

⇔

v →∞

lim M v = 0. v →∞

Teorema 3.2.2, (Sturler, 2003)

lim M v = 0

⇔

v →∞

µ M

<1

dengan µ M adalah spectral radius matrik M .

Bukti. v (a) lim M = 0 v →∞

⇒

µ M

< 1.

Diambil sebarang λ , nilai eigen matriks M . Dari persamaan nilai eigen, yaitu Mu

=

λ u , didapat M v u

Karena lim M v v →∞

= 0,

=

λ v u .

maka 0 = lim || M v u ||

=

lim | λ | v || u ||

0

v →∞

v →∞

Karena | λ |< 1 , jadi µ M

=

lim | λ |v || u || v →∞

⇔

| λ |< 1 .

< 1.

27

(b)

lim M v = 0 .

µ M < 1 ⇒

v →∞

Untuk setiap ε > 0 dan matrik M ,

terdapat norm (sebut

|| . || a ), sedemikian

sehingga || M || a ≤ µ M + ε . Diketahui µ M

< 1,

dibentuk ε = (1 − µ M ) / 2 , jadi µ M + ε < 1 . v

|| M || a

≤

v

|| M || a ≤ ( µ M + ε )

lim || M v || a

≤

lim || M v || a

=

v →∞

lim ( µ M + ε ) v v →∞

v →∞

v

=

0

0.

Karena setiap norm (dalam ruang vektor berdimensi hingga) saling ekuivalen, v jadi lim || M ||

=

v →∞

0 , untuk sebarang norm || . ||, berakibat lim M v = 0. • v →∞

Dari persamaan (3.2.9), Teorema 3.2.1 dan Teorema 3.2.2, diperoleh lim e ( v )

=

v→∞

Karena e ( v )

=

u − u(

v)

,

lim M v e ( 0)

=

v →∞

0

⇔ µ M < 1

selanjutnya diperoleh lim u ( v→∞

v)

=

u

⇔ µ M < 1 .

Jadi metode iterasi konvergen jika dan hanya jika nilai µ M

<1

(nilai spectral

radius matriks iterasinya kurang dari satu).

28

3.3 Skema Iterasi Metode SOR (Bailey, 2003) Diberikan sistem persamaan linier (3.3.1)

Au=b

dengan A matriks berukuran mxm, b adalah matriks berukuran mx1 , dan u adalah matriks yang dicari berukuran mx1. Matriks A diasumsikan dapat dinyatakan sebagai (3.3.2)

A = D - L –U

dengan D matriks diagonal, L

matriks segitiga bawah dengan komponen

diagonalnya nol ( strictly lower triangular ),

U matriks

segitiga atas dengan

komponen diagonalnya nol ( strictly upper triangular ). Persamaan (3.3.2) disubstitusi pada persamaan (3.3.1), diperoleh ( D - L – U )u = b.

(3.3.3)

Persamaan (3.3.3) dikalikan dengan bilangan

ω

(bilangan

ω

selanjutnya disebut

parameter SOR), diperoleh ω

( D − L − U )u

= ω

(3.3.4)

b

kedua ruas ditambahkan Du, menjadi ω

( D − L − U )u + Du

= ω

b + Du

selanjutnya ditulis menjadi ( D − ω L)u = (ω U + (1 − ω ) D )u

+ ω

b

sehingga diperoleh skema iterasi SOR, yaitu

29

( D − ω L)u ( v)

= ω U +

(

(1 − ω ) D )u ( v

1)

−

+ ω b

,

v

= 1,2,3,...

(3.3.5)

atau u

( v)

=

( D − ω L ) 1 (ω U + (1 − ω ) D ) u ( v −

1)

−

+

( D − ω L)

1

−

ω b

atau u

dengan

u

(0)

(v)

= M ω u

+

c, v

= 1,2,3,...

(3.3.6)

diberikan sebarang, M ω

Matriks

( v −1)

M ω

=

( D − ω L) 1 (ω U + (1 − ω ) D ) , −

c

=

( D − ω L)

1

−

ω b

.

disebut matriks iterasi metode SOR. Adapun mengenai nilai

parameter SOR ( ω ) akan dibicarakan pada bagian berikut ini.

3.4 Parameter Relaksasi Metode SOR (O’ Brien, 1986) Teorema 3.4.1

Diberikan sistem persamaan linier Au

Jika M adalah matriks

=

b.

iterasi

ω

dan dapat dinyatakan sebagai

metode SOR dari matriks A yang berukuran A=I–L–U ,

dengan I matriks identitas,

segitiga bawah dengan komponen diagonal nol, serta matriks

mxm

L

matriks

U adalah

matriks

segitiga atas dengan komponen diagonal nol, dan M ω

=

( I − ω L) 1 (ω U + (1 − ω ) I ) −

(3.4.1)

30

maka (3.4.2)

µ ω ≥ ω − 1 dengan µ ω adalah spectral radius dari matriks M ω .

Bukti.

Dibentuk ϕ (λ ) = det(λ I − M ω ) . Karena M ω

=

(3.4.3)

− ( I − ω L) 1 (ω U + (1 − ω ) I ) , maka

ϕ (λ ) = det(λ I − ( I − ω L) −1 (ω U + (1 − ω ) I )) =

det( I − ω L) −1 (( I − ω L)λ I − (ω U + (1 − ω ) I )))

=

det( I − ω L) −1 ((λ I − ω Lλ ) − (ω U + (1 − ω ) I )))

=

1 det( I − ω L) − ((λ ( I − ω L) − (ω U + (1 − ω ) I )))

=

det( I − ω L) −1 . det(λ ( I − ω L) − (ω U + (1 − ω ) I ))) .

1 Karena det( I − ω L) = 1, det( I − ω L) −

=

1 det ( I − ω L)

= 1,

maka

ϕ (λ ) = det(λ ( I − ω L) − (ω U + (1 − ω ) I ) . Dilain pihak, misalkan λ i (ω ) adalah pembuat nol ϕ ( λ )

(3.4.4) dengan

i=1,2,3 ,...,m, maka ϕ (λ ) dapat dinyatakan sebagai

ϕ (λ ) = (λ − λ 1 (ω ) )(λ − λ 2 (ω ))Λ (λ − λ m (ω )) .

(3.4.5)

Jika diberikan λ = 0 diperoleh

ϕ (0) = (− λ 1 (ω ) )(− λ 2 (ω ) )Λ (− λ m (ω ))

(3.4.6)

31

atau m

ϕ (0) = (−1)

m

∏ λ (ω ) i

(3.4.7)

.

i =1

Jika ë = 0 diberikan pada persamaan ( 3.4.4 ), diperoleh

ϕ (0) = det((ω − 1) I − ω U ) .

(3.4.8)

Diketahui matriks I dan U masing-masing berukuran m x m dengan bentuk

1   I =    

1 1

Ο

   ,   1 

0   U =     

η 1, m 

η 1,2 Λ 0 η 2,3



η 2, m  η 3, m  .

0

Ο 0

   

Jadi

(ω − 1) I − ω U

(ω − 1)  0    (ω − 1)    − (ω − 1) =     Ο     (ω − 1)  

ωη1,2 0

Λ

ωη1 ,m 



ωη 2,3 Λ ωη 2 ,m  0 ωη 3,4 Λ ωη 3 ,m  

Ο 0

  

- ωη1 ,m  (ω − 1) - ωη1,2 Λ   (ω − 1) - ωη 2,3 Λ - ωη 2 ,m   (ω − 1) - ωη3,4 - ωη3 ,m  =      Ο   (ω − 1)  

diperoleh det((ω − 1) I − ω U ) = (ω − 1) m .

32

Dari persamaan (3.4.8) didapat

ϕ (0) = det((ω − 1) I − ω U ) = (ω − 1) m .

(3.4.9)

Dari persamaan (3.4.9 ) dan persamaan (3.4.7) diperoleh m

m

ϕ (0) = (ω − 1) = (−1)

m

∏ λ i

.

(3.4.10)

i =1

Diketahui µ ω adalah spectral radius matriks M ω , maka m

ϕ (0) = (ω − 1) m = (−1) m ∏ λ i ≤ ( µ ω ) m i =1

sehingga

(ω − 1)m

m

≤ (µ ω )

jadi

ω − 1 ≤ µ ω . •

(3.4.11)

Karena metode iterasi (termasuk metode SOR) konvergen jika dan hanya jika nilai spectral radius matriks iterasinya kurang dari 1, maka pemilihan parameter SOR harus memenuhi,

ω − 1 ≤ µ ω < 1 atau

ω − 1 < 1 jadi 0 < ω < 2 .

33

3.5

Menyelesaikan Persamaan Differensial Poisson dan Laplace 2D Menggunakan Metode Beda Hingga dan Metode Iterasi SOR Diberikan persamaan differensial Poisson 2D

2

∂ u ( x, y ) ∂ x

2

2

+

∂ u ( x, y) ∂ y

2

= f ( x, y )

(3.5.1)

dengan domain Ω = {( x, y ) | 0 ≤ x ≤ p, 0 ≤ y ≤ q}. Batas-batas domain dapat bertipe Dirichlet, Neumann, atau Robin. Penyelesaian persamaan (3.5.1) dengan metode beda hingga, dimulai dengan mempartisi domain seperti pada Gambar 3.2.

Gambar 3.2

Domain dan partisinya

34

Persamaan (3.5.1) selanjutnya ditulis menjadi

(

∂ 2u ∂ 2u ) ( ) = + ∂ x 2 i , j ∂ y 2 i , j

(3.5.2)

f i , j

jika digunakan rumus pendekatan

(

u i −1, j − 2u i , j + u i +1, j ∂ 2u ≈ ) , i j h2 ∂ x 2

(

u i , j −1 − ∂ 2u ) ≈ 2 i , j ∂ y

dan 2u i , j

+

u i , j +1

2

k

maka persamaan (3.5.2) dapat didekati dengan skema u i −1, j

− 2u i , j + u i +1, j h

2

+

u i , j −1

− 2u i , j + u i , j +1

= f i , j

2

k

atau u i −1, j u i , j

=

+ u i +1, j h2

+

u i . j −1

+ u i , j +1

k 2 1 

 1 2 2 + 2   h k 

− f i, j .

(3.5.3)

Persamaan (3.5.3) merupakan skema untuk mencari u pada titik-titik grid yang terletak pada bagian dalam domain, jadi berindeks i=2,3,...,I-1 j=2,3,...,J-1. Adapun untuk titik-titik grid yang terletak pada batas domain dalam hal ini u i , j dengan salah satu atau kedua indeksnya adalah i=1, i=I, j=1, j=1, j=J, dibutuhkan modifikasi persamaan (3.5.3) sesuai syarat batas yang diberikan ( Robin atau Neumann). Pada batas bertipe Dirichlet , nilai u

telah diketahui, sehingga tidak

35

dibutuhkan skema numerik untuk mencari nilai u pada batas tersebut. Pada sudut batas yang dibentuk oleh dua batas bertipe Dirichlet , nilai u diasumsikan sama dengan nilai u rata-ratanya. Jika sudut batas dibentuk oleh batas bertipe Dirichlet dan batas lainnya bertipe Robin atau Neumann maka nilai u pada sudut batas tersebut diasumsikan mengikuti nilai u dari batas Dirichlet . Nilai-nilai u pada batas bertipe Robin atau Neumann belum diketahui, oleh karena itu dibutuhkan skema numerik untuk mencari nilai u pada batas-batas tersebut. Diketahui bentuk umum syarat batas merupakan tipe Robin, yaitu P u(x,y) + Q

∂u ( x, y ) ∂n

= G

dengan n adalah vektor arah normal satuan (vektor satuan yang tegak lurus pada kurva batas).

n

Y

n ∂u ∂n

=−

∂u ∂n

=

∂u ∂ y

n

Domain

∂u

∂u

∂ x

∂n

n

∂u ∂n

Gambar 3.3

=−

∂u

=

∂u ∂ x

X

∂ y

Vektor normal satuan

36

Jika suatu batas tidak bertipe Dirichlet (Q• 0) maka bentuk syarat batas Robin dapat dinyatakan sebagai

∂u ( x, y ) ∂n

= α u ( x, y ) + β

(3.5.4)

dengan α = − P / Q , β = G / Q . Jika α = 0 , maka persamaan (3.5.4) menjadi syarat batas bertipe Neumann. Skema numerik untuk mencari nilai u pada suatu batas bertipe Robin atau Neumann pada prinsipnya adalah modifikasi persamaan (3.5.3) sesuai syarat batas

yang diberikan. Hal ini akan dicontohkan pada bagian berikut.

∂u

Y

∂ y

=

α 2 u

+

β 2

q

∂u ∂ x

∂u

= α 3u + β 3

∂ x

0

X ∂u ∂ y

Gambar 3.4

= α 4 u + β 4

=

α 1 u

+

β 1

p

Domain dengan syarat batas Robin atau Neumann

37

Misalkan domain pada Gambar 3.4 dengan batas batas kiri bertipe Robin, yaitu

∂u( x, y ) = α 3 u( x, y) + β 3 , ∂ x

x

= 0,

0 ≤ y ≤ q

(3.5.5)

dengan α 3 , β 3 bilangan konstan. Persamaan (3.5.5) dapat ditulis menjadi (

∂u ) = α 3u i , j + β 3 , ∂ x i, j

i

= 1,

= 1 ,2,..., J .

(3.5.6)

i=1; j=1,2,...,J.

(3.5.7)

j

Jika digunakan rumus pendekatan (

u i +1, j − u i −1, j ∂u ) i , j ≈ 2h ∂ x

maka persamaan (3.5.6) dapat didekati dengan skema u i +1, j

− ui −1, j 2h

= α 3ui , j + β 3 ;

Jika i=1 digunakan pada skema PD Poisson (3.5.3) dan persamaan (3.5.7) maka akan dijumpai u 0, j dengan j=1,2,...,J. Sedangkan diketahui bahwa indeks terkecil untuk i adalah 1, ini berarti u 0, j adalah titik-titik fiktif . Oleh karena itu u0,j tidak dapat digunakan secara langsung. Persamaan (3.5.3) untuk i=1, diperoleh

u 0, j u1, j

=

+ u 2, j h

2

+

u1. j −1

+ u1, j +1 2

k 1

 1  2 2 + 2   h k 

− f 1, j (3.5.8)

38

Persamaan (3.5.7) untuk i=1, diperoleh u 2, j − u 0 , j

2h

= α 3u1, j + β 3

atau u 0, j = u 2 , j − 2h β 3 − 2hα 3u1, j .

(3.5.9)

Persamaan (3.5.9) disubstitusi pada persamaan persamaan (3.5.8), diperoleh 2u 2, j − 2h β 3 u1, j =

h

2

+

u1. j −1 + u1, j +1 2

k 1

− f 1, j .

1 ⎞ ⎛ 1 2⎜ 2 + 2 ⎟ + 2 2hα 3 k ⎠ h ⎝ h

Persamaan (3.5.10) merupakan skema untuk mencari

u1, j ,

(3.5.10)

j=2,3,..., J-1.

Sedangkan skema untuk mencari u1,1 dan u1, J belum dapat ditentukan, karena melibatkan syarat batas lain yang membentuk sudut-sudut tersebut. Misalkan domain pada Gambar 3.4, batas bagian atas domain domain juga bertipe Robin,yaitu

∂u ( x, y ) = α 2 u ( x, y ) + β 2 , 0 ≤ x ≤ p, y = q ∂ y

(3.5.11)

dengan α 2 , β 2 bilangan konstan. Persamaan (3.5.11) dapat ditulis menjadi (

∂u ) i , j = α 2 u i , j + β 2 , i=1,2,...,I, j=J. ∂ y

(3.5.12)

Jika digunakan rumus pendekatan (

u i , j +1 − u i , j −1 ∂u ) i , j ≈ ∂ y 2k

39

maka persamaan (3.5.12) dapat didekati dengan u i , j +1 − u i , j −1

2k

= α 2 u i , j + β 2 ; i=1,2,...,I; j=J

atau u i , J +1 = u i ,J −1 + 2k β 2 + α 2 ,

i=1,2,...,I .

(3.5.13)

Diketahui indeks terbesar untuk j adalah J. Jadi ui ,J +1 , i=1,2,3,...,I adalah titik-titik fiktif. Skema PD Poisson (3.5.3) untuk j = J adalah u i −1,J + u i +1,J u i,J =

h

2

+

u i .J −1 + u i ,J +1 2

k 1 ⎞

− f i ,J (3.5.14)

⎛ 1 2⎜ 2 + 2 ⎟ k ⎠ ⎝ h

Persamaan (3.5.13) disubstitusi ke persamaan (3.5.14) diperoleh skema untuk u i ,J yaitu u i −1, J − u i +1, J u i ,J =

h

2

+

2u i .J −1 + 2k β 2 k 2 1

1 ⎞ ⎛ 1 2⎜ 2 + 2 ⎟ − 2 2k α 2 k ⎠ k ⎝ h

− f i , J

(3.5.15)

dengan i=2,3,...,I-1. Untuk u1,J dan u I,J ( u pada sudut batas domain) masih memerlukan informasi tambahan dari syarat batas yang lain yang membentuk sudut-sudut tersebut. Karena diketahui domain pada Gambar 3.4, batas kiri dan batas kanan domain bertipe Robin, maka skema u pada sudut kiri atas ( u1,J ) diperoleh dengan mensubstitusi persamaan (3.5.9) dan (3.5.13) ke persamaan (3.5.3), diperoleh

40

2u 2 ,J − 2h β 3 h

u1,J =

2

+

2u i .J −1 + 2k β 2 2

k

− f 1,J

1 ⎞ 1 1 ⎛ 1 2⎜ 2 + 2 ⎟ + 2 2hα 3 − 2 2k α 2 k ⎠ h k ⎝ h

.

(3.5.16)

Khusus untuk sudut batas yang dibentuk oleh batas-batas Dirichlet , nilai u pada sudut diasumsikan sama dengan nilai u rata-ratanya. Jika sudut batas dibentuk oleh batas bertipe Dirichlet dan batas yang lainnya bertipe Robin atau Neumann maka nilai u pada sudut tersebut diasumsikan mengikuti nilai u dari

batas bertipe Dirichlet . Setelah semua skema untuk u tersedia, maka untuk u yang belum diketahui nilainya diberikan sebarang nilai awal u ( 0 ) , selanjutnya diiterasi menggunakan skema SOR, yaitu (u i , j ) v = ω .u i , j + (1 − ω ).(u i , j ) v −1

)

(3.5.17)

dengan v nomor iterasi ( v = 1,2,3,...), u i , j adalah ui,j dari skema titik grid pada

)

bagian dalam domain (persamaan (3.5.3)) dan skema batas Robin atau Neumann. Parameter SOR ω dipilih 0 < ω < 2 . Proses iterasi terus dilakukan sampai memenuhi kriteria konvergensi yang diberikan, yaitu

(u ) − (u ) −

v 1

v

i , j

i , j

< ε , ∀i, j

(3.5.18)

dengan nilai toleransi kekonvergenan, biasanya dipilih nilai 0 < ε < 1 . Setelah ketaksamaan (3.5.18) dipenuhi, maka iterasi dihentikan. Akhirnya bilangan asli v0 yang dimaksud oleh barisan Cauchy (Definisi 2.5.2), identik dengan nomor iterasi terakhir ( v =1,2,3,..., v0 ).

41

3.6 Mengatur Nilai Awal Iterasi 3.6.1

Mengatur Nilai Awal Iterasi pada Batas Domain

Nilai-nilai u pada batas domain bertipe Dirichlet dapat digunakan untuk ( 0) mendapatkan nilai awal u pada batas lainnya yang bertipe Robin atau

Neumman. Sebagai contoh, misalkan syarat batas bagian kiri domain pada

Gambar 3.5 bertipe Dirichlet , yaitu u1, j

=

g j , j

= 1,2, Λ

,J

(3.6.1.1)

dan batas bagian kanan domain bertipe Robin, yaitu (

∂u ∂ x

) I , j = α 4 u + β 4 , j=1,2,...,J

(3.6.1.2)

dengan α 4 , β 4 bilangan konstan. Jika α 4 = 0 , batas bagian kanan domain bertipe Neumann.

Gambar 3.5.

Mengatur tebakan awal

u

pada batas bertipe Robin

42

Skema untuk mendapatkan tebakan awal

untuk u I, j , j=1,2,...,J dapat

diperoleh dengan menggunakan rumus forward difference (

u i +1, j − u i , j ∂u ) i , j ≈ h ∂ x

,

sehingga persamaan (3.6.1.2) dapat didekati dengan

u I, j

− g1, j p

= α 4u I, j + β 4 , j=1,2,...,J

atau

u I, j

=

g1, j

+ pβ 4

1 − pα 4

.

(3.6.1.3)

Persamaan (3.6.1.3) selanjutnya digunakan untuk mendapatkan nilai awal u I, j , yaitu (0)

u I, j

=

g1, j

+ aβ 4

1 − aα 4

, j=1,2,3...,J.

(3.6.1.4)

Nilai awal iterasi untuk u pada sudut batas, pada tulisan ini menggunakan nilai u

3.6.2

(0)

rata-rata dari kedua batas yang membentuk sudut tersebut.

Skema untuk Menentukan Nilai Awal Iterasi pada Titik Grid Bagian Dalam Domain

Titik grid yang berada pada bagian dalam domain berindeks i=2,3,...I-1, j=2,3,...,J-1, menggunakan skema beda hingga PD Poisson (3.5.3), yaitu u i −1, j u i , j

=

+ u i+1, j h2

+

u i . j −1

+ u i, j +1

k 2 1 

− f i , j

 1 +   h 2 k 2 

2

dengan i=2,3,...I-1, j=2,3,...,J-1.

43

Jika h=k , didapat

u i , j

=

ui

1, j +

−

ui

1, j +

+

u i . j

1 +

−

u i , j

1 −

+

2

h f i , j

.

4

(3.6.2.1)

Berdasarkan persamaan (3.6.2.1) dibentuk skema pertama untuk menentukan nilai awal iterasi u

(0)

pada titik grid bagian dalam domain, yaitu

u

( 0) i , j =

u i(01), j −

+

u i( 01), j +

+

u i(.0 j) 1

+

−

u i(,0 j) 1

−

+

h 2 f i , j

4

.

(3.6.2.2)

dengan i=2,3,…,I-1, j=2,3,…,J-1, dan skema kedua, yaitu

u

( 0) i , j =

( 0) + −1, j −1

ui

( 0) + +1, j −1

ui

(0) + −1. j +1

ui

( 0) − +1, j +1

ui

2

2h f i , j

4

. (3.6.2.3)

dengan i=2,3,…,I-1, j=2,3,…,J-1.

(0)

ui,j+1

ui-1 ,j

(0)

ui-1 ,j

(0)

(0)

ui+1,j+1

(0)

ui,j

(0)

ui+1,j

ui,j

h

2h

(0)

ui,j-1

Skema I

Gambar 3.6.

(0)

2

(0)

(0)

ui-1 ,j-1

ui+1,j-1 Skema II

Skema untuk menetukan nilai awal iterasi

44

3.6.3 Urutan Langkah Mengatur Nilai Awal Iterasi

Diasumsikan nilai awal iterasi selanjutnya akan ditentukan

u

(0)

u

(0)

pada batas domain telah diperoleh,

pada titik grid pada bagian dalam domain

dengan urutan langkah sebagai berikut 1) Langkah pertama (level grid ke-1) adalah mencari

u

( 0)

pada bagian tengah

domain menggunakan skema pertama atau dapat pula menggunakan skema kedua.

Gambar 3.7 Sketsa mencari tebakan awal pada level grid ke- 1

menggunakan skema I atau skema II 2) Langkah ke-2 (level grid ke-2) adalah

menghitung

u

(0)

lainnya

menggunakan skema II, dilanjutkan dengan menggunakan skema I.0

(a)

45

(b) Gambar 3.8 Sketsa mencari tebakan awal pada level grid ke- 2 menggunakan (a) skema II, (b) skema I

3) Langkah ke-3 (level grid ke- 3) dan seterusnya adalah identik dengan langkah ke-2.

Jika dilakukan

dihasilkan nilai awal

u

( 0)

sampai level grid ke-L maka akan L

2

sebanyak (2 +1) .

3.7 Bentuk Matriks PD Poisson dan Laplace 2D Diberikan PD Poisson 2D, yaitu 2

(

∂ u ∂ x

2

2

) i , j + (

∂ u ∂ y

2

) i , j = f i , j

(3.7.1)

dengan bentuk domain Ω ∆ = {( x i , y j ) | 0 ≤ x i ≤ p,0 ≤ y j ≤ q; i = 1,2, Λ , I; j = 1,2, Λ , J } .

dan syarat atas dapat bertipe Dirichlet , Neumann atau Robin. Bentuk domain seperti terlihat pada gambar berikut.

46

Gambar 3.9.

Domain dan partisinya

Jika digunakan rumus pendekatan 2

(

∂ u ∂ x

2

) i , j ≈

u i −1, j − 2u i , j + u i +1, j h2

dan 2

(

∂ u ∂ y

) ≈ 2 i , j

u i , j −1 − 2u i , j + u i , j +1 k 2

maka PD Poisson (3.7.1) dapat didekati dengan skema u i −1, j − 2u i , j + u i +1, j h2

+

u i , j −1 − 2u i , j + u i , j +1 k 2

= f i , j

atau

47

(−

1

)ui , j k 2

1

−

(

+ −

1 h

)ui 2

1, j

−

+

2(

1 h2

+

1

)ui , j k 2

(

+ −

1 h

)u i 2

1, j

+

(

+ −

1 k 2

)u i, j

1

+

f

= − i j . ,

(3.7.2)

Untuk h = k , didapat (−1)u i , j

( 1)u i

1 + −

−

1, j +

−

4ui , j

( 1)ui

+ −

( 1)u i , j

1, j + −

+

h 2 f i , j .

1 =−

+

(3.7.3)

Persamaan (3.7.2) dan (3.7.3) berlaku untuk i=2,3,...I-1; j=2,3,...,J-1. Sedangkan untuk i=1, i=I, j=1 dan j=J merupakan indeks-indeks pada batas domain, jadi persamaan (3.7.2) dan (3.7.3) harus disesuaikan dengan syarat batas yang diberikan. Bentuk sistem persamaan linier Au=b pada masing-masing subbab berikut ini menggunakan h=k .

3.7.1

PD Poisson dan Laplace 2D dengan Syarat Batas Domain Seluruhnya Bertipe Dirichlet Syarat batas domain seluruhnya bertipe Dirichlet , berarti semua nilai u

pada batas domain yaitu u1, j , u I , j , u i ,1 , u i ,J i=1,2,...,I; j=1,2,...,J telah diketahui sehingga tidak dibutuhkan skema numerik untuk mencarinya. Nilai-nilai u yang dicari hanya terletak pada bagian dalam domain (titik interior ), yaitu u i , j dengan i=2,3,...,I-1, j=2,3,...,J-1. Jadi masalah ini hanya

menggunakan skema numerik (3.7.3) tanpa tambahan skema untuk mencari u pada batas domain.

48

Persamaan-persaman linier yang dihasilkan oleh persamaan (3.7.3) dapat dilihat pada Tabel 3.1. Tabel 3.1.

j

j=2

i=1

4u 2,2

i=2

(−1)u 2, 2

…

( 1)u 3, 2

+ −

4u 3, 2

+

(−1)u I

−

2, 2

+

+

i=2

(−1)u3, 2

+ −

…

i=1

=−

( 1)u 4,2

+ −

+ −

4u I-1,2

(−1)u 2, 2

…

h 2 f 2, 2

( 1)u 2,3

+ −

( 1)u 3, 4

u 2,1

+

u1,2

2

h f 3, 2

=−

4u 2,3

( 1)u I-1,3

( 1)u3,3

+ −

( 1)u 2,3

+

4u 3,3

2

h f I

+ −

=−

1, 2

−

=−

( 1)u 4,3

+ −

+ −

+

uI

+

u 3,1

( 1)u3,4

u I, 2

+

1,1

−

h 2 f 2,3

( 1)u 2, 4

+ −

+

u1,3 h 2 f 3,3

=−

…

(−1)u I

1, 2 +

−

(−1)u I

−

2, 3 +

4u I-1,3

( 1)u I-1, 4

+ −

2

h f I

=−

1, 3 +

−

u I ,3

…

(−1)u 2,J

−

2 +

4u 2,J

( 1)u3,J

1 + −

−

2

h f 2,J-1

1 =−

−

i=2

(−1)u3, J - 2 + (−1)u2, J -1 + 4u3, J -1 + (−1)u4, J -1

…

…

i=I-1

(−1)u I

j=J-1

+

…

i=1

i=I-1 ...

Persamaan linier

i

i=I-1

j=3

Persaman-persamaan linier

1, J - 2

−

( 1)u I

+ −

−

2, J -1

+

4u I-1,J-1

+

+

u1, J

1

−

2

h f 3,J -1 + u3,J

=−

( 1)u I,J

+ −

u 2,J

1

−

2

h f I

=−

1, J −1

−

+

u I-1, J

+

u I, J-1

49

Persamaan-persamaan linier pada Tabel 3.1, dapat dinyatakan dalam bentuk Au=b

(3.7.1.1)

dengan 4 - 1 Λ - 1 4 - 1   -1 4 -1  Ο   A =     -1   -1   -1

-1

    -1        4 -1  -1 4 -1   - 1 4 

-1

-1

u 2,2    u 3,2  Μ    u I−1,2  u   2,3  u 3,3  u =   Μ  u   I−1,3  Μ    u 3,J−1  Μ    u I−1,J−1

 − h 2 f 2 , 2 + u 2 , 1 + u 1 , 2    2  − h f 3 , 2 + u 3 , 1     Μ   − h 2 f I − 1 , 2 + u I − 1 , 1 + u I , 2     − h 2 f 2 , 3 + u 1 , 3  .   − h 2 f 3 , 3   b =    Μ   − h 2 f  + u I −1 , 3 I,3    Μ    2  − h f 3 , J -1 + u 3 , J     Μ   − h 2 f I − 1 , J − 1 + u I -1 , J + u I , J -1 

50

Matriks u berukuran mx1 dengan m=(I-2)x (J-2), matriks A berukuran mxm, dan matriks b berukuran mx1. Hal yang menarik adalah matriks A dapat dinyatakan dalam bentuk block tridiagonal yang mempunyai kesesuaian dengan banyaknya nilai-nilai u yang

dicari pada arah sumbu X dan pada arah sumbu Y. Hal ini diilustasikan pada Gambar 3.10.

Gambar 3.10.

Ilustrasi letak nilai-nilai u pada domain

Bentuk matriks A pada persamaan (3.7.1.1) dapat ditulis menjadi matriks block tridiagonal, yaitu T - I - I T - I   - I T - I  Ο   A =        

           - I T - I  - I T - I   - I T 

51

dengan 4 - 1 - 1 4 - 1   -1 4 -1  Ο   T =        

-1

           4 -1  -1 4 -1   - 1 4 

matriks T berukuran m1 x m1 dengan m1=I-2 ( m1 sama dengan banyaknya u yang dicari pada arah sumbu X ), matriks I adalah matriks identitas yang berukuran sama dengan matriks T , matriks A terdiri dari susunan matriks sebanyak m2 x m2 dengan m2 = J - 2 ( m2 sama dengan banyaknya u yang dicari searah sumbu Y). Masing-masing matriks yang menyusun matriks A berukuran m1 x m1. Jadi matriks A berukuran m x m dengan m=m1 m2=(I-2)(J-2). Bentuk block tridiagonal tersebut memudahkan untuk mengidentifikasi sifat-sifat matriks A (sifat simetri, dominan diagonal, dan lain-lain), disamping itu memudahkan pemrograman untuk mencari nilai spectral radius matriks iterasinya. Nilai spectral radius

matriks iterasi tidak digunakan untuk

menyelesaikan sistem persamaan linier, tetapi nilai spectral radiusnya digunakan untuk menjamin metode iterasi konvergen atau tidak. Contoh untuk masalah PD Poisson dan Laplace 2D dengan syarat batas seluruhnya bertipe Dirichlet dapat dilihat pada Bab IV melalui Contoh 4.1 dan Contoh 4.2.

52

3.7.2

PD Poisson dan Laplace 2D dengan Syarat Batas Terdiri dari Tipe Dirichlet

dan Tipe Neumman

Diberikan PD Poisson ∂ u ( x, y ) 2

∂ x

2

∂ u ( x, y ) 2

+

∂ y

2

= f ( x, y )

dengan domain Ω = {( x, y ) | 0 ≤ x ≤ p,0 ≤ y ≤ q} dan syarat batas u( x,0) = g 1 ( x )

 ∂u  = β   2  ∂ x  1, y

u(0 ,y) =

g 3 ( y )

 ∂u     = β 4  ∂ y  ,1 x

dengan diberikan fungsi

g1 ( x ), g 3 ( y)

Ilustrasi letak nilai nilai

u

, dan β 1 , β 2 bilangan konstan.

yang dicari pada masalah ini dapat dilihat pada

Gambar 3.11.

Gambar 3.11.

Ilustrasi letak nilai-nilai

u

pada domain

53

Sistem persamaan linier Au=b untuk kasus ini, setelah diterapkan skema

u

pada batas-batas domain maka akan diperoleh matriks A dengan bentuk

A

T I - I T - I   - I T - I  Ο   =        

- I

T

- I

           - I  - I  T  - 2 I T 

dengan 4 - 1 - 1 4 - 1   -1 4 -1  Ο   T =        

matriks

T

berukuran

m1

x

m1

dengan

yang dicari pada arah sumbu X ),

m1

-1

           4 -1  -1 4 -1   -2 4 

= I-1 ( m1 sama dengan banyaknya

u

matriks I adalah matriks identitas yang

berukuran sama dengan matriks T , dan matriks A tersusun dari matriks sebanyak m2

x

m2

dengan

m2

= J-1 ( m2 sama dengan banyaknya

u

yang dicari searah

sumbu Y). Masing-masing matriks yang menyusun matriks A berukuran m1 x Jadi matriks A berukuran m x m dengan

m1.

m= m1 m2 = (I-1)(J-1).

Contoh PD Poisson dengan batas domain terdiri dari tipe Dirichlet dan tipe Neumman diberikan pada Bab IV melalui Contoh 4.3.

54

3.7.3

PD Poisson dan Laplace 2D dengan Syarat Batas Domain Seluruhnya Bertipe Neumann

Diberikan PD Poisson

∂ 2 u( x, y ) ∂ 2 u( x, y ) + = f ( x, y ) ∂ x 2 ∂ y 2 dengan domain

Ω = {( x, y) | 0 ≤ x ≤ p,0 ≤ y ≤ q}

dan syarat batas

 ∂u     = β 1  ∂ y  x,0  ∂u     = β 2  ∂ y  x,q  ∂u  = β   3  ∂ x  0, y  ∂u     = β 4 ∂ y   p, y Sistem persamaan linier Au=b untuk kasus ini, setelah diterapkan skema u pada batas-batas domain selanjutnya akan diperoleh matriks A dengan bentuk

T - I      A =        

- 2 I T - I

- I

- I

T

Ο

- I

T

- I

           - I  T - I   - 2 I T 

55

dengan

4 - 2 - 1 4 - 1   - 1 4 -1  Ο   T =        

-1

           4 -1  - 1 4 -1   -2 4 

matriks T berukuran m1 x m1 dengan m1 = I ( m1 sama dengan banyaknya u yang dicari pada arah sumbu X ), matriks I adalah matriks identitas yang berukuran sama dengan matriks T , sedangkan susunan block tridiagonal matriks A terdiri dari matriks sebanyak m2 x m2 dengan m2 = J ( m2 sama dengan banyaknya u yang dicari pada arah sumbu Y). Masing-masing matriks yang menyusun matriks A berukuran m1 x m1. Jadi matriks A berukuran m x m dengan m= m1 m2 = I.J.

Contoh PD Poisson dengan batas domain seluruhnya bertipe Neumann diberikan pada Bab IV melalui Contoh 4.4.

56

3.7.4

PD Poisson dan Laplace 2D dengan Syarat Batas Terdapat Tipe Robin

Diberikan PD Poisson 2D

∂ 2 u ( x, y ) ∂ 2 u( x, y ) + = f ( x, y ) ∂ x 2 ∂ y 2 dengan domain

Ω = {( x, y ) | 0 ≤ x ≤ p,0 ≤ y ≤ q}

dan syarat batas u(x,0)= g1 ( x )

 ∂u     = β 2 ∂ y   x,1  ∂u  = α u +   3  ∂ x  0, y  ∂u  = α u +   4  ∂ x  1, y

β 3 β 4

dengan β 2 , β 3 , β 4 , α 3 , α 4 bilangan konstan. Ilustrasi letak nilai-nilai u yang akan dicari untuk masalah ini dapat dilihat pada Gambar 3.12.


57

Sistem persamaan linier Au=b untuk kasus ini, setelah diterapkan skema u pada batas-batas domain maka akan diperoleh matriks A dengan bentuk T - 2 I  - I T - I   T - I - I  Ο   A =        

           T - I - I  T - I  - I  T  - 2 I

dengan Φ - 2 - 1 4 - 1   - 1 4 -1  Ο   T =        

-1

           4 -1  -1 4 -1   - 2 Ψ 

dengan Φ = 4 + 2hα 3 , Ψ = 4 − 2hα 4 , matriks T berukuran m1xm1 dengan m1 = I ( m1 sama dengan banyaknya u yang dicari pada arah sumbu X ), matriks I adalah matriks identitas yang berukuran sama dengan matriks T , sedangkan susunan block tridiagonal matriks A terdiri dari matriks sebanyak m2xm2 dengan m2 = J-1

( m2 sama dengan banyaknya u yang dicari searah sumbu Y). Masing-masing matriks yang menyusun matriks A berukuran m1 x m1. Jadi matriks A berukuran mxm dengan m= m1 m2 = I (J-1). Contoh untuk masalah ini diberikan pada Bab

IV melalui Contoh 4.5 dan Contoh 4.6.

58

3.8 Sifat Skema Beda Hingga PD Poisson 2D dan Laplce 2D

Secara umum skema beda hingga persamaan differensial Poisson dan Laplace 2D pada tulisan ini dapat dinyatakan dalam bentuk a1u i , j

1 +

−

a2ui

1, j +

−

a 3u i , j

+

a 4ui

1, j +

+

a 5u i , j

f i , j

1 = −

+

(3.8.1)

a5 ui ,j+1

a2 ui-1, j

a3 ui,j

a4 ui+1, j

a1 ui,j-1

Gambar 3.13. Ilustrasi letak nilai-nilai u dan koefisiennya

dengan a1 , a2 , a3 , a4 , a5 adalah koefisien yang mempunyai sifat sebagai berikut I) Jika syarat batas domain seluruhnya bertipe Dirichlet maka semua sifat berikut dipenuhi (a) a1 <0, a2 <0, a3 >0, a4 <0, a5 <0. (b) Jika persamaan (3.8.1) dinyatakan dalam bentuk matriks Au=b maka komponen diagonal matriks A adalah a3 . (c) Bersifat simetri, yaitu a1 = a5, a2 = a4 .

59

(d) Bersifat dominan diagonal , yaitu |a3|

•

|a1| + |a2| + |a4| + |a5|

dengan bentuk ketaksamaannya (strong dominan diagonal ) dipenuhi oleh paling sedikit satu titik grid. (e) Jika himpunan persamaan linier yang dihasilkan oleh persamaan (3.8.1) lebih dari satu persamaan linier maka tidak ada persamaan linier yang dapat diselesaikan secara terpisah dari persamaan lainnya (irreducibility). II)

Jika batas domain seluruhnya bertipe Neumann maka sifat I(a), I(b), dan I(e) dipenuhi, sedangkan

sifat dominan diagonal (Id) hanya dipenuhi

bentuk kesamaannya saja, yaitu |a3| = |a1| + |a2| + |a4| + |a5|

III)

Jika pada batas domain hanya terdiri dari tipe Neumann dan tipe Dirichlet maka sifat I(a), I(b), I(d), dan I(e) dipenuhi, sedangkan sifat I(c) tidak dipenuhi.

IV)

Jika pada syarat batas terdapat tipe Robin dengan bentuk

∂u ∂n

= α u + β ,

α ≠ 0 , maka sifat I(b) dan I(e) dipenuhi. Karena a3 bergantung pada nilai α maka nilai a3 dapat bernilai nol atau negatif, jadi sifat I(a) dapat tidak dipenuhi. Sedangkan sifat dominan diagonal I(d) dapat pula tidak dipenuhi tergantung pada tanda (positif atau negatif) nilai α pada syarat batas Robin.

60

3.9. Urutan Penyelesaian Secara umum urutan penyelesaian persamaan differensial Poisson dan Laplace dalam tulisan ini dapat dilihat pada bagan berikut.

Masalah syarat batas Persamaan Differensial Poisson dan La lace 2D Rumus Beda Hingga

Skema Numerik

Sistem Persamaan Linier (SPL )

Metode Iterasi SOR dengan Pengaturan nilai awal iterasi

Solusi SPL

Solusi SPL

≈

Solusi Pers. Differensial Poisson dan Laplace 2D

Gambar 3.14. Skema penyelesaian masalah

61

BAB IV PENGGUNAAN METODE DAN APLIKASINYA

Masalah syarat batas PD Poisson dan Laplace 2D pada contoh-contoh berikut diambil domain

Ω = {( x, y ) | 0 ≤ x ≤ 1, 0 ≤ y ≤ 1} . Jika masalah syarat batas tersebut didefinisikan pada domain

Ω1 = {( x, y ) | a ≤ x ≤ b, c ≤ y ≤ d } maka dapat ditransformasi mejadi

γ.∂2 ũ / ∂ x*2 + Ψ.∂2ũ / ∂ y* 2

= f ( x*,y*)

yang didefinisikan pada

Ω 2 = {( x*, y*) | 0 ≤ x* ≤ p, 0 ≤ y* ≤ q} dengan u ( x, y ) =

ũ

( x*, y*)

γ = ( p /(b − a )) 2 , Ψ = (q /(d − c)) 2 . Pembicaran mengenai transformasi koordinat ini selanjutnya akan dibicarakan pada akhir bab ini.

62

4.1 PD Poisson 2D dengan Syarat Batas Seluruhnya Bertipe

Dirichlet

Contoh 4.1 (Abayomi, 2003)

Diberikan PD Poisson

∂ 2u ( x, y ) ∂ 2u ( x, y ) + =4 ∂ x 2 ∂ y 2

(4.1.1)

dengan domain Ω = {( x, y ) | 0 ≤ x ≤ 1,0 ≤ y ≤ 1} dan syarat batas

⎧ x 2 , y = 0 ⎪ 2 ⎪1 + y , x = 1 u=⎨ 2 ⎪1 + x , y = 1 ⎪ y 2 , x = 0 ⎩ . Persoalan ini telah diketahui mempunyai solusi eksak u ( x, y ) = x 2 + y 2 , dalam bentuk grafik contour tampak seperti pada Gambar 4.1.

Gambar 4.1. Grafik contour u ( x, y ) = x 2

+ y 2

Solusi eksak tersebut selanjutnya akan dibandingkan dengan solusi numerik dari masing-masing perlakuan berikut ini.

63

4.1.1 PD Poisson 2D dengan Nilai Awal Iterasi u ( 0)

=0

Persoalan ini diselesaikan dengan grid IxJ=33x33, jadi h=1/(I-1)=1/32, dan k =1/(J-1)=1/32. Nilai awal iterasi dipilih u ( 0) = 0 , parameter SOR ω = 1 .2 , nilai toleransi

ε

= 0. 01 . Dari Lampiran 1 diperoleh bahwa spectral radius matriks

iterasi untuk masalah ini adalah

ω

= 0.9880 . Karena

ω

<1 maka metode iterasi

SOR dijamin konvergen. Seperti terlihat pada Lampiran 2, proses iterasi konvergen, dibutuhkan 51 iterasi untuk mencapai kekonvergenan. Grafik solusi numeriknya tampak pada Gambar 4.2.


64

Solusi numerik, solusi eksak dan error global (e) untuk masalah ini dapat dilihat pada Lampiran 2A, sebagian diperlihatkan pada Tabel 4.1.

Tabel 4.1.

No

1 2 3 4 ... 640 641 642 643 644 645 646 647 648 649 650 ... 1085 1086 1087 1088 1089

x

Perbedaan solusi numerik dan solusi eksak

y

Solusi numerik Solusi eksak 2 2 u numerik u eksak = x + y

Error global e = u numerik − u eksak

0.00000 0.00000 0.00000 0.00000

0.00000 0.03125 0.06250 0.09375

0.00000000 0.00097656 0.00390625 0.00878906

0.00000000 0.00097656 0.00390625 0.00878906

0.00000000 0.00000000 0.00000000 0.00000000

0.59375 0.59375 0.59375 0.59375 0.59375 0.59375 0.59375 0.59375 0.59375 0.59375 0.59375

0.37500 0.40625 0.43750 0.46875 0.50000 0.53125 0.56250 0.59375 0.62500 0.65625 0.68750

0.06989351 0.06913434 0.07290868 0.08170180 0.09607642 0.11665160 0.14407205 0.17896975 0.22192028 0.27339754 0.33373044

0.49316406 0.51757813 0.54394531 0.57226563 0.60253906 0.63476563 0.66894531 0.70507813 0.74316406 0.78320313 0.82519531

-0.42327055 -0.44844378 -0.47103663 -0.49056383 -0.50646264 -0.51811403 -0.52487326

1.00000 1.00000 1.00000 1.00000 1.00000

0.87500 0.90625 0.93750 0.96875 1.00000

1.76562500 1.82128906 1.87890625 1.93847656 2.00000000

1.76562500 1.82128906 1.87890625 1.93847656 2.00000000

-0.52610838

-0.52124378 -0.50980559 -0.49146487 0.00000000 0.00000000 -0.00000000 0.00000000 0.00000000

Keseluruhan data error global, apabila dinyatakan dalam bentuk grafik tampak seperti pada Gambar 4.3. Nilai absolut terbesar dari error global yaitu maks | e( x, y ) | =0.52610838 terjadi pada koordinat x=0.59375, y=0.59375. Jika error global tersebut dibandingkan dengan solusi eksak pada koordinat tersebut, yaitu u = 0.70507813 maka error global e = -0.52610838 dianggap besar.

65

Gambar 4.3.

Grafik error global

Error global ( e) biasanya biasanya dapat dapat diperkecil diperkecil dengan dengan toleransi

ε

yang lebih kecil (misalnya

ε

=

10 6 ), tetapi −

cara memilih memilih nilai sebagai konsekuensi

banyaknya iterasi untuk mencapai kekonvergenan kekonvergenan umumnya umumnya

akan meningkat. meningkat.

Dilain pihak, banyaknya iterasi untuk mencapai kekonvergenan (

v

0

) dipengaruhi

pula oleh parameter SOR ( ω ), sebagai contoh dapat dilihat pada Gambar 4.4. 1200

1000

( ) v ) 0 n ( i s a r e t I k a y n a B

800

600

400

200

0 0

Gambar 4.4.

0. 5

1

1.5

2

Parameter SOR (ω

2. 5

)

Pengaruh parameter SOR ( ù ) terhadap banyak iterasi untuk mencapai kekonvergenan ( v0) menggunakan toleransi å =10 =10

-6

66

4.1.2 PD Poisson 2D 2D dengan Nilai Awal Iterasi u ( 0) Dimodifikasi

Nilai u (0)

=

ω

=

1.2 pada uji coba sebelumnya menggunakan tebakan awal

0 dianggap kurang baik karena dibutuhkan iterasi yang banyak untuk

mencapai kekonvergenan, yaitu dengan

ε

=0.01

dibutuhkan 51 iterasi untuk -6

konvergen dan data grafik pada Gambar 4.4 bahwa untuk å =10 =10 dibutuhkan 667 iterasi untuk mencapai kekonvergenan. Tebakan awal iterasi selanjutnya dimodifikasi dengan memanfaatkan nilainilai u pada syarat batas Dirichlet menggunakan urutan langkah seperti yang dikemukakan pada subbab 3.6.2, ternyata proses iterasi konvergen konvergen pada iterasi iterasi ke- 1, nilai maks | e( x, y ) | = 2.2204 x 10

-16

terjadi pada koordinat x = 0.90625,

y = 0.96875 (Lampiran 2B). Hasil Hasil ini tentunya jauh lebih baik.

Grafik solusi numerik numerik dan error global global pada uji coba ini dapat dilihat pada Gambar 4.5 dan 4.6.


67

Gambar 4.6.

Grafik error global dengan parameter SOR

ω

=

1.2

Hal yang menarik adalah error global dapat tidak tampak. Hal ini terjadi jika digunakan

ω

=

1 (dapat dilihat pada Lampiran 2C).

Nilai absolut terbesar dari error global pada

ω

=

1 tetap memperlihatkan

nilai nol sekalipun digunakan grid lebih banyak yang memungkinkan munculnya error pembulatan yang besar. Sebagai contoh pada Lampiran 2D dipilih grid IxJ=2049x2049. Grafik error globalnya seperti terlihat pada Gambar 4.7.

Gambar 4.7.

Grafik error global dengan parameter SOR

ω

=

1

68

4.2 PD Laplace 2D dengan Syarat Batas Seluruhnya Bertipe Dirichlet Contoh 4.2 (Abayomi, 2003) Diberikan PD Laplace 2D

∂ 2 u ( x, y ) ∂ 2 u ( x, y ) + =0 ∂ x 2 ∂ y 2 dengan domain Ω = {( x, y ) | 0 ≤ x ≤ 1,0 ≤ y ≤ 1} dan syarat batas seluruhnya bertipe Dirichlet yaitu

 x 2 + x, y = 0  2 2 + 3 y − y − 2 y, x = 1 u= 2  x − 3 y + x − 3, y = 1 − y 2 − 2 y, x = 0 

2

2

Telah diketahui solusi eksak untuk masalah ini adalah u(x,y) = x –3xy- y + x -2y, dalam grafik contour tampak seperti pada Gambar 4.8.

Gambar 4.8.

Grafik

2

2

u(x,y) = x –3xy- y + x - 2y

69

Solusi eksak tersebut selanjutnya akan dibandingkan dengan solusi numerik yang diperoleh dari masing-masing uji coba berikut ini.

4.2.1 PD Laplace 2D dengan Nilai Awal Iterasi

u

( 0)

=

0

Uji coba ini menggunakan grid IxJ=33x33, jadi h=k =1/32, parameter SOR dipilih ω 1.2 , diperoleh nilai =

adalah µ ω

=

0.9880 <1

spectral radius matrik iterasi untuk masalah ini

(Lampiran 1A),

jadi

metode iterasi SOR dijamin

konvergen. (0) Selanjutnya tebakan awal u dipilih seperti yang umum digunakan yaitu u

dipilih nilai toleransi ε

=

=

0 , dan

0.01 .

Output komputer pada Lampiran 3A memperlihatkan bahwa masalah ini membutuhkan 91 iterasi untuk mencapai kekonvergenan dan nilai absolut error global

terbesar

yaitu

maks | e( x, y ) | 0.61238212 =

terjadi

pada

koordinat

x=0.46875, y=0.56250. Grafik solusi numerik dan error global tampak pada Gambar

4.9 dan 4.10.

. Gambar 4.9. Grafik solusi numerik

70

Gambar 4.10.

Grafik error global

Nilai absolut terbesar dari error global yaitu maks | e( x, y ) | 0.61238212 =

terjadi pada koordinat x=0.46875, y=0.56250. Diketahui solusi eksak masalah ini 2 2 adalah u(x,y) = x –3xy- y + x - 2y, jadi untuk x=0.46875, y=0.56250 menghasilkan

u(x,y)

≈

-1.405517570. Jika dibandingkan dengan solusi eksak ini tentunya nilai

error global tersebut dianggap besar. Seperti yang umum diketahui bahwa untuk memperkecil error global biasanya dipilih nilai toleransi ε yang lebih kecil misalnya ε

=

6

10 , atau sekaligus −

dengan memperbanyak partisi domain agar dihasilkan ukuran kisi ( h dan k ) yang relatif lebih kecil, tetapi sebagai konsekuensinya, banyaknya iterasi untuk mencapai kekonvergenan umumnya akan meningkat. Dilain pihak, banyaknya iterasi untuk mencapai kekonvergenan dapat pula dipengaruhi oleh parameter

ω

, sebagai contoh

dapat dilihat pada Gambar 4.11.

71

1200

1000

(v 0 ) ) n ( i s a r e t I

k a y n a B

800

600

400

200

0 0

0.5

1

1. 5

2

2.5

( )

P a ra m e t e r S O R ω

Gambar 4.11. Hubungan parameter SOR terhadap banyak iterasi untuk mencapai kekonvergenan menggunakan

4.2.2 1 PD Laplace 2D dengan Nilai Awal Iterasi

u

(0)

ε

=

10

6

−

Dimodifikasi

Tebakan awal iterasi selanjutnya dimodifikasi sesuai urutan langkah yang dikemukakan pada subbab 3.6.3, ternyata proses iterasi dapat konvergen pada iterasi ke-1, (Lampiran 3B). Grafik solusi numeriknya tampak pada Gambar 4.12.


72

Output program komputer pada Lampiran 3B memperlihatkan bahwa maks | e( x, y ) |= 4.4409 × 10

-16

. Adapun grafik error global untuk masalah ini dapat

dilihat pada gambar berikut.

Gambar 4.13.

Grafik error global dengan

ω =

1.2 , grid 33x33

Hal yang menarik seperti pada contoh 4.1, bahwa jika tebakan awal dimodifikasi dan digunakan

ω =

1 , ternyata nilai error global sama dengan nol,

(Lampiran 3C). Hal ini tidak berubah sekalipun grid diperbanyak, (Lampiran 3D).

Gambar 4.14.

Grafik error global dengan

ω =

1 , grid 2049x2049

73

4. 3 PD Poisson 2D dengan Domain Mempunyai Batas bertipe Dirichlet dan Batas Lainnya Bertipe Neumann Contoh 4.3 Diberikan PD Poisson

∂ 2 u ( x, y ) ∂ 2 u ( x, y ) + =4 ∂ x 2 ∂ y 2 dengan domain

Ω = {( x, y) | 0 ≤ x ≤ 1,0 ≤ y ≤ 1}

dan syarat batas 2

u(0 ,y) = y

 ∂u  =    ∂ x  1, y

2 2

u( x,0) = x

 ∂u     =  ∂ y  x,1

2 2

2

yang telah diketahui mempunyai solusi eksak u( x,y) = x + y . Persoalan ini diselesaikan secara numerik dengan menggunakan grid IxJ=33x33, parameter SOR dipilih ù =1, diperoleh nilai spectral radius matrik iterasi µ ω

= 0.9976<1,

(Lampiran 4A). Karena µ ω <1, maka metode iterasi SOR dijamin

konvergen. Tebakan awal

u

( 0)

tidak menggunakan nilai nol untuk seluruh titik grid,

tetapi dimodifikasi dengan memanfaatkan nilai-nilai

u

pada batas bertipe

Dirichlet ,

-6

nilai toleransi dipilih å =10 . Seperti terlihat pada Lampiran 4B, proses iterasi konvergen pada iterasi ke- 1, dan error global menunjukkan nilai nol. Grafik solusi numerik tampak seperti pada Gambar 4.15.

74


4.4 Contoh Kasus, Metode SOR Menghasilkan Solusi yang Tidak Tunggal Contoh 4. 4 Diberikan PD Laplace 2D

∂ 2 u ( x, y ) ∂ 2 u( x, y ) + =0 ∂ x 2 ∂ y 2 dengan domain

Ω = {( x, y) | 0 ≤ x ≤ 1,0 ≤ y ≤ 1}

syarat batasnya seluruhnya bertipe Neumann, yaitu

 ∂u     = 1  ∂ y  x ,0  ∂u     = 1 y ∂   x,1  ∂u  = 1    ∂ x  0, y  ∂u     = 1 y ∂   1, y

75

Persoalan ini diselesaikan dengan menggunakan grid IxJ=5x5. Seperti terlihat pada Lampiran 5A bahwa jika digunakan parameter SOR ù =1.5 maka spectral radius matrik iterasinya adalah µ ω µ ω

=

=

0.5000, sedangkan untuk ù =1.333 dihasilkan

0.7286, jadi metode iterasi SOR dengan menggunakan nilai-nilai parameter

tersebut dijamin konvergen. Nilai toleransi dipilih ε

=

0.01 . Tebakan awal iterasi digunakan nilai yang

berbeda untuk masing-masing parameter ω

=

1.5 , kedua

u

( 0)

=

SOR, yaitu pertama

u

( 0)

=

0 dengan

10 dengan ω 1.333 . Seperti terlihat pada Lampiran 5B dan =

5C proses iterasi dapat konvergen. Solusi numerik masing-masing perlakuan ini sekalipun memperlihatkan pola grafik yang sama seperti terlihat pada Gambar 4.16 dan 4.17, tetapi interval nilai

u

menunjukkan perbedaaan yang mencolok. Hal ini

menjadi petunjuk bahwa metode SOR untuk masalah ini menghasilkan solusi yang tidak tunggal

Gambar 4.16.

Solusi numerik dengan

u

( 0)

=

0, ω 1.5 =

76

Gambar 4.17. Solusi numerik dengan u

( 0)

=

10, ω

=

1. 333

4.5 Persamaan Differensial Laplace dengan Batas Domain Terdapat Tipe Dirichlet, Neumman dan Robin Solusi eksak pada contoh-contoh berikut ini tidak tersedia, jadi error yang dipergunakan bukanlah error global, tetapi error relatif dari bentuk sistem persamaan linier Au=b, yaitu R

(v )

=

Au

(v)

−

b

dengan v adalah nomor iterasi v=1,2,3,… .

Contoh 4.5 Diberikan PD Laplace 2D ∂ u ( x, y ) 2

∂ x

2

∂ u ( x, y ) 2

+

∂ y

dengan domain Ω = {( x, y ) | 0 ≤ x ≤ 1,0 ≤

2

=0

y ≤ 1}

77

dan syarat batas u( x,0)=

1

 ∂u     = 2 y ∂   x ,1  ∂u  = 3 u + 1    ∂ x  0, y  ∂u  = −4 u + 1    ∂ x  1, y Masalah ini diselesaikan dengan menggunakan grid IxJ=33x33, jadi diperoleh h=k = 1/32. Parameter SOR dipilih ω = 1.2 , diperoleh nilai spectral radius matrik

iterasi yaitu µ ω

=

0.99551 , (Lampiran 6A). Karena µ ω

<

1 maka metode

SOR

( 0) dijamin konvergen. Nilai awal iterasi u dimodifikasi dengan memanfaatkan nilai u

6

−

yang diketahui pada syarat batas Dirichlet , nilai toleransi dipilih ε = 10 . Proses iterasi konvergen pada iterasi ke- 258,

(Lampiran 7).

Grafik solusi numeriknya

tampak pada Gambar 4.18.

Gambar 4.18.

Solusi numerik

78

Realitas solusi numerik tersebut dapat dilihat melalui sketsa nilai-nilai u sesuai letaknya pada domain, seperti pada Gambar 4.19.

Gambar 4.19.

Sketsa solusi numerik

Terlihat dari Gambar 4.19 bahwa nilai-nilai u pada titik interior mengikuti skema beda hingga PD Laplace 2D, sedangkan nilai-nilai u pada batas domain mengikuti

syarat batas yang diberikan. Adapun error relatif untuk masalah ini seperti terlihat pada Gambar 4.20. Nilai absolut terbesar dari error relatifnya adalah 9.2821 x 10

-4

terjadi pada koordinat x=0.000, y=0.71875, (Lampiran 7).

X

Y

Gambar 4.20.

Grafik error relatif

79

Apabila nilai absolut error relatif terbesar dari tiap-tiap iterasi di plot sesuai nomor iterasinya maka diperoleh grafik yang memperlihatkan fenomena error relatif yang teredam seiring dengan meningkatnya jumlah iterasi. Terlihat pada Gambar 4.21 bahwa nilai absolut error relatif terbesar dari tiap-tiap iterasi setelah sekitar 10 iterasi pertama, mulai mendekati nol.

Gambar 4.21.

Error teredam

Contoh 4.6

Diberikan PD Laplace 2D dengan domain

Ω = {( x, y ) | 0 ≤ x ≤ 1,0 ≤ y ≤ 1}

dengan syarat batas bertipe Dirichlet, Neumann, dan Robin, yaitu u( x,0) = 1

 ∂u     = 2  ∂ y  x,1  ∂u  = -3u+1    ∂ x  0, y  ∂u     = 4u +1 y ∂   1, y

80

Masalah ini diselesaikan dengan menggunakan grid 33x33, parameter SOR ω 1.2 , diperoleh spectral radius matriks iterasi SOR, yaitu µ ω = 1.0101, =

(Lampiran 6B). Karena µ ω >1, secara teoritis metode SOR tidak dapat konvergen. Nilai toleransi dipilih ε

=

10 6 , proses iterasi untuk masalah ini dapat −

dilihat pada Lampiran 8. Banyaknya titik grid yang konvergen cenderung tidak bertambah sekalipun iterasi diperbanyak. Apabila nilai

absolut error relatif

terbesar dari tiap-tiap iterasi di plot sesuai nomor iterasinya maka diperoleh grafik seperti pada Gambar 4.22 .

. Gambar 4.22.

Error relatif tidak teredam

Terlihat pada Gambar 4.22 bahwa setelah 10 iterasi pertama nilai absolut error relatif terbesar semakin meningkat bila jumlah iterasi diperbanyak. Error relatif terlihat tidak dapat diredam. Bentuk grafik mulai cekung keatas, yang berarti error relatif mengalami percepatan tumbuh. Hal ini menunjukkan gejala bahwa metode iterasi SOR tidak dapat konvergen.

81

4.6 Konduksi Panas Steady State (Constantinides, 1987)

Persamaan konduksi panas steady state pada suatu lempengan logam merupakan PD Poisson 2D, yaitu ∂ u ( x, y ) 2

∂ x

dengan fungsi sumber

2

∂ u ( x, y ) 2

+

∂ y

2

=−

Q ' ( x, y )

(4.6.1)

K

f(x,y) = -Q’( x, y)/ K( x, y), K adalah konduktifitas termal

logam tersebut, dan Q’(x,y) menyatakan besarnya panas yang dibangkitkan oleh sumber panas internal pada logam itu. Bila logam tersebut tidak memiliki sumber panas internal berarti Q’(x,y)=0, maka persamaan (4.6.1) menjadi persamaan differensial Laplace 2D.

Contoh 4.7

Sebuah lempengan tipis terbuat dari logam berukuran 1 ft x 1 ft . Pada o

o

batas kiri dan atasnya diberikan suhu tetap masing-masing 500 F dan 1000 F sedangkan batas-batas lainnya diisolasi sempurna. Penyekatan dilakukan sedemikian sehingga suhu lingkungan tidak mempengaruhi suhu lempengan tersebut. Pada masalah ini akan ditentukan suhu steady state pada lempengan. Karena logam tersebut tidak mempunyai sumber panas internal, jadi diselesaikan menggunakan PD Laplace 2D, yaitu

∂ u ( x, y ) 2

∂ x

2

∂ u ( x, y ) 2

+

∂ y

2

=0

82

dengan domain

Ω = {( x, y ) | 0 ≤ x ≤ 1,0 ≤ y ≤ 1}

dan syarat batas

 ∂u     = 0 ∂ y   x,0 o

u( x,1) = 1000 F u(0 ,y) = 500

o

F

 ∂u  = 0    ∂ x  1, y

Input dan output program komputer untuk menyelesaikan masalah ini

dapat dilihat pada Lampiran 9. Solusi numerik berupa suhu yang terjadi pada lempengan logam tersebut tampak seperti pada gambar berikut.

Gambar 4.23.

Suhu pada lempengan logam

Realitas solusi numerik tersebut dapat dilihat dengan menempatkan nilainila suhu sesuai letaknya pada domain (lempengan logam) seperti tampak pada sketsa berikut.

83

750 = (500 + 1000) / 2

1000

750.00000

1000.00000

…

1000.00000

1000.00000

1000.00000

1000.00000

500.00000

750.00000

…

986.60542

986.76193

986.87259

986.96045

500.00000

651.17418

…

973.56972

973.78989

973.92112

973.96471

500.00000

604.69670

…

960.46870

960.79615

960.99132

961.05617

500.00000

579.25303

…

947.50324

947.93466

948.19186

948.27732

…

…

…

…

…

…

…

∂u

500.00000

513.23806

…

749.99997

751.33493

752.13583

754.27159

∂ x

500.00000

513.12740

…

748.66501

749.99997

750.80091

751.06789

500.00000

513.06146

…

747.86410

749.19902

749.99997

750.26695

500.00000

513.03955

…

747.59715

748.93204

749.73299

749.99997

500

∂u ∂ y

=0

ui , j =

ui −1, j + ui + 1, j + ui. j −1 + ui , j +1

749.19902 ≈

Gambar 4.24.

=0

4 747. 86410 + 749.99997 + 748.93204 + 749.99997 4

Sketsa suhu pada lempengan logam

Terlihat pada Gambar 4.24 bahwa suhu lempengan logam tersebut pada o

o

batas bagian kiri 500 F dan batas bagian atas 1000 F, sedangkan suhu pada batas bagian kanan dan bawah cenderung mengikuti suhu pada titik-titik grid terdekat sesuai arah normal pada masing-masing batas tersebut. Keadaan ini sesuai dengan syarat batas yang diberikan. Adapun suhu-suhu pada titik-titik grid yang berada pada bagian dalam domain (titik interior ) mengikuti skema beda hingga PD Laplace 2D.

Contoh 4.8

Sebuah lempengan tipis berukuran 1 ft x 1 ft , diketahui terbuat dari logam o

o dengan titik lebur 1200 F dan konduktifitas termalnya K = 10 BTU /(h.ft. F).

Lempengan logam dialiri arus listrik sehingga terbentuk sumber panas internal pada lempengan tersebut dengan panas yang dibangkitkan merata sebesar

84

3

Q’=200000 BTU /(h.ft ). Setiap batas lempengan tersebut dipertahankan bersuhu o

tetap 100 F. Penyekatan dilakukan sedemikian sehingga suhu pada lempengan logam tidak dipengaruhi oleh lingkungan. Melalui perlakuan ini akan ditentukan suhu steady state yang terjadi pada lempengan logam tersebut. Karena suhu lebur logam diketahui, maka bentuk lelehan pada lempengan dapat diperoleh dengan cara mem plot koordinat-koordinat yang suhunya melebihi suhu lebur. Lempengan logam ini mempunyai sumber panas internal, jadi digunakan PD Poisson 2D, yaitu

∂ u ( x, y ) 2

∂ x

2

∂ u ( x, y ) 2

+

∂ y

2

= f ( x, y )

dengan f ( x, y ) = −Q' / K = −20000 , domain Ω = {( x, y ) | 0 ≤ x ≤ 1,0 ≤ y ≤ 1} , syarat batas u( x,0) = 100 u( x,1) = 100 u(0, y) = 100 u(1, y) = 100.

Masalah ini diselesaikan dengan menggunakan grid IxJ=33x33. Kemudian domain dipartisi lebih banyak agar diperoleh aproksimasi bentuk lelehan yang lebih baik, digunakan grid IxJ= 1025x1025, (Lampiran 10). Solusi numerik berupa suhu-suhu pada lempengan logam dan aproksimasi bentuk lelehannya tampak pada gambar berikut.

85

Gambar 4.25.

Aproksimasi suhu dan bentuk lelehan

Terlihat pada Gambar 4.25 bahwa terdapat koordinat-koordinat pada o

lempengan logam yang suhunya melebihi suhu lebur logam (1200 F), mulai dari bagian tengah logam sampai pada batas lelehan. Batas lelehan bersuhu sama o

dengan suhu lebur logam (1200 F ). Grafik error relatif masalah ini dapat dilihat pada Gambar 4.26. Error relatif untuk masalah ini mempunyai satuan o

o

2

2

F / ft . Grafik 4.26 memperlihatkan

bahwa error relatif berkisar antara -0.1 F / ft sampai 0.0035

o

2

F / ft , yang berarti

o 2 nilai absolut error relatif terbesar sekitar 0.1 F / ft . Hal ini menunjukkan akurasi o 2 numerik yang baik karena perbedaaan kurang dari 0.1 F/ ft untuk masalah suhu

umumnya tidaklah terlalu dipersoalkan.

86

Gambar 4.26.

4.7

Error relatif suhu pada lempengan logam logam

Transformasi Koordinat

Diberikan PD Poisson ∂ u ( x, y ) 2

∂ x

2

∂ u ( x, y ) 2

+

∂ y

2

= f ( x, y )

(4.7.1)

dengan domain Ω 1 = {( x, y ) | a ≤ x ≤ b, c ≤ y ≤ d }

Domain Ω1 selanjutnya ditransformasi ke koordinat yang baru, yaitu Ω 2 = {( x*, y*) | 0 ≤ x* ≤ p, 0 ≤ y* ≤ q} .

Fungsi transformasi dipilih memenuhi persamaan garis lurus, yaitu x* = g1 ( x) = m 1 x + c1 .

87

Karena x = a, x* = g(a g(a) = 0, dan x = b, b, x* = g(b g(b) = p = p,, diperoleh

⎛ p ⎞ x − ⎛ p ⎞a ⎟ ⎜ ⎟ ⎝ b − a ⎠ ⎝ b − a ⎠

x* = ⎜

(4.7.2)

∂ x * p = . ∂ x b − a

(4.7.3)

jadi

Cara untuk mendapatkan mendapatkan persamaan (4.7.2) dan (4.7.3) digunakan pula untuk fungsi transformasi

y* = g 2 ( y) y) = m 2 y + c 2 sehingga diperoleh

⎛ q ⎞ y − ⎛ q ⎞c ⎟ ⎜ ⎟ ⎝ d − c ⎠ ⎝ d − c ⎠

Diketahui

x*, y*) y*) ũ( x*,

y* = ⎜

(4.7.4)

∂ y * q = . ∂ y d − c

(4.7.5)

= u( x, x y), , y), jadi

∂u/∂x

=

x*. ∂ x*/ x*/ ∂ x ∂ũ/ ∂ x*.

=

x*. p p/( /(b b-a) + ∂ũ/ ∂ x*.

= p/( p/(b b-a).∂ũ/ ∂ x* misal H =

+ ∂ũ/ ∂ y*. y*. ∂ y*/ y*/ ∂ y 0 (4.7.6)

∂u/∂ x

88

jadi

∂2u/∂ x2 = ∂ H /∂ x = ∂ H /∂ x*. x*. ∂ x*/ x*/∂ x + ∂ H /∂ y*. y*. ∂ y*/ y*/∂ y = ∂ H /∂ x*. x*. p p/( /(b b-a) +

0

= p/( p/(b-a b-a)) . ∂ H /∂ x* x* = p/( p/(b-a b-a)) . p . p/( /(b-a b-a)) ∂2ũ/∂ x*2 diperoleh 2 2 2 ∂2u/ ∂ x2 = ( p/(b-a)) p/(b-a)) . ∂ ũ/∂ x* x*

(4.7.7)

Dengan cara serupa (cara yang digunakan untuk mendapatkan persamaan (4.7.6) dan (4.7.7)) diperoleh (q/(d /(d -c)) . ∂ũ/∂ y* ∂u/∂ y = (q

(4.7.8)

2 2 2 (q/(d /(d -c)) . ∂ ũ/∂ y* y* ∂2u/∂ y2 = (q

(4.7.9)

Persamaan (4.7.7) dan (4.7.9) disubstitusi ke persamaan (4.7.1) diperoleh persamaan differensial Poisson yang telah ditransformasi, yaitu

2 2 2 γ.∂2 ũ / ∂ x x* * + Ψ.∂ ũ / ∂ y

dengan

γ = ( p /(b − a)) 2 , Ψ

= f ( x*, x* y*) , y*)

(4.7.10)

= (q /( d − c)) 2

89

Mengenai syarat batas, misalkan syarat batas sebelah atas domain sebelum ditransformasi bertipe Robin, yaitu

∂u/∂ y = αu + β.

(4.7.11)

Persamaan (4.7.8) disubstitusi ke persamaan (4.7.11), diperoleh (q/(d -c)) . ∂ũ/ ∂ y* = αu + β atau

∂ũ/∂ y*= (α ũ + β )(d-c)/q

(4.7.12)

Persamaan (4.7.12) merupakan transformasi syarat batas (4.7.11). Transformasi koordinat diatas diilustrasikan pada gambar berikut.

∂u/∂y = αu + β

Y d

Y*

∂ũ/∂ y*= (α ũ + β )(d-c)/q

q

c

a

∂2u/ ∂x2

b +

X

∂2u/ ∂y2 = f(x,y)

0

p

γ.∂2 ũ / ∂ x*2+ Ψ.∂2ũ / ∂ y*2 dengan: γ=(p/(b-a))2, Ψ

Gambar 4.27.

X* = f(x*,y*) (q/(d-c))2

=

Transformasi koordinat

90

BAB V PENUTUP

5.1 Kesimpulan Metode iterasi untuk menyelesaikan sistem persamaan linier konvergen jika dan hanya jika spectral radius matrik iterasinya ( iterative matrix) bernilai kurang dari satu. Tebakan awal iterasi SOR yang diatur dengan memanfaatkan batas domain bertipe Dirichlet dapat membantu meningkatkan akurasi numerik serta memperkecil banyaknya iterasi untuk mencapai kekonvergenan. Sistem persamaan linier yang dihasilkan skema beda hingga PD Poisson dan Laplace 2D dengan syarat batas seluruhnya bertipe Dirichlet , atau hanya terdapat tipe Dirichlet dan Neummann, bersifat dominan diagonal . Apabila pada batas domain terdapat tipe Robin yang tidak termasuk tipe Dirichlet atau Neumann maka sifat dominan diagonal dapat tidak dipenuhi.

Metode iterasi SOR dapat tidak konvergen jika digunakan untuk menyelesaikan sistem persamaan linier yang dihasilkan oleh skema beda hingga PD poisson dan Laplace 2D dengan syarat batas domain terdapat tipe Robin (yang bukan termasuk tipe Dirichlet atau Neumann).

5.2 Masalah lanjutan Penentuan parameter SOR yang menghasilkan banyaknya iterasi paling kecil untuk mencapai kekonvergenan dan kasus solusi numerik

yang tidak

tunggal, menjadi hal menarik untuk dikaji lebih lanjut.

91

LAMPIRAN 1 Mencari spectral

radius matrik

iterasi untuk contoh 4.1 dan contoh 4.2

»

Matrik

I=33 J=33

Ukuran_Matrik_E:

Sistim persamaan linier Au=b Ukuran matrik A : 961 x 961 Parameter SOR w =

1.2000

A = 4 -1 0 0 0 0 0

7 baris dan 7 kolom pertama -1 0 0 0 0 4 -1 0 0 0 -1 4 -1 0 0 0 -1 4 -1 0 0 0 -1 4 -1 0 0 0 -1 4 0 0 0 0 -1

0 0 0 0 0 -1 4

4 -1 0 0 0 0 0

... 7 baris dan 7 kolom terakhir -1 0 0 0 0 4 -1 0 0 0 -1 4 -1 0 0 0 -1 4 -1 0 0 0 -1 4 -1 0 0 0 -1 4 0 0 0 0 -1

0 0 0 0 0 -1 4

Matrik

D

Ukuran_Matrik_D:

961 x

961

D = 4 0 0 0 0 0 0

4 0 0 0 0 0 0

7 baris dan 7 kolom pertama 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0 ... 7 baris dan 7 kolom terakhir 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0

0 0 0 0 0 0 4

E 961 x

961

E = 7 baris dan 7 kolom pertama 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 ... 7 baris dan 7 kolom terakhir 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1

Ukuran matrik

0 0 0 0 0 0 0

0 0 0 0 0 0 0

F: 961 x 961

F = 0 0 0 0 0 0 0 ... 0 0 0 0 0 0 0

7 baris 1 0 0 0 0 0 0

dan 7 kolom 0 0 1 0 0 1 0 0 0 0 0 0 0 0

pertama 0 0 0 0 0 0 1 0 0 1 0 0 0 0

0 0 0 0 0 1 0

7 baris 1 0 0 0 0 0 0

dan 7 kolom 0 0 1 0 0 1 0 0 0 0 0 0 0 0

terakhir 0 0 0 0 0 0 1 0 0 1 0 0 0 0

0 0 0 0 0 1 0

0 0 0 0 0 0 4

92

Ukuran matrik iterasi M : 961X961

M = inv(D-w*E)*((1-w)*D+w*F) M = 7 baris dan 7 kolom pertama -0.2000 0.3000 0 0 -0.0600 -0.1100 0.3000 0 -0.0180 -0.0330 -0.1100 0.3000 -0.0054 -0.0099 -0.0330 -0.1100 -0.0016 -0.0030 -0.0099 -0.0330 -0.0005 -0.0009 -0.0030 -0.0099 -0.0001 -0.0003 -0.0009 -0.0030 ... 7 baris dan 7 kolom terakhir -0.0200 0.3000 0 0 0.0210 -0.0200 0.3000 0 0.0144 0.0210 -0.0200 0.3000 0.0068 0.0144 0.0210 -0.0200 0.0028 0.0068 0.0144 0.0210 0.0010 0 .0028 0. 0068 0.0 144 0.0004 0 .0010 0. 0028 0.0 068

0 0 0 0.3000 -0.1100 -0.0330 -0.0099

0 0 0 0 0.3000 -0.1100 -0.0330

0 0 0 0 0 0.3000 -0.1100

0 0 0 0.3000 -0.0200 0.0 210 0.0 144

0 0 0 0 0.3000 -0.0200 0.02 10

0 0 0 0 0 0.3000 -0.0200

Spectral radius matrik M -----------------------spectral_radius = max(abs(eig(M)) = 0.9856 selesai »

» I=33 J=33 Sistim persamaan linier Au=b Ukuran matrik A : 961 x 961 Parameter SOR w = 1 ... Spectral radius matrik M -----------------------spectral_radius = max(abs(eig(M)) = 0.9912 selesai »

.

93

LAMPIRAN 2 Proses iterasi untuk contoh 4.1 2A) Nilai awal Iterasi

u

(0)

=0

» Contoh 4.1 =========== PD Poisson 2D dengan syarat batas Dirichlet Nilai Awal Iterasi u = 0 Bentuk: Uxx + Uyy = f f =

4

Dimensi Domain Segi Empat -------------------------lebar = ?1 tinggi = ?1 -------------------------I:indeks maksimum searah sb-x. i=1,2,...,I. J:indeks maksimum searah sb-y. j=1,2,...,J. Parameter SOR : w = 1.2000 Toleransi kekonvergenan : epsilon =

--------------Iterasi ke 1 Iterasi ke 2 Iterasi ke 3 Iterasi ke 4 Iterasi ke 5 Iterasi ke 6 Iterasi ke 7 Iterasi ke 8 Iterasi ke 9 Iterasi ke 10 Iterasi ke 11 Iterasi ke 12 Iterasi ke 13 Iterasi ke 14 Iterasi ke 15 Iterasi ke 16 Iterasi ke 17 Iterasi ke 18 Iterasi ke 19 Iterasi ke 20 Iterasi ke 21 Iterasi ke 22 Iterasi ke 23 Iterasi ke 24 ...

Proses banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak banyak

I = ?33 J = ?33

0.0100

iterasi SOR ---------------U konvergen 886 // 1089 U konvergen 825 // 1089 U konvergen 786 // 1089 U konvergen 751 // 1089 U konvergen 720 // 1089 U konvergen 731 // 1089 U konvergen 709 // 1089 U konvergen 721 // 1089 U konvergen 721 // 1089 U konvergen 700 // 1089 U konvergen 710 // 1089 U konvergen 734 // 1089 U konvergen 735 // 1089 U konvergen 735 // 1089 U konvergen 727 // 1089 U konvergen 738 // 1089 U konvergen 753 // 1089 U konvergen 767 // 1089 U konvergen 779 // 1089 U konvergen 787 // 1089 U konvergen 797 // 1089 U konvergen 808 // 1089 U konvergen 823 // 1089 U konvergen 837 // 1089

Iterasi ke 48 banyak u konvergen 1070 // 1089 Iterasi ke 49 banyak u konvergen 1076 // 1089 Iterasi ke 50 banyak u konvergen 1085 // 1089 Iterasi ke 51 banyak u konvergen 1089 // 1089 --------------------------------------------------Hasil dalam bentuk tabel dan grafik

94

TABEL PERBEDAAN SOLUSI NUMERIK DENGAN SOLUSI EKSAK --------------------------------------------------------------------------no

x

y

Solusi numerik Solusi eksak Error global u_numerik u_eksak e = u_numerik - u_eksak --------------------------------------------------------------------------| 1 | 0.00000 | 0.00000 | 0. 00000000 | 0. 00000000 | 0 .00000000 | | 2 | 0.00000 | 0.03125 | 0. 00097656 | 0. 00097656 | 0 .00000000 | | 3 | 0.00000 | 0.06250 | 0. 00390625 | 0. 00390625 | 0 .00000000 | | 4 | 0.00000 | 0.09375 | 0. 00878906 | 0. 00878906 | 0 .00000000 | | 5 | 0.00000 | 0.12500 | 0. 01562500 | 0. 01562500 | 0 .00000000 | | 6 | 0.00000 | 0.15625 | 0. 02441406 | 0. 02441406 | 0 .00000000 | | 7 | 0.00000 | 0.18750 | 0. 03515625 | 0. 03515625 | 0 .00000000 | | 8 | 0.00000 | 0.21875 | 0. 04785156 | 0. 04785156 | 0 .00000000 | | 9 | 0.00000 | 0.25000 | 0. 06250000 | 0. 06250000 | 0 .00000000 | | 10 | 0.00000 | 0.28125 | 0.07910156 | 0.07910156 | 0.00000000 | | 11 | 0.00000 | 0.31250 | 0.09765625 | 0.09765625 | 0.00000000 | | 12 | 0.00000 | 0.34375 | 0.11816406 | 0.11816406 | 0.00000000 | | 13 | 0.00000 | 0.37500 | 0.14062500 | 0.14062500 | 0.00000000 | | 14 | 0.00000 | 0.40625 | 0.16503906 | 0.16503906 | 0.00000000 | | 15 | 0.00000 | 0.43750 | 0.19140625 | 0.19140625 | 0.00000000 | | 16 | 0.00000 | 0.46875 | 0.21972656 | 0.21972656 | 0.00000000 | | 17 | 0.00000 | 0.50000 | 0.25000000 | 0.25000000 | 0.00000000 | | 18 | 0.00000 | 0.53125 | 0.28222656 | 0.28222656 | 0.00000000 | | 19 | 0.00000 | 0.56250 | 0.31640625 | 0.31640625 | 0.00000000 | ..... | 640 | 0.59375 | 0.37500 | 0.06989351 | 0.49316406 | -0.42327055 | | 641 | 0.59375 | 0.40625 | 0.06913434 | 0.51757813 | -0.44844378 | | 642 | 0.59375 | 0.43750 | 0.07290868 | 0.54394531 | -0.47103663 | | 643 | 0.59375 | 0.46875 | 0.08170180 | 0.57226563 | -0.49056383 | | 644 | 0.59375 | 0.50000 | 0.09607642 | 0.60253906 | -0.50646264 | | 645 | 0.59375 | 0.53125 | 0.11665160 | 0.63476563 | -0.51811403 | | 646 | 0.59375 | 0.56250 | 0.14407205 | 0.66894531 | -0.52487326 | 0.59375 | 0.59375 | 0.17896975 | 0.70507813 | -0.52610838 | | 647 | | 648 | 0.59375 | 0.62500 | 0.22192028 | 0.74316406 | -0.52124378 | | 649 | 0.59375 | 0.65625 | 0.27339754 | 0.78320313 | -0.50980559 | | 650 | 0.59375 | 0.68750 | 0.33373044 | 0.82519531 | -0.49146487 | | 651 | 0.59375 | 0.71875 | 0.40306570 | 0.86914063 | -0.46607492 | | 652 | 0.59375 | 0.75000 | 0.48134004 | 0.91503906 | -0.43369902 | | 653 | 0.59375 | 0.78125 | 0.56826454 | 0.96289063 | -0.39462608 | | 654 | 0.59375 | 0.81250 | 0.66332282 | 1.01269531 | -0.34937250 | | 655 | 0.59375 | 0.84375 | 0.76578337 | 1.06445313 | -0.29866975 | ..... | 1084 | 1.00000 | 0.84375 | 1.71191406 | 1.71191406 | 0.00000000 | | 1085 | 1.00000 | 0.87500 | 1.76562500 | 1.76562500 | 0.00000000 | | 1086 | 1.00000 | 0.90625 | 1.82128906 | 1.82128906 | 0.00000000 | | 1087 | 1.00000 | 0.93750 | 1.87890625 | 1.87890625 | -0.00000000 | | 1088 | 1.00000 | 0.96875 | 1.93847656 | 1.93847656 | 0.00000000 | | 1089 | 1.00000 | 1.00000 | 2.00000000 | 2.00000000 | 0.00000000 | ---------------------------------------------------------------------------

Nilai absolut terbesar dari error global : maks|e(x,y)|= terjadi pada koordinat x = 0.59375 y = 0.59375 dengan nilai error global e(x,y)= -0.52610838

0.52610838

Selesai »

95

2B) Nilai Awal Iterasi

u

( 0)

dimodifikasi

» PD Poisson 2D dengan syarat batas Dirichlet Nilai Awal Iterasi u di modifikasi Bentuk: uxx + uyy = f f = 4 Domain segi empat berukuran : ---------------------------lebar = 1 tinggi = 1 Level grid. L =5 ------------------------Banyak Titik Grid : 33 X 33 =1089 w = 1.2000 Toleransi kekonvergenan : epsilon = 0.0100

----------------Proses iterasi SOR---------------------Iterasi ke 1 banyak U konvergen 1089 // 1089 -------------------------------------------------------Hasil dalam bentuk grafik -------------------------------------------------------Error absolut terbesar: maks|e(x,y)|= 2.2204e-016 terjadi pada koordinat x= 0.90625 y= 0.96875 dengan nilai e(x,y)= 2.2204e-016

2C) Nilai Awal Iterasi

u

( 0)

dimodifikasi dan

ω

=

1

» PD Poisson 2D dengan syarat batas Dirichlet Nilai Awal Iterasi u di modifikasi Bentuk: uxx + uyy = f f = 4 Domain segi empat berukuran : ---------------------------lebar = 1 tinggi = 1 Level grid. L =5 ------------------------Banyak Titik Grid : 33 X 33 =1089 w = 1.0000 Toleransi kekonvergenan : epsilon = 0.0100

----------------Proses iterasi SOR---------------------Iterasi ke 1 banyak U konvergen 1089 // 1089 -------------------------------------------------------Hasil dalam bentuk grafik -------------------------------------------------------Error absolut terbesar: maks|e(x,y)|= 0 terjadi pada koordinat x= 0 y= 0 dengan nilai e(x,y)= 0

96

2D) Nilai Awal Iterasi

u

( 0)

dimodifikasi, grid diperbanyak dan

ω

=

1

PD Poisson 2D dengan syarat batas Dirichlet Nilai Awal Iterasi u di modifikasi Bentuk: Uxx + Uyy = f f = 4 Domain segi empat berukuran : ---------------------------lebar = 1 tinggi = 1 Level grid. L =11 ------------------------Banyak Titik Grid : 2049 X 2049 = 4198401 w = 1.0000 Toleransi kekonvergenan : epsilon = 0.0100

-------------------Proses iterasi SOR-----------------------Iterasi ke 1 banyak U konvergen 4198401 //4198401 ------------------------------------------------------------Hasil dalam bentuk grafik ------------------------------------------------------------Error absolut terbesar: maks|e(x,y)|= 0 terjadi pada koordinat x= 0 y= 0 dengan nilai e(x,y)= 0 Selesai »

97

LAMPIRAN 3 Proses iterasi untuk contoh 4.2 3A) Nilai Awal Iterasi

u

( 0)

=0

» Contoh 4.2 =========== PD Laplace 2D dengan syarat batas Dirichlet Nilai Awal Iterasi u = 0 Bentuk: Uxx + Uyy = f f = 4 Ukuran Domain : -------------lebar = 1 tinggi = 1 ---------------I:indeks maksimum pilahan pada arah sb-x. i=1,2,...,I. J:indeks maksimum pilahan pada arah sb-y. j=1,2,...,J. epsilon = w =

I=33 J=33

0.0100

1.2000

-----------Proses iterasi SOR-------------------Iterasi ke 1 banyak U interior yang konvergen 650 // 961 Iterasi ke 2 banyak U interior yang konvergen 566 // 961 Iterasi ke 3 banyak U interior yang konvergen 504 // 961 Iterasi ke 4 banyak U interior yang konvergen 458 // 961 Iterasi ke 5 banyak U interior yang konvergen 425 // 961 ... Iterasi ke 74 banyak U interior yang konvergen 823 // 961 Iterasi ke 75 banyak U interior yang konvergen 827 // 961 Iterasi ke 76 banyak U interior yang konvergen 837 // 961 Iterasi ke 77 banyak U interior yang konvergen 847 // 961 Iterasi ke 78 banyak U interior yang konvergen 854 // 961 Iterasi ke 79 banyak U interior yang konvergen 862 // 961 Iterasi ke 80 banyak U interior yang konvergen 871 // 961 Iterasi ke 81 banyak U interior yang konvergen 882 // 961 Iterasi ke 82 banyak U interior yang konvergen 888 // 961 Iterasi ke 83 banyak U interior yang konvergen 895 // 961 Iterasi ke 84 banyak U interior yang konvergen 906 // 961 Iterasi ke 85 banyak U interior yang konvergen 915 // 961 Iterasi ke 86 banyak U interior yang konvergen 922 // 961 Iterasi ke 87 banyak U interior yang konvergen 933 // 961 Iterasi ke 88 banyak U interior yang konvergen 941 // 961 Iterasi ke 89 banyak U interior yang konvergen 950 // 961 Iterasi ke 90 banyak U interior yang konvergen 960 // 961 Iterasi ke 91 banyak U interior yang konvergen 961 // 961 ---------------------------------------------------------------Hasil dalam bentuk tabel dan grafik

98

TABEL PERBEDAAN SOLUSI NUMERIK DENGAN SOLUSI EKSAK --------------------------------------------------------------------------no

x

y

Solusi numerik Solusi eksak Error global u_numerik u_eksak e = u_numerik - u_eksak --------------------------------------------------------------------------| 1 | 0.00000 | 0.00000 | 0. 00000000 | 0. 00000000 | 0 .00000000 | | 2 | 0.00000 | 0.03125 | -0.06347656 | -0.06347656 | 0.00000000 | | 3 | 0.00000 | 0.06250 | -0.12890625 | -0.12890625 | 0.00000000 | | 4 | 0.00000 | 0.09375 | -0.19628906 | -0.19628906 | 0.00000000 | ... | 511 | 0.46875 | 0.46875 | -0.54168296 | -1.12792969 | 0.58624673 | | 512 | 0.46875 | 0.50000 | -0.66334940 | -1.26464844 | 0.60129904 | | 513 | 0.46875 | 0.53125 | -0.79315660 | -1.40332031 | 0.61016371 | | 514 | 0.46875 | 0.56250 | -0.93156319 | -1.54394531 | 0.61238212 | | 515 | 0.46875 | 0.59375 | -1.07892181 | -1.68652344 | 0.60760163 | | 516 | 0.46875 | 0.62500 | -1.23546050 | -1.83105469 | 0.59559419 | | 517 | 0.46875 | 0.65625 | -1.40126783 | -1.97753906 | 0.57627123 | | 518 | 0.46875 | 0.68750 | -1.57628268 | -2.12597656 | 0.54969389 | | 519 | 0.46875 | 0.71875 | -1.76028902 | -2.27636719 | 0.51607817 | | 520 | 0.46875 | 0.75000 | -1.95291623 | -2.42871094 | 0.47579471 | | 521 | 0.46875 | 0.78125 | -2.15364483 | -2.58300781 | 0.42936298 | | 522 | 0.46875 | 0.81250 | -2.36181750 | -2.73925781 | 0.37744031 | | 523 | 0.46875 | 0.84375 | -2.57665498 | -2.89746094 | 0.32080596 | | 524 | 0.46875 | 0.87500 | -2.79727608 | -3.05761719 | 0.26034110 | | 525 | 0.46875 | 0.90625 | -3.02272111 | -3.21972656 | 0.19700545 | | 526 | 0.46875 | 0.93750 | -3.25197759 | -3.38378906 | 0.13181147 | | 527 | 0.46875 | 0.96875 | -3.48400728 | -3.54980469 | 0.06579740 | | 528 | 0.46875 | 1.00000 | -3.71777344 | -3.71777344 | 0.00000000 | | 529 | 0.50000 | 0.00000 | 0.75000000 | 0.75000000 | 0.00000000 | | 530 | 0.50000 | 0.03125 | 0.68530529 | 0.63964844 | 0.04565685 | | 1079 | 1.00000 | 0.68750 | -1.91015625 | -1.91015625 | 0.00000000 | | 1080 | 1.00000 | 0.71875 | -2.11035156 | -2.11035156 | 0.00000000 | | 1081 | 1.00000 | 0.75000 | -2.31250000 | -2.31250000 | 0.00000000 | | 1082 | 1.00000 | 0.78125 | -2.51660156 | -2.51660156 | 0.00000000 | | 1083 | 1.00000 | 0.81250 | -2.72265625 | -2.72265625 | 0.00000000 | | 1084 | 1.00000 | 0.84375 | -2.93066406 | -2.93066406 | 0.00000000 | | 1085 | 1.00000 | 0.87500 | -3.14062500 | -3.14062500 | 0.00000000 | | 1086 | 1.00000 | 0.90625 | -3.35253906 | -3.35253906 | 0.00000000 | | 1087 | 1.00000 | 0.93750 | -3.56640625 | -3.56640625 | 0.00000000 | | 1088 | 1.00000 | 0.96875 | -3.78222656 | -3.78222656 | 0.00000000 | | 1089 | 1.00000 | 1.00000 | -4.00000000 | -4.00000000 | 0.00000000 | --------------------------------------------------------------------------Error absolut terbesar : maks|e(x,y)|= terjadi pada koordinat x= 0.46875 y= 0.56250 dengan nilai e(x,y)= 0.61238212

0.61238212

Selesai »

99

3B) Nilai Awal Iterasi

u

(0)

dimodifikasi

» PD Laplace 2D dengan syarat batas Dirichlet Nilai Awal Iterasi u dimodifikasi Domain segi empat berukuran : ---------------------------lebar = ? 1 tinggi = ? 1 level grid 1,2,...L. L = ?5 ---------------------------Banyak Titik Grid : 33 X 33 = 1089 epsilon = 0.0100 w =1.2000 ---------------------------Proses iterasi SOR----------------------Iterasi ke 1 Banyak titik grid yang konvergen : 1089 // 1089 -------------------------------------------------------------------Hasil dalam bentuk grafik --------------------------Error absolut terbesar: maks|e(x,y)|= 4.4409e-016 terjadi pada koordinat x= 0.37500 y= 0.96875 dengan nilai e(x,y)= 4.4409e-016

3C) Nilai Awal Iterasi

u

( 0)

dimodifikasi dan

ω

=

1

» PD Laplace 2D dengan syarat batas Dirichlet Nilai Awal Iterasi u dimodifikasi Domain segi empat berukuran : ---------------------------lebar = ? 1 tinggi = ? 1 level grid 1,2,...L. L = ?5 ---------------------------Banyak Titik Grid : 33 X 33 = 1089 epsilon = 0.0100 w =1.0000 ---------------------------Proses iterasi SOR----------------------Iterasi ke 1 Banyak titik grid yang konvergen : 1089 // 1089 -------------------------------------------------------------------Hasil dalam bentuk grafik --------------------------Error absolut terbesar: maks|e(x,y)|= 0 terjadi pada koordinat x= 0 y= 0 dengan nilai e(x,y)= 0

100

3D) Nilai Awal Iterasi

u

( 0)

dimodifikasi, grid diperbanyak dan

ω

=

1

» PD Laplace 2D dengan syarat batas Dirichlet Nilai Awal Iterasi u dimodifikasi Domain segi empat berukuran : ---------------------------lebar = ? 1 tinggi = ? 1 level grid 1,2,...L. L = ?11 ---------------------------Banyak Titik Grid : 2049 X 2049 =4198401 epsilon = 0.0100 w = 1.0000 -----------------------------Proses iterasi SOR-------------------------------Iterasi ke 1 Banyak titik grid yang konvergen : 4198401 // 4198401 ------------------------------------------------------------------------------Hasil dalam bentuk grafik -------------------------Error absolut terbesar : maks|e(x,y)|= 0 terjadi pada koordinat x= 0 y= 0 dengan nilai e(x,y)= 0 Selesai »

101

LAMPIRAN 4 Contoh 4.3 4A) Mencari spectral radius matrik iterasi » Ukuran Matrik A,D,E,F,dan M masing-masing :

1024

X

1024

A = 33 baris dan 33 kolom pertama Columns 1 through 18

4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Columns 19 through 33 0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

0 0 0 0 0 0 0 0 0

-1 0 0 0 0 0 0 0 0

102

0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 -2 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ...

0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4

33 baris dan 33 kolom terakhir Columns 1 through 18 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

103

Columns 19 through 33 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 -2 4 Matrik D,E,F ... Parameter SOR : w=1 Matrik iterasi M=inv(D-w*E)*((1-w)*D+w*F) M = 0 0 0 0 0 0 0 ...

0.2500 0.0625 0.0156 0.0039 0.0010 0.0002 0.0001

0.1875 0.0781 0.0273 0.0088 0.0027 0.0008 0.0005

7 baris 0 0.2500 0.0625 0.0156 0.0039 0.0010 0.0002

0.2500 0.1875 0.0781 0.0273 0.0088 0.0027 0.0016

dan 7 kolom pertama 0 0 0 0 0 0 0.2500 0 0 0.0625 0.2500 0 0.0156 0.0625 0.2500 0.0039 0.0156 0.0625 0.0010 0.0039 0.0156

7 baris dan 7 kolom terakhir 0 0 0 0.2500 0 0 0.1875 0.2500 0 0.0781 0.1875 0.2500 0.0273 0.0781 0.1875 0.0088 0.0273 0.0781 0.0054 0.0176 0.0547

0 0 0 0 0 0.2500 0.0625

0 0 0 0 0.2500 0.1875 0.1563

0 0 0 0 0 0.2500 0.2500

Spectral radiusmatrik M = max(abs(eig(M)) spectral_radius = 0.9976

104

4B) Proses iterasi PD Poisson 2D dengan dua tipe syarat batas yaitu Dirichlet dan Neumann Nilai Awal Iterasi u di modifikasi Bentuk: uxx + uyy = f f = 4 Domain segi empat berukuran : ---------------------------lebar = 1 tinggi = 1 Level grid. L =5 ------------------------Banyak Titik Grid : 33 X 33 =1089 w = 1.0000 Toleransi kekonvergenan : epsilon = 0.000001

----------------Proses iterasi SOR---------------------Iterasi ke 1 banyak U konvergen 1089 // 1089 -------------------------------------------------------Hasil dalam bentuk grafik -------------------------------------------------------Error absolut terbesar: maks|e(x,y)|= 0 terjadi pada koordinat x= 0 y= 0 dengan nilai e(x,y)= 0 Selesai »

.

105

LAMPIRAN 5 Contoh 4.4 5A) Mencari spectral radius matrik iterasi SOR dengan Ukuran Matrik A,D,E,F dan M masing-masing : 25 x 25 A= Columns 1 through 11 4 -2 0 0 0 -1 0 0 0 0 0 -1 4 -1 0 0 0 -1 0 0 0 0 0 -1 4 -1 0 0 0 -1 0 0 0 0 0 -1 4 -1 0 0 0 -1 0 0 0 0 0 -2 4 0 0 0 0 -1 0 -1 0 0 0 0 4 -2 0 0 0 -1 0 -1 0 0 0 -1 4 -1 0 0 0 0 0 -1 0 0 0 -1 4 -1 0 0 0 0 0 -1 0 0 0 -1 4 -1 0 0 0 0 0 -1 0 0 0 -2 4 0 0 0 0 0 0 -1 0 0 0 0 4 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Columns 12 through 22 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 -2 0 0 0 -1 0 0 0 0 0 0 4 -1 0 0 0 -1 0 0 0 0 0 -1 4 -1 0 0 0 -1 0 0 0 0 0 -1 4 -1 0 0 0 -1 0 0 0 0 0 -2 4 0 0 0 0 -1 0 0 0 0 0 0 4 -2 0 0 0 -1 0 -1 0 0 0 -1 4 -1 0 0 0 -1 0 -1 0 0 0 -1 4 -1 0 0 0 0 0 -1 0 0 0 -1 4 -1 0 0 0 0 0 -1 0 0 0 -2 4 0 0 0 0 0 0 -1 0 0 0 0 4 -2 0 0 0 0 0 -1 0 0 0 -1 4 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0

ω

=

1.5 dan

ω

=

1.333

Columns 23 through 25 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 -1 0 0 0 -1 0 0 4 -1 0 -1 4 -1 0 -2 4

Matrik D,E,F ...

106

Parameter SOR : w =1.5 -----------------------------Matrik iterasi : M=inv(D-w*E) *((1-w)*D+w *F) M = -0.5000 -0.1875 -0.0703 -0.0264 -0.0198 -0.1875 -0.1406

0.7500 -0.2188 -0.0820 -0.0308 -0.0231 0.2813 0.0234

7 baris dan 7 kolom pertama 0 0 0 0.3750 0 0 -0.3594 0.3750 0 -0.1348 -0.3594 0.3750 -0.1011 -0.2695 -0.2188 0 0 0 0.1406 0 0

0.3750 0.1406 0.0527 0.0198 0.0148 -0.3594 -0.0820

0 0.3750 0.1406 0.0527 0.0396 0.7500 -0.0781

0.3750 -0.0781 0 0 0 0.1406 0.0762

7 baris dan 7 kolom terakhir 0.0198 0.0527 0.1406 0.0148 0 .0396 0 .1055 -0.3594 0.7500 0 -0.0820 -0.0781 0.3750 -0.0110 0.0234 -0.2188 0.0033 0.0286 -0.0293 0.0080 0.0363 0.0176

0.3750 0. 2813 0 0 0.3750 -0.2188 -0.0586

0 0.3 750 0 0 0 0.3750 -0.0781

... -0.2188 -0.0586 0 0 0.1406 -0.0293 -0.0439

----------------------------------------------------------------------Spectral_radius_matrik_M = max(abs(eig(M)) = 0.5000

Parameter SOR : w=1.333 ---------------------Matrik iterasi : M=inv(D-w*E)*((1-w)*D+w*F) M = -0.3330 -0.1110 -0.0370 -0.0123 -0.0082 -0.1110 -0.0740

0.6665 -0.1109 -0.0370 -0.0123 -0.0082 0.2221 0.0371

7 baris dan 7 kolom pertama 0 0 0 0.3332 0 0 -0.2219 0.3332 0 -0.0740 -0.2219 0.3332 -0.0493 -0.1479 -0.1109 0 0 0 0.1111 0 0

0.3332 0.1111 0.0370 0.0123 0.0082 -0.2219 -0.0370

0 0.3332 0.1111 0.0370 0.0247 0.6665 0.0002

... 7 baris dan 7 kolom terakhir -0.1109 0.3332 0.0123 0.0370 0.1111 0.3332 0 0.0001 0.0002 0.0082 0.0247 0.0740 0.2221 0.3332 0 0 -0.2219 0.6665 0 0 0 0 0 -0.0370 0.0002 0.3332 0 0 0.1111 0 0.0000 0.0371 -0.1109 0.3332 0 0.0001 0.1111 0.0041 0 .0247 0 .0001 -0.1109 0.3332 0.0001 0.0741 0.0055 0.0247 0.0247 0.0001 0.0002 ----------------------------------------------------------------------Spectral_radius_ matrik_M = max(abs(eig(M)) = 0.7286 >>

.

107

5B) Proses iterasi dengan Nilai Awal Iterasi u » BENTUK UMUM : Uxx + Uyy = f ------------------------------Jenis Persamaan Differensial: 1. Laplace 2. Poisson * pilih 1-2 : ? 1 UKURAN DOMAIN -------------lebar = ? 1 tinggi = ? 1 I:indeks maks.pilahan pd arah sb-x. i=1,2,.,I. J:indeks maks.pilahan pd arah sb-y. j=1,2,.,J.

SYARAT BATAS -----------Syarat batas sisi ke 1 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dy = 1 Syarat batas sisi ke 2 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dy = 1 Syarat batas sisi ke 3 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dx = 1

( 0)

=

0 dan

ω

=

1.5

I = ? 5 J = ? 5

Syarat batas sisi ke 4 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dx = 1 -------------------------------------Ganti atribut iterasi (y/t) = ? y -------------------------------------parameter SOR: w = ? 1.5 tolerans : epsilon =?0.01 Nilai Awal Iterasi u = ?0 ------------Proses iterasi SOR------------------Iterasi ke 1 banyak U konvergen 7 // 25 Iterasi ke 2 banyak U konvergen 14 // 25 Iterasi ke 3 banyak U konvergen 15 // 25 Iterasi ke 4 banyak U konvergen 19 // 25 Iterasi ke 5 banyak U konvergen 22 // 25 Iterasi ke 6 banyak U konvergen 25 // 25 ------------------------------------------------Hasil dalam bentuk grafik Selesai »

5C) Proses iterasi dengan Nilai Awal Iterasi u » BENTUK UMUM : Uxx + Uyy = f ------------------------------Jenis Persamaan Differensial: 1. Laplace 2. Poisson * pilih 1-2 : ? 1 UKURAN DOMAIN -------------lebar = ? 1 tinggi = ? 1 I:indeks maks.pilahan pd arah sb-x. i=1,2,.,I. J:indeks maks.pilahan pd arah sb-y. j=1,2,.,J.

( 0)

=

(0)

10 (u i , j

=

10) dan

ω

=

1. 5

I = ? 5 J = ? 5

108

SYARAT BATAS -----------Syarat batas sisi ke 1 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dy = 1 Syarat batas sisi ke 2 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dy = 1 Syarat batas sisi ke 3 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dx = 1

Syarat batas sisi ke 4 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dx = 1 -------------------------------------Ganti atribut iterasi (y/t) = ? y -------------------------------------parameter SOR: w = ? 1.333 tolerans : epsilon =?0.01 Nilai Awal Iterasi u = ?10

-----------Proses iterasi SOR-------------------Iterasi ke 1 banyak U konvergen 8 // 25 Iterasi ke 2 banyak U konvergen 14 // 25 Iterasi ke 3 banyak U konvergen 18 // 25 Iterasi ke 4 banyak U konvergen 21 // 25 Iterasi ke 5 banyak U konvergen 23 // 25 Iterasi ke 6 banyak U konvergen 25 // 25 ------------------------------------------------Hasil dalam bentuk grafik Selesai »

.

.

109

LAMPIRAN 6

Mencari spectral radius matrik iterasi untuk contoh 4.5 dan contoh 4.6

6A) S pectral radius matrik iterasi untuk contoh 4.5

>> alfa3 = 3 alfa4 = -4 Ukuran Matrik A,D,E,F, dan M masing-masing :

1056 X 1056

A = 4.1875 -1.0000 0 0 0 0 0

-2.0000 4.0000 -1.0000 0 0 0 0

4.0000 -1.0000 0 0 0 0 0

-1.0000 4.0000 -1.0000 0 0 0 0

7 baris dan 7 kolom pertama 0 0 0 -1.0000 0 0 4.0000 -1.0000 0 -1.0000 4.0000 -1.0000 0 -1.0000 4.0000 0 0 -1.0000 0 0 0 ... 7 baris dan 7 kolom terakhir 0 0 0 -1.0000 0 0 4.0000 -1.0000 0 -1.0000 4.0000 -1.0000 0 -1.0000 4.0000 0 0 -1.0000 0 0 0

0 0 0 0 -1.0000 4.0000 -1.0000

0 0 0 0 0 -1.0000 4.0000

0 0 0 0 -1.0000 4.0000 -2.0000

0 0 0 0 0 -1.0000 4.2500

Matrik D,E,F ... Parameter SOR : w=1.2 Matrik iterasi M=inv(D-w*E)*((1-w)*D+w*F) M = -0.2000 -0.0600 -0.0180 -0.0054 -0.0016 -0.0005 -0.0001

0.5731 -0.0281 -0.0084 -0.0025 -0.0008 -0.0002 -0.0001

7 baris dan 7 kolom pertama 0 0 0 0.3000 0 0 -0.1100 0.3000 0 -0.0330 -0.1100 0.3000 -0.0099 -0.0330 -0.1100 -0.0030 -0.0099 -0.0330 -0.0009 -0.0030 -0.0099

0 0 0 0 0.3000 -0.1100 -0.0330

0 0 0 0 0 0.3000 -0.1100

0.3000 0.0700 0.0750 0.0387 0.0165 0.0064 0.0044


0 0 0 0 0.3000 0.0700 0.1352

0 0 0 0 0 0.3000 0.1289

... 0.0700 0.0750 0.0387 0.0165 0.0064 0.0024 0.0016

110

KOMPONEN MATRIK ITERASI M -----------------------------Baris Kolom no i j M(i,j) -------------------------------| 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12 | 13 | 14 | 15 | 16 | 17 | 18 | 19 | 20 | 21 | 22 | 23 | 24 | 25 | 26 | 27 | 28 | 29 | 30 | 31 | 32 | 33 | 34 | 35 | 36 | 37 | 38 | 39 | 40 | 41 | 42 | 43 | 44 | 45 | 46 | 47 | 48 | 49 ...

| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | |

1 1 1 1 1 1 1 2 2 2 2 2 2 2 3 3 3 3 3 3 3 4 4 4 4 4 4 4 5 5 5 5 5 5 5 6 6 6 6 6 6 6 7 7 7 7 7 7 7

| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | |

1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7

|-0.20000 | 0.57313 | 0.00000 | 0.00000 | 0.00000 | 0.00000 | 0.00000 |-0.06000 |-0.02806 | 0.30000 | 0.00000 | 0.00000 | 0.00000 | 0.00000 |-0.01800 |-0.00842 |-0.11000 | 0.30000 | 0.00000 | 0.00000 | 0.00000 |-0.00540 |-0.00253 |-0.03300 |-0.11000 | 0.30000 | 0.00000 | 0.00000 |-0.00162 |-0.00076 |-0.00990 |-0.03300 |-0.11000 | 0.30000 | 0.00000 |-0.00049 |-0.00023 |-0.00297 |-0.00990 |-0.03300 |-0.11000 | 0.30000 |-0.00015 |-0.00007 |-0.00089 |-0.00297 |-0.00990 |-0.03300 |-0.11000

| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | |

|1115086 | 1050 | 1054 | 0.00000 | |1115087 | 1050 | 1055 | 0.00000 | |1115088 | 1050 | 1056 | 0.00000 | |1115089 | 1051 | 1049 | 0.03870 | |1115090 | 1051 | 1050 | 0.07500 | |1115091 | 1051 | 1051 | 0.07000 | |1115092 | 1051 | 1052 | 0.30000 | |1115093 | 1051 | 1053 | 0.00000 | |1115094 | 1051 | 1054 | 0.00000 | |1115095 | 1051 | 1055 | 0.00000 | |1115096 | 1051 | 1056 | 0.00000 | |1115097 | 1052 | 1049 | 0.01647 | |1115098 | 1052 | 1050 | 0.03870 | |1115099 | 1052 | 1051 | 0.07500 | |1115100 | 1052 | 1052 | 0.07000 | |1115101 | 1052 | 1053 | 0.30000 | |1115102 | 1052 | 1054 | 0.00000 | |1115103 | 1052 | 1055 | 0.00000 | |1115104 | 1052 | 1056 | 0.00000 | |1115105 | 1053 | 1049 | 0.00640 | |1115106 | 1053 | 1050 | 0.01647 | |1115107 | 1053 | 1051 | 0.03870 | |1115108 | 1053 | 1052 | 0.07500 | |1115109 | 1053 | 1053 | 0.07000 | |1115110 | 1053 | 1054 | 0.30000 | |1115111 | 1053 | 1055 | 0.00000 | |1115112 | 1053 | 1056 | 0.00000 | |1115113 | 1054 | 1049 | 0.00236 | |1115114 | 1054 | 1050 | 0.00640 | |1115115 | 1054 | 1051 | 0.01647 | |1115116 | 1054 | 1052 | 0.03870 | |1115117 | 1054 | 1053 | 0.07500 | |1115118 | 1054 | 1054 | 0.07000 | |1115119 | 1054 | 1055 | 0.30000 | |1115120 | 1054 | 1056 | 0.00000 | |1115121 | 1055 | 1049 | 0.00084 | |1115122 | 1055 | 1050 | 0.00236 | |1115123 | 1055 | 1051 | 0.00640 | |1115124 | 1055 | 1052 | 0.01647 | |1115125 | 1055 | 1053 | 0.03870 | |1115126 | 1055 | 1054 | 0.07500 | |1115127 | 1055 | 1055 | 0.07000 | |1115128 | 1055 | 1056 | 0.30000 | |1115129 | 1056 | 1049 | 0.00054 | |1115130 | 1056 | 1050 | 0.00156 | |1115131 | 1056 | 1051 | 0.00439 | |1115132 | 1056 | 1052 | 0.01188 | |1115133 | 1056 | 1053 | 0.03046 | |1115134 | 1056 | 1054 | 0.07105 | |1115135 | 1056 | 1055 | 0.13520 | |1115136 | 1056 | 1056 | 0.12886 | ----------------------------------Spectral radius matrik M = max(abs(eig(M)) spectral_radius = 0.9951

111

6B) S pectral radius matrik iterasi untuk contoh 4.6 >> alfa3 = -3 alfa4 = 4 Ukuran Matrik A,D,E,F dan M masing-masing:

1056X1056

A = 3.8125 -1.0000 0 0 0 0 0

-2.0000 4.0000 -1.0000 0 0 0 0

7 baris dan 7 kolom pertama 0 0 0 -1.0000 0 0 4.0000 -1.0000 0 -1.0000 4.0000 -1.0000 0 -1.0000 4.0000 0 0 -1.0000 0 0 0

0 0 0 0 -1.0000 4.0000 -1.0000

0 0 0 0 0 -1.0000 4.0000

0 0 0 0 -1.0000 4.0000 -2.0000

0 0 0 0 0 -1.0000 3.7500

...

4.0000 -1.0000 0 0 0 0 0

-1.0000 4.0000 -1.0000 0 0 0 0

7 baris dan 7 kolom terakhir 0 0 0 -1.0000 0 0 4.0000 -1.0000 0 -1.0000 4.0000 -1.0000 0 -1.0000 4.0000 0 0 -1.0000 0 0 0

Matrik D,E,F ... Parameter SOR : w=1.2 Matrik iterasi M=inv(D-w*E)*((1-w)*D+w*F) M = -0.2000 -0.0600 -0.0180 -0.0054 -0.0016 -0.0005 -0.0001

0.6295 -0.0111 -0.0033 -0.0010 -0.0003 -0.0001 -0.0000

7 baris dan 7 kolom pertama 0 0 0 0.3000 0 0 -0.1100 0.3000 0 -0.0330 -0.1100 0.3000 -0.0099 -0.0330 -0.1100 -0.0030 -0.0099 -0.0330 -0.0009 -0.0030 -0.0099

0 0 0 0 0.3000 -0.1100 -0.0330

0 0 0 0 0 0.3000 -0.1100

0.3000 0.0700 0.0750 0.0387 0.0165 0.0064 0.0051


0 0 0 0 0.3000 0.0700 0.1677

0 0 0 0 0 0.3000 0.1968

... 0.0700 0.0750 0.0387 0.0165 0.0064 0.0024 0.0018

112

KOMPONEN MATRIK ITERASI M -----------------------------Baris Kolom no i j M(i,j) -------------------------------| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | ...

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49

| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | |

1 1 1 1 1 1 1 2 2 2 2 2 2 2 3 3 3 3 3 3 3 4 4 4 4 4 4 4 5 5 5 5 5 5 5 6 6 6 6 6 6 6 7 7 7 7 7 7 7

| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | |

1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7 1 2 3 4 5 6 7

|-0.20000 | 0.62951 | 0.00000 | 0.00000 | 0.00000 | 0.00000 | 0.00000 |-0.06000 |-0.01115 | 0.30000 | 0.00000 | 0.00000 | 0.00000 | 0.00000 |-0.01800 |-0.00334 |-0.11000 | 0.30000 | 0.00000 | 0.00000 | 0.00000 |-0.00540 |-0.00100 |-0.03300 |-0.11000 | 0.30000 | 0.00000 | 0.00000 |-0.00162 |-0.00030 |-0.00990 |-0.03300 |-0.11000 | 0.30000 | 0.00000 |-0.00049 |-0.00009 |-0.00297 |-0.00990 |-0.03300 |-0.11000 | 0.30000 |-0.00015 |-0.00003 |-0.00089 |-0.00297 |-0.00990 |-0.03300 |-0.11000

| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | |

|1115085 | 1050 | 1053 | 0.00000 | |1115086 | 1050 | 1054 | 0.00000 | |1115087 | 1050 | 1055 | 0.00000 | |1115088 | 1050 | 1056 | 0.00000 | |1115089 | 1051 | 1049 | 0.03870 | |1115090 | 1051 | 1050 | 0.07500 | |1115091 | 1051 | 1051 | 0.07000 | |1115092 | 1051 | 1052 | 0.30000 | |1115093 | 1051 | 1053 | 0.00000 | |1115094 | 1051 | 1054 | 0.00000 | |1115095 | 1051 | 1055 | 0.00000 | |1115096 | 1051 | 1056 | 0.00000 | |1115097 | 1052 | 1049 | 0.01647 | |1115098 | 1052 | 1050 | 0.03870 | |1115099 | 1052 | 1051 | 0.07500 | |1115100 | 1052 | 1052 | 0.07000 | |1115101 | 1052 | 1053 | 0.30000 | |1115102 | 1052 | 1054 | 0.00000 | |1115103 | 1052 | 1055 | 0.00000 | |1115104 | 1052 | 1056 | 0.00000 | |1115105 | 1053 | 1049 | 0.00640 | |1115106 | 1053 | 1050 | 0.01647 | |1115107 | 1053 | 1051 | 0.03870 | |1115108 | 1053 | 1052 | 0.07500 | |1115109 | 1053 | 1053 | 0.07000 | |1115110 | 1053 | 1054 | 0.30000 | |1115111 | 1053 | 1055 | 0.00000 | |1115112 | 1053 | 1056 | 0.00000 | |1115113 | 1054 | 1049 | 0.00236 | |1115114 | 1054 | 1050 | 0.00640 | |1115115 | 1054 | 1051 | 0.01647 | |1115116 | 1054 | 1052 | 0.03870 | |1115117 | 1054 | 1053 | 0.07500 | |1115118 | 1054 | 1054 | 0.07000 | |1115119 | 1054 | 1055 | 0.30000 | |1115120 | 1054 | 1056 | 0.00000 | |1115121 | 1055 | 1049 | 0.00084 | |1115122 | 1055 | 1050 | 0.00236 | |1115123 | 1055 | 1051 | 0.00640 | |1115124 | 1055 | 1052 | 0.01647 | |1115125 | 1055 | 1053 | 0.03870 | |1115126 | 1055 | 1054 | 0.07500 | |1115127 | 1055 | 1055 | 0.07000 | |1115128 | 1055 | 1056 | 0.30000 | |1115129 | 1056 | 1049 | 0.00063 | |1115130 | 1056 | 1050 | 0.00181 | |1115131 | 1056 | 1051 | 0.00509 | |1115132 | 1056 | 1052 | 0.01386 | |1115133 | 1056 | 1053 | 0.03583 | |1115134 | 1056 | 1054 | 0.08486 | |1115135 | 1056 | 1055 | 0.16768 | |1115136 | 1056 | 1056 | 0.19680 | --------------------------------------Spectral radius matrik M=max(abs(eig(M)) spectral_radius = 1.0101 >>

113

LAMPIRAN 7

Proses iterasi untuk contoh 4.5 » BENTUK UMUM : Uxx + Uyy = f -----------------------------Jenis Persamaan Differensial: 1. Laplace 2. Poisson * pilih 1-2 : ? 1 UKURAN DOMAIN -------------lebar = ? 1 tinggi = ? 1 I:indeks maks.pilahan pd arah sb-x. i=1,2,.,I. J:indeks maks.pilahan pd arah sb-y. j=1,2,.,J.

SYARAT BATAS -----------Syarat batas sisi ke 1 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 u = 1 Syarat batas sisi ke 2 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dy = 2 Syarat batas sisi ke 3 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 3 dU/dx=alfa.u+betta alfa =3 betta = 1

I = ? 33 J = ? 33

---------------------------------------Ganti atribut iterasi (y/t) = ? t ---------------------------------------Parameter SOR : w=1.2 Nilai toleransi =0.000001 -----------Proses iterasi SOR-------------------Iterasi ke 1 banyak U konvergen 33 // 1089 Iterasi ke 2 banyak U konvergen 33 // 1089 Iterasi ke 3 banyak U konvergen 33 // 1089 ... Iterasi ke 262 banyak U konvergen 798 // 1089 Iterasi ke 263 banyak U konvergen 848 // 1089 Iterasi ke 264 banyak U konvergen 903 // 1089 Iterasi ke 265 banyak U konvergen 971 // 1089 Iterasi ke 266 banyak U konvergen 1037 // 1089 Iterasi ke 267 banyak U konvergen 1089 // 1089 -------------------------------------------------Hasil dalam bentuk tabel dan grafik --------------------------------------------------

Syarat batas sisi ke 4 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 3 dU/dx=alfa.u+betta alfa =-4 betta = 1

114

Nilai u sesuai letaknya pada Domain ------------------------------------0.63832

0.72247

0.79764

0.86579

0.92791 ...

1.12587

1.08054

1.02961

0.97310

0.91296

0.58276

0.66446

0.73865

0.80630

0.86812 ...

1.06601

1.02084

0.96991

0.91241

0.84629

0.53556

0.61396

0.68621

0.75265

0.81366 ...

1.01143

0.96690

0.91677

0.86034

0.79647

0.49448

0.56961

0.63958

0.70441

0.76427 ...

0.96200

0.91857

0.86993

0.81571

0.75551

0.45832

0.53044

0.59808

0.66115

0.71967 ...

0.91746

0.87543

0.82867

0.77706

0.72054

0.42637

0.49575

0.56115

0.62243

0.67953 ...

0.87747

0.83701

0.79229

0.74330

0.69018

0.39812

0.46503

0.52834

0.58789

0.64357 ...

0.84169

0.80286

0.76016

0.71368

0.66361

0.37318

0.43789

0.49930

0.55722

0.61153 ...

0.80981

0.77258

0.73181

0.68765

0.64030

0.35130

0.41407

0.47375

0.53017

0.58318 ...

0.78154

0.74583

0.70687

0.66480

0.61986

0.33225

0.39331

0.45147

0.50653

0.55834 ...

0.75663

0.72234

0.68502

0.64484

0.60202

0.31588

0.37547

0.43228

0.48613

0.53686 ...

0.73487

0.70187

0.66602

0.62751

0.58654

0.30207

0.36040

0.41606

0.46884

0.51860 ...

0.71610

0.68425

0.64970

0.61264

0.57326

0.29072

0.34802

0.40270

0.45457

0.50349 ...

0.70016

0.66932

0.63590

0.60008

0.56206

0.28178

0.33826

0.39214

0.44326

0.49146 ...

0.68694

0.65697

0.62450

0.58972

0.55281

0.27522

0.33108

0.38436

0.43487

0.48247 ...

0.67636

0.64710

0.61542

0.58148

0.54547

0.27103

0.32649

0.37934

0.42939

0.47651 ...

0.66834

0.63967

0.60860

0.57531

0.53997

0.26925

0.32451

0.37710

0.42685

0.47360 ...

0.66286

0.63463

0.60401

0.57117

0.53630

0.26993

0.32521

0.37772

0.42729

0.47378 ...

0.65989

0.63197

0.60165

0.56908

0.53444

0.27315

0.32867

0.38128

0.43081

0.47713 ...

0.65945

0.63172

0.60153

0.56904

0.53444

0.27906

0.33504

0.38793

0.43755

0.48378 ...

0.66158

0.63391

0.60371

0.57112

0.53632

0.28782

0.34451

0.39784

0.44766

0.49386 ...

0.66633

0.63864

0.60829

0.57542

0.54020

0.29968

0.35732

0.41128

0.46139

0.50759 ...

0.67382

0.64602

0.61539

0.58206

0.54618

0.31496

0.37382

0.42856

0.47903

0.52522 ...

0.68418

0.65622

0.62521

0.59123

0.55445

0.33406

0.39443

0.45010

0.50096

0.54706 ...

0.69759

0.66948

0.63801

0.60322

0.56525

0.35755

0.41975

0.47647

0.52766

0.57349 ...

0.71431

0.68609

0.65413

0.61837

0.57894

0.38618

0.45054

0.50836

0.55971

0.60496 ...

0.73461

0.70644

0.67404

0.63721

0.59598

0.42100

0.48787

0.54671

0.59785

0.64200 ...

0.75884

0.73104

0.69838

0.66044

0.61708

0.46353

0.53324

0.59277

0.64299

0.68518 ...

0.78740

0.76049

0.72800

0.68908

0.64324

0.51603

0.58879

0.64815

0.69614

0.73506 ...

0.82066

0.79553

0.76404

0.72464

0.67603

0.58228

0.65772

0.71490

0.75836

0.79201 ...

0.85892

0.83692

0.80799

0.76943

0.71808

0.66933

0.74493

0.79534

0.83041

0.85602 ...

0.90216

0.88523

0.86157

0.82699

0.77447

0.79318

0.85731

0.89114

0.91192

0.92612 ...

0.94974

0.94026

0.92609

0.90250

0.85694

1.00000

1.00000

1.00000

1.00000

1.00000 ...

1.00000

1.00000

1.00000

1.00000

1.00000

TABEL SOLUSI NUMERIK DENGAN ERROR RELATIF --------------------------------------------------------------Solusi numerik Error relatif no x y u R --------------------------------------------------------------| | | | | | | | | | | | | | | | | |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18

| | | | | | | | | | | | | | | | | |

0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000

| | | | | | | | | | | | | | | | | |

0.00000 0.03125 0.06250 0.09375 0.12500 0.15625 0.18750 0.21875 0.25000 0.28125 0.31250 0.34375 0.37500 0.40625 0.43750 0.46875 0.50000 0.53125

| | | | | | | | | | | | | | | | | |

1.000000 00 0.79318447 0.66933100 0.58228233 0.51603062 0.46352551 0.42100398 0.38618037 0.35754758 0.33405765 0.31495539 0.29968336 0.28782361 0.27906034 0.27315493 0.26992888 0.26925190 0.27103364

| | | | | | | | | | | | | | | | | |

0.00000 000 -0.00014127 -0.00021465 -0.00028445 -0.00035062 -0.00041311 -0.00047190 -0.00052696 -0.00057830 -0.00062592 -0.00066984 -0.00071008 -0.00074668 -0.00077970 -0.00080918 -0.00083519 -0.00085779 -0.00087708

| | | | | | | | | | | | | | | | | |

115

| 19 | 0.00000 | 0.56250 | 0.27521792 | -0.00089313 | | 20 | 0.00000 | 0.59375 | 0.28177894 | -0.00090604 | | 21 | 0.00000 | 0.62500 | 0.29071913 | -0.00091590 | | 22 | 0.00000 | 0.65625 | 0.30206828 | -0.00092281 | | 23 | 0.00000 | 0.68750 | 0.31588411 | -0.00092688 | | 24 | 0.00000 | 0.71875 | 0.33225428 | -0.00092821 | | 25 | 0.00000 | 0.75000 | 0.35130022 | -0.00092693 | | 26 | 0.00000 | 0.78125 | 0.37318342 | -0.00092313 | | 27 | 0.00000 | 0.81250 | 0.39811527 | -0.00091694 | | 28 | 0.00000 | 0.84375 | 0.42637270 | -0.00090848 | | 29 | 0.00000 | 0.87500 | 0.45832408 | -0.00089786 | | 30 | 0.00000 | 0.90625 | 0.49447506 | -0.00088521 | | 31 | 0.00000 | 0.93750 | 0.53556024 | -0.00087064 | | 32 | 0.00000 | 0.96875 | 0.58275900 | -0.00085426 | | 33 | 0.00000 | 1.00000 | 0.63831888 | -0.00016341 | | 34 | 0.03125 | 0.00000 | 1.00000000 | 0.00000000 | | 35 | 0.03125 | 0.03125 | 0.85731442 | -0.00008250 | | 36 | 0.03125 | 0.06250 | 0.74492827 | -0.00009539 | | 37 | 0.03125 | 0.09375 | 0.65772269 | -0.00010749 | | 38 | 0.03125 | 0.12500 | 0.58878501 | -0.00011879 | | 39 | 0.03125 | 0.15625 | 0.53323904 | -0.00012929 | | 40 | 0.03125 | 0.18750 | 0.48787392 | -0.00013900 | | 41 | 0.03125 | 0.21875 | 0.45053911 | -0.00014792 | | 42 | 0.03125 | 0.25000 | 0.41974596 | -0.00015606 | | 43 | 0.03125 | 0.28125 | 0.39443141 | -0.00016343 | | 44 | 0.03125 | 0.31250 | 0.37381703 | -0.00017003 | | 45 | 0.03125 | 0.34375 | 0.35732219 | -0.00017589 | | 46 | 0.03125 | 0.37500 | 0.34450847 | -0.00018102 | | 47 | 0.03125 | 0.40625 | 0.33504293 | -0.00018542 | | 48 | 0.03125 | 0.43750 | 0.32867313 | -0.00018913 | | 49 | 0.03125 | 0.46875 | 0.32520977 | -0.00019216 | ... | 1075 | 1.00000 | 0.56250 | 0.54546914 | 0.00025982 | | 1076 | 1.00000 | 0.59375 | 0.55281430 | 0.00026447 | | 1077 | 1.00000 | 0.62500 | 0.56205538 | 0.00026797 | | 1078 | 1.00000 | 0.65625 | 0.57326494 | 0.00027031 | | 1079 | 1.00000 | 0.68750 | 0.58654256 | 0.00027150 | | 1080 | 1.00000 | 0.71875 | 0.60201926 | 0.00027151 | | 1081 | 1.00000 | 0.75000 | 0.61986470 | 0.00027034 | | 1082 | 1.00000 | 0.78125 | 0.64029826 | 0.00026795 | | 1083 | 1.00000 | 0.81250 | 0.66360689 | 0.00026430 | | 1084 | 1.00000 | 0.84375 | 0.69017513 | 0.00025938 | | 1085 | 1.00000 | 0.87500 | 0.72054062 | 0.00025329 | | 1086 | 1.00000 | 0.90625 | 0.75551165 | 0.00024675 | | 1087 | 1.00000 | 0.93750 | 0.79646646 | 0.00024328 | | 1088 | 1.00000 | 0.96875 | 0.84629359 | 0.00025961 | | 1089 | 1.00000 | 1.00000 | 0.91296165 | 0.00015760 | -------------------------------------------------------------Keterangan : maksimum absolut Error relatif : maks|R(x,y)|= terjadi pada koordinat x= 0.00000 y= 0.71875 dengan besar suhu U(x,y)= 0.3323 dan error relatif R(x,y)= -9.2821e-004

9.2821e-004

Selesai »

116

LAMPIRAN 8

Proses iterasi untuk contoh 4.6 » BENTUK UMUM : Uxx + Uyy = f ------------------------------Jenis Persamaan Differensial: 1. Laplace 2. Poisson * pilih 1-2 : ? 1 UKURAN DOMAIN -------------lebar = ? 1 tinggi = ? 1 I:indeks maks.pilahan pd arah sb-x. i=1,2,.,I. J:indeks maks.pilahan pd arah sb-y. j=1,2,.,J.

SYARAT BATAS -----------Syarat batas sisi ke 1 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 u = 1 Syarat batas sisi ke 2 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dy = 2 Syarat batas sisi ke 3 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 3 dU/dx=alfa.u+betta alfa =-3 betta = 1 Syarat batas sisi ke 4 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 3 dU/dx=alfa.u+betta alfa =4 betta = 1

I = ? 33 J = ? 33

-----------------------------------------Ganti atribut iterasi (y/t) = ? t -----------------------------------------Parameter SOR : w=1.2 Nilai toleransi =0.000001 -----------Proses iterasi SOR-------------------Iterasi ke 1 banyak U konvergen 33 // 1089 Iterasi ke 2 banyak U konvergen 33 // 1089 Iterasi ke 3 banyak U konvergen 33 // 1089 Iterasi ke 4 banyak U konvergen 33 // 1089 .... Iterasi ke 186 banyak U konvergen 33 // 1089 Iterasi ke 187 banyak U konvergen 33 // 1089 Iterasi ke 188 banyak U konvergen 33 // 1089 Iterasi ke 189 banyak U konvergen 33 // 1089 Iterasi ke 190 banyak U konvergen 33 // 1089 Iterasi ke 191 banyak U konvergen 33 // 1089 Iterasi ke 192 banyak U konvergen 33 // 1089 Iterasi ke 193 banyak U konvergen 33 // 1089 Iterasi ke 194 banyak U konvergen 33 // 1089 Iterasi ke 195 banyak U konvergen 33 // 1089 Iterasi ke 196 banyak U konvergen 33 // 1089 Iterasi ke 197 banyak U konvergen 33 // 1089 Iterasi ke 198 banyak U konvergen 33 // 1089 Iterasi ke 199 banyak U konvergen 33 // 1089 Iterasi ke 200 banyak U konvergen 33 // 1089 Hentikan program (y/t) = ? y

Selesai »

117

LAMPIRAN 9

Contoh 4.7. Konduksi panas pada Lempengen logam yang tidak mempunyai sumber panas

internal

» Persamaan Panas steady state: Uxx + Uyy = f ---------------------------------------------Jenis: 1. Laplace (tidak ada sumber panas internal) 2. Poisson (ada sumber panas internal) * pilih 1-2 : ? 1 UKURAN PLAT LOGAM -----------------lebar = ? 1 tinggi = ? 1 Level grid. L =5 Banyak Titik Grid : 33 X 33

= 1089

* Suhu leleh logam : ? 1200 -------------SYARAT BATAS -----------Syarat batas sisi ke 1 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dy = 0 Syarat batas sisi ke 2 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 u ? = 1000 Syarat batas sisi ke 3 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 u ? = 500 Syarat batas sisi ke 4 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 2 du/dx = 0 w=1.8311 Nilai toleransi =0.000001

118

----------------------------------------------ingin mengganti atribut iterasi (y/t) = ? t -------------------------------------------------------------Proses iterasi SOR --------------Iterasi ke 1 banyak U konvergen 65 // 1089 Iterasi ke 2 banyak U konvergen 65 // 1089 ... Iterasi ke 558 banyak U konvergen 1034 // 1089 Iterasi ke 559 banyak U konvergen 1052 // 1089 Iterasi ke 560 banyak U konvergen 1073 // 1089 Iterasi ke 561 banyak U konvergen 1087 // 1089 Iterasi ke 562 banyak U konvergen 1089 // 1089 ------------------------------------------------Hasil dalam bentuk tabel dan grafik --------------------------------------------------BESARNYA SUHU SESUAI LETAKNYA PADA LEMPENGAN LOGAM ---------------------------------------------------------------------------------750.00000 1000.00000 1000.00000 1000.00000 1000.00000 ... 1000.00000 1000.00000 1000.00000 1000.00000 1000.00000 500.00000 750.00000 848.82582 895.30330 920.74697 ...

986.60542 986.76193 986.87259 986.93854 986.96045

500.00000 651.17418 750.00000 811.64040 851.25830 ...

973.25835 973.56972 973.78989 973.92112 973.96471

500.00000 604.69670 688.35960 750.00000 794.65391 ...

960.00571 960.46870 960.79615 960.99132 961.05617

500.00000 579.25303 648.74170 705.34609 750.00000 ...

946.89337 947.50324 947.93466 948.19186 948.27732

500.00000 563.57371 622.00809 672.64264 715.04876 ...

933.96555 934.71622 935.24741 935.56414 935.66939

500.00000 553.03373 603.07432 648.16762 687.55240 ...

921.26447 922.14869 922.77460 923.14790 923.27196

500.00000 545.48689 589.08784 629.40112 665.67570 ...

908.82992 909.83948 910.55440 910.98090 911.12266

500.00000 539.82600 578.38903 614.67333 648.03164 ...

896.69896 897.82492 898.62261 899.09863 899.25687

500.00000 535.42807 569.96894 602.87153 633.60239 ...

884.90572 886.13864 887.01250 887.53414 887.70757

500.00000 531.91736 563.18712 593.24145 621.64641 ...

873.48127 874.81139 875.75462 876.31785 876.50514

500.00000 529.05427 557.62071 585.26076 611.62128 ...

862.45352 863.87103 864.87675 865.47749 865.67730

500.00000 526.67899 552.98071 578.55959 603.12648 ...

851.84725 853.34246 854.40386 855.03808 855.24905

500.00000 524.68099 549.06355 572.87040 595.86265 ...

841.68422 843.24771 844.35816 845.02191 845.24274

500.00000 522.98142 545.72209 567.99581 589.60319 ...

831.98322 833.60599 834.75914 835.44866 835.67811

500.00000 521.52261 542.84758 563.78757 584.17437 ...

822.76030 824.43389 825.62377 826.33548 826.57236

500.00000 520.26143 540.35806 560.13251 579.44139 ...

814.02894 815.74551 816.96656 817.69715 817.94035

500.00000 519.16507 538.19071 556.94302 575.29844 ...

805.80027 807.55267 808.79979 809.54622 809.79473

500.00000 518.20812 536.29671 554.15041 571.66166 ...

798.08333 799.86509 801.13370 801.89322 802.14613

500.00000 517.37068 534.63762 551.70023 568.46396 ...

790.88527 792.69068 793.97668 794.74685 795.00333

500.00000 516.63699 533.18286 549.54894 565.65126 ...

784.21162 786.03569 787.33551 788.11415 788.37350

500.00000 515.99442 531.90789 547.66140 563.17965 ...

778.06652 779.90495 781.21549 782.00077 782.26235

500.00000 515.43280 530.79289 546.00912 561.01332 ...

772.45292 774.30209 775.62075 776.41107 776.67437

500.00000 514.94387 529.82175 544.56888 559.12294 ...

767.37283 769.22975 770.55436 771.34841 771.61297

500.00000 514.52095 528.98136 543.32171 557.48448 ...

762.82746 764.68972 766.01851 766.81522 767.08070

500.00000 514.15857 528.26102 542.25210 556.07823 ...

758.81744 760.68316 762.01475 762.81328 763.07938

500.00000 513.85230 527.65206 541.34745 554.88809 ...

755.34296 757.21072 758.54406 759.34376 759.61027

500.00000 513.59858 527.14746 540.59755 553.90103 ...

752.40390 754.27270 755.60703 756.40741 756.67417

500.00000 513.39457 526.74165 539.99428 553.10662 ...

749.99997 751.86915 753.20395 754.00469 754.27159

500.00000 513.23806 526.43027 539.53129 552.49675 ...

748.13079 749.99997 751.33493 752.13583 752.40278

500.00000 513.12740 526.21010 539.20384 552.06532 ...

746.79598 748.66501 749.99997 750.80091 751.06789

500.00000 513.06146 526.07887 539.00866 551.80812 ...

745.99524 747.86410 749.19902 749.99997 750.26695

500.00000 513.03955 526.03528 538.94382 551.72266 ...

745.72835 747.59715 748.93204 749.73299 749.99997

119

TABEL SUHU PADA LEMPENGAN LOGAM ---------------------------------------------------------------------Suhu Error no x y u relatif Keterangan ---------------------------------------------------------------------| 1 | 0.00000 | 0.00000 |500.00000000 | 0.00000000 | ----- | | 2 | 0.00000 | 0.03125 |500.00000000 | 0.00000000 | ----- | | 3 | 0.00000 | 0.06250 |500.00000000 | 0.00000000 | ----- | | 4 | 0.00000 | 0.09375 |500.00000000 | 0.00000000 | ----- | ... | 1062 | 1.00000 | 0.15625 |756.67416988 | -0.00072865 | ----- | | 1063 | 1.00000 | 0.18750 |759.61026579 | -0.00071088 | ----- | | 1064 | 1.00000 | 0.21875 |763.07938110 | -0.00069177 | ----- | | 1065 | 1.00000 | 0.25000 |767.08069977 | -0.00068603 | ----- | | 1066 | 1.00000 | 0.28125 |771.61296972 | -0.00066405 | ----- | | 1067 | 1.00000 | 0.31250 |776.67436707 | -0.00064097 | ----- | | 1068 | 1.00000 | 0.34375 |782.26235047 | -0.00061689 | ----- | | 1069 | 1.00000 | 0.37500 |788.37349775 | -0.00059189 | ----- | | 1070 | 1.00000 | 0.40625 |795.00333036 | -0.00056607 | ----- | | 1071 | 1.00000 | 0.43750 |802.14612714 | -0.00053951 | ----- | | 1072 | 1.00000 | 0.46875 |809.79473045 | -0.00051230 | ----- | | 1073 | 1.00000 | 0.50000 |817.94034818 | -0.00048454 | ----- | | 1074 | 1.00000 | 0.53125 |826.57235614 | -0.00045630 | ----- | | 1075 | 1.00000 | 0.56250 |835.67810594 | -0.00042767 | ----- | | 1076 | 1.00000 | 0.59375 |845.24274420 | -0.00039875 | ----- | | 1077 | 1.00000 | 0.62500 |855.24904969 | -0.00036960 | ----- | | 1078 | 1.00000 | 0.65625 |865.67729533 | -0.00034031 | ----- | | 1079 | 1.00000 | 0.68750 |876.50514234 | -0.00031097 | ----- | | 1080 | 1.00000 | 0.71875 |887.70757356 | -0.00028165 | ----- | | 1081 | 1.00000 | 0.75000 |899.25687296 | -0.00025241 | ----- | | 1082 | 1.00000 | 0.78125 |911.12265683 | -0.00022335 | ----- | | 1083 | 1.00000 | 0.81250 |923.27196159 | -0.00019452 | ----- | | 1084 | 1.00000 | 0.84375 |935.66939073 | -0.00016599 | ----- | | 1085 | 1.00000 | 0.87500 |948.27732176 | -0.00013783 | ----- | | 1086 | 1.00000 | 0.90625 |961.05617128 | -0.00011010 | ----- | | 1087 | 1.00000 | 0.93750 |973.96471396 | -0.00008286 | ----- | | 1088 | 1.00000 | 0.96875 |986.96044829 | -0.00005617 | ----- | | 1089 | 1.00000 | 1.00000 |1000.00000000| 0.00000000 | ----- | ---------------------------------------------------------------------Keterangan ----- : tidak meleleh maksimum absolut Error relatif : maks|R(x,y)|= terjadi pada koordinat x= 0.78125 y= 0.00000 dengan besar suhu : u(x,y) = 736.9206 dan error relatif : R(x,y) = 0.0012

0.0012

Selesai »

120

LAMPIRAN 10

10A) Mengaproksimasi bentuk lelehan logam menggunakan grid IxJ = 33 x 33 » Persamaan Panas steady state: Uxx + Uyy = f ---------------------------------------------Jenis: 1. Laplace (tidak ada sumber panas internal) 2. Poisson (ada sumber panas internal) * pilih 1-2 : ? 2 f = ? -20000 UKURAN PLAT LOGAM -----------------lebar = ? 1 tinggi = ? 1 Level grid. L =5 Banyak Titik Grid : 33 X 33

= 1089

* Suhu leleh logam : ? 1200 -------------SYARAT BATAS -----------Syarat batas sisi ke 1 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 u ? = 100 Syarat batas sisi ke 2 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 u ? = 100 Syarat batas sisi ke 3 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 u ? = 100 Syarat batas sisi ke 4 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 u ? = 100 w=1.8311 Nilai toleransi =0.000001 ----------------------------------------------ingin mengganti atribut iterasi (y/t) = ? t ----------------------------------------------...

121

TABEL SUHU PADA LEMPENGAN LOGAM ---------------------------------------------------------------------Suhu Error no x y u relatif Keterangan ---------------------------------------------------------------------... | | | | ... | | | | | | | | | | | | | | | | | | | | ... | | | | | | | | | | | | | | | | | | | ... | | | | | | | | | |

29 30 487 488

| | | |

0.00000 0.00000 0.43750 0.43750

| | | |

0.05469 0.05664 0.75000 0.78125

|100.00000000 |100.00000000 |1231.54936794 |1141.45840458

| | | |

0.00000000 0.00000000 0.00003070 0.00001168

| | | |

----- | ----- | meleleh | ----- |

502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521

| | | | | | | | | | | | | | | | | | | |

0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875 0.46875

| | | | | | | | | | | | | | | | | | | |

0.18750 0.21875 0.25000 0.28125 0.31250 0.34375 0.37500 0.40625 0.43750 0.46875 0.50000 0.53125 0.56250 0.59375 0.62500 0.65625 0.68750 0.71875 0.75000 0.78125

|1046.92468660 |1151.12874321 |1242.26040914 |1321.01419169 |1388.00793086 |1443.77994477 |1488.78718566 |1523.40415200 |1547.92232634 |1562.54994355 |1567.41193573 |1562.54994327 |1547.92232585 |1523.40415131 |1488.78718478 |1443.77994379 |1388.00792984 |1321.01419067 |1242.26040812 |1151.12874224

| | | | | | | | | | | | | | | | | | | |

-0.00001032 0.00000044 -0.00008682 -0.00003363 -0.00011423 -0.00013594 -0.00014293 -0.00006410 -0.00003574 -0.00006615 -0.00001000 0.00003538 0.00000130 0.00000232 0.00001603 0.00000540 0.00001169 0.00001305 0.00003675 0.00000536

| | | | | | | | | | | | | | | | | | | |

----- | ----- | meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh ----- |

536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554

| | | | | | | | | | | | | | | | | | |

0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000 0.50000

| | | | | | | | | | | | | | | | | | |

0.21875 0.25000 0.28125 0.31250 0.34375 0.37500 0.40625 0.43750 0.46875 0.50000 0.53125 0.56250 0.59375 0.62500 0.65625 0.68750 0.71875 0.75000 0.78125

|1154.34022156 |1245.81808264 |1324.88004085 |1392.14244747 |1448.14263750 |1493.33696310 |1528.09959371 |1552.72185796 |1567.41193573 |1572.29474812 |1567.41193549 |1552.72185752 |1528.09959307 |1493.33696226 |1448.14263655 |1392.14244649 |1324.88003985 |1245.81808166 |1154.34022061

| | | | | | | | | | | | | | | | | | |

-0.00000892 -0.00001743 -0.00005660 -0.00001206 -0.00011238 -0.00011908 -0.00004546 0.00001283 -0.00001000 -0.00004089 0.00004264 0.00000679 -0.00001392 -0.00000068 -0.00000725 0.00001186 0.00002668 0.00000521 0.00000902

| | | | | | | | | | | | | | | | | | |

----- | meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh ----- |

569 570 571 572 573 574 575 576 577 578

| | | | | | | | | |

0.53125 0.53125 0.53125 0.53125 0.53125 0.53125 0.53125 0.53125 0.53125 0.53125

| | | | | | | | | |

0.21875 0.25000 0.28125 0.31250 0.34375 0.37500 0.40625 0.43750 0.46875 0.50000

|1151.12874312 |1242.26040899 |1321.01419153 |1388.00793068 |1443.77994455 |1488.78718542 |1523.40415174 |1547.92232606 |1562.54994327 |1567.41193549

| | | | | | | | | |

-0.00007754 0.00000059 -0.00000344 -0.00007267 -0.00005431 -0.00007219 0.00002130 0.00003022 0.00003538 0.00004264

| | | | | | | | | |

----- | meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh meleleh

| | | | | | | | | | | | | | | | |

| | | | | | | | | | | | | | | | |

| | | | | | | | |

122

|

579 |

0.53125 |

0.53125 |1562.54994308 | -0.00000449 | meleleh |

| |

580 | 581 |

0.53125 | 0.53125 |

0.56250 |1547.92232567 | 0.59375 |1523.40415115 |

0.00002008 | meleleh | 0.00001435 | meleleh |

|

582 |

0.53125 |

0.62500 |1488.78718465 |

0.00000648 | meleleh |

|

583 |

0.53125 |

0.65625 |1443.77994367 |

0.00001559 | meleleh |

|

584 |

0.53125 |

0.68750 |1388.00792974 |

0.00000893 | meleleh |

|

585 |

0.53125 |

0.71875 |1321.01419059 |

0.00001832 | meleleh |

| |

586 | 587 |

0.53125 | 0.53125 |

0.75000 |1242.26040807 | 0.00000146 | meleleh | 0.78125 |1151.12874220 | -0.00000941 | ----- |

|

601 |

0.56250 |

0.18750 |1038.40349727 | -0.00008504 | ----- |

|

602 |

0.56250 |

0.21875 |1141.45840534 | -0.00000797 | ----- |

|

603 |

0.56250 |

0.25000 |1231.54936869 | -0.00002846 | meleleh |

|

604 |

0.56250 |

0.28125 |1309.37713560 | -0.00000993 | meleleh |

|

605 |

0.56250 |

0.31250 |1375.56388908 | -0.00004369 | meleleh |

| |

606 | 607 |

0.56250 | 0.56250 |

0.34375 |1430.65077455 | -0.00006205 | meleleh | 0.37500 |1475.09643221 | -0.00004113 | meleleh |

|

608 |

0.56250 |

0.40625 |1509.27625177 | -0.00001606 | meleleh |

|

609 |

0.56250 |

0.43750 |1533.48210131 |

0.00000616 | meleleh |

|

610 |

0.56250 |

0.46875 |1547.92232585 |

0.00000130 | meleleh |

|

611 |

0.56250 |

0.50000 |1552.72185752 |

0.00000679 | meleleh |

|

612 |

0.56250 |

0.53125 |1547.92232567 |

0.00002008 | meleleh |

|

613 |

0.56250 |

0.56250 |1533.48210095 |

0.00001260 | meleleh |

|

614 |

0.56250 |

0.59375 |1509.27625124 |

0.00001553 | meleleh |

|

615 |

0.56250 |

0.62500 |1475.09643152 |

0.00001198 | meleleh |

|

616 |

0.56250 |

0.65625 |1430.65077374 |

0.00002646 | meleleh |

|

617 |

0.56250 |

0.68750 |1375.56388823 |

0.00001697 | meleleh |

|

618 |

0.56250 |

0.71875 |1309.37713473 |

0.00001124 | meleleh |

|

619 |

0.56250 |

0.75000 |1231.54936783 |

0.00000913 | meleleh |

|

620 |

0.56250 |

0.78125 |1141.45840449 |

0.00000461 | ----- |

|

635 |

0.59375 |

0.21875 |1125.22076228 | -0.00002178 | ----- |

|

636 |

0.59375 |

0.25000 |1213.57027480 | -0.00000825 | meleleh |

|

637 |

0.59375 |

0.28125 |1289.84984307 | -0.00003001 | meleleh |

| |

638 | 639 |

0.59375 | 0.59375 |

0.31250 |1354.68846546 | 0.00001908 | meleleh | 0.34375 |1408.63158228 | -0.00004442 | meleleh |

|

640 |

0.59375 |

0.37500 |1452.14026708 | -0.00004937 | meleleh |

|

641 |

0.59375 |

0.40625 |1485.59107182 |

0.00000351 | meleleh |

| |

642 | 643 |

0.59375 | 0.59375 |

0.43750 |1509.27625155 | 0.46875 |1523.40415131 |

0.00001574 | meleleh | 0.00000232 | meleleh |

|

644 |

0.59375 |

0.50000 |1528.09959307 | -0.00001392 | meleleh |

|

645 |

0.59375 |

0.53125 |1523.40415115 |

0.00001435 | meleleh |

|

646 |

0.59375 |

0.56250 |1509.27625124 |

0.00001553 | meleleh |

|

647 |

0.59375 |

0.59375 |1485.59107135 | -0.00000026 | meleleh |

|

648 |

0.59375 |

0.62500 |1452.14026647 |

0.00000299 | meleleh |

|

649 |

0.59375 |

0.65625 |1408.63158157 |

0.00002033 | meleleh |

|

650 |

0.59375 |

0.68750 |1354.68846471 | -0.00000212 | meleleh |

|

651 |

0.59375 |

0.71875 |1289.84984229 |

|

652 |

0.59375 |

0.75000 |1213.57027402 | -0.00000255 | meleleh |

|

653 |

0.59375 |

0.78125 |1125.22076150 | -0.00001336 | ----- |

|

669 |

0.62500 |

0.25000 |1188.12987517 | -0.00003184 | ----- |

|

670 |

0.62500 |

0.28125 |1262.23224641 |

| |

671 | 672 |

0.62500 | 0.62500 |

0.31250 |1325.17729741 | -0.00001119 | meleleh | 0.34375 |1377.51557200 | -0.00001952 | meleleh |

|

673 |

0.62500 |

0.37500 |1419.71073196 | -0.00003166 | meleleh |

|

674 |

0.62500 |

0.40625 |1452.14026688 | -0.00002516 | meleleh |

|

675 |

0.62500 |

0.43750 |1475.09643179 |

0.00004223 | meleleh |

|

676 |

0.62500 |

0.46875 |1488.78718478 |

0.00001603 | meleleh |

|

677 |

0.62500 |

0.50000 |1493.33696226 | -0.00000068 | meleleh |

|

678 |

0.62500 |

0.53125 |1488.78718465 |

0.00000648 | meleleh |

|

679 |

0.62500 |

0.56250 |1475.09643152 |

0.00001198 | meleleh |

|

680 |

0.62500 |

0.59375 |1452.14026647 |

0.00000299 | meleleh |

|

681 |

0.62500 |

0.62500 |1419.71073143 | -0.00002040 | meleleh |

|

682 |

0.62500 |

0.65625 |1377.51557138 |

...

...

0.00000177 | meleleh |

... 0.00001125 | meleleh |

0.00000892 | meleleh |

123

|

683 |

0.62500 |

0.68750 |1325.17729674 | -0.00000980 | meleleh |

| |

684 | 685 |

0.62500 | 0.62500 |

0.71875 |1262.23224572 | -0.00000206 | meleleh | 0.75000 |1188.12987446 | -0.00001469 | ----- |

|

701 |

0.65625 |

0.21875 |1072.23199573 | -0.00001717 | ----- |

|

702 |

0.65625 |

0.25000 |1154.95314989 |

|

703 |

0.65625 |

0.28125 |1226.24072000 | -0.00003448 | meleleh |

| |

704 | 705 |

0.65625 | 0.65625 |

0.31250 |1286.74165576 | -0.00003669 | meleleh | 0.34375 |1337.01142632 | -0.00005414 | meleleh |

|

706 |

0.65625 |

0.37500 |1377.51557184 | -0.00003168 | meleleh |

|

707 |

0.65625 |

0.40625 |1408.63158192 | -0.00000327 | meleleh |

|

708 |

0.65625 |

0.43750 |1430.65077397 |

0.00003021 | meleleh |

|

709 |

0.65625 |

0.46875 |1443.77994379 |

0.00000540 | meleleh |

|

710 |

0.65625 |

0.50000 |1448.14263655 | -0.00000725 | meleleh |

|

711 |

0.65625 |

0.53125 |1443.77994367 |

0.00001559 | meleleh |

| |

712 | 713 |

0.65625 | 0.65625 |

0.56250 |1430.65077374 | 0.59375 |1408.63158157 |

0.00002646 | meleleh | 0.00002033 | meleleh |

|

714 |

0.65625 |

0.62500 |1377.51557138 |

0.00000892 | meleleh |

|

715 |

0.65625 |

0.65625 |1337.01142577 | -0.00002037 | meleleh |

|

716 |

0.65625 |

0.68750 |1286.74165516 |

|

717 |

0.65625 |

0.71875 |1226.24071938 | -0.00000980 | meleleh |

|

718 |

0.65625 |

0.75000 |1154.95314926 |

|

736 |

0.68750 |

0.28125 |1181.50457791 | -0.00001049 | ----- |

|

737 |

0.68750 |

0.31250 |1239.00592928 | -0.00002956 | meleleh |

|

738 |

0.68750 |

0.34375 |1286.74165563 | -0.00002203 | meleleh |

|

739 |

0.68750 |

0.37500 |1325.17729713 | -0.00002139 | meleleh |

|

740 |

0.68750 |

0.40625 |1354.68846500 |

|

741 |

0.68750 |

0.43750 |1375.56388842 | -0.00000191 | meleleh |

|

742 |

0.68750 |

0.46875 |1388.00792984 |

0.00001169 | meleleh |

|

743 |

0.68750 |

0.50000 |1392.14244649 |

0.00001186 | meleleh |

|

744 |

0.68750 |

0.53125 |1388.00792974 |

0.00000893 | meleleh |

|

745 |

0.68750 |

0.56250 |1375.56388823 |

0.00001697 | meleleh |

|

746 |

0.68750 |

0.59375 |1354.68846471 | -0.00000212 | meleleh |

| |

747 | 748 |

0.68750 | 0.68750 |

0.62500 |1325.17729674 | -0.00000980 | meleleh | 0.65625 |1286.74165516 | 0.00001269 | meleleh |

|

749 |

0.68750 |

0.68750 |1239.00592876 | -0.00000072 | meleleh |

| |

770 | 771 |

0.71875 | 0.71875 |

0.31250 |1181.50457780 | -0.00000752 | ----- | 0.34375 |1226.24071976 | -0.00000815 | meleleh |

|

772 |

0.71875 |

0.37500 |1262.23224604 | -0.00001757 | meleleh |

|

773 |

0.71875 |

0.40625 |1289.84984253 |

0.00001086 | meleleh |

|

774 |

0.71875 |

0.43750 |1309.37713488 |

0.00003094 | meleleh |

|

775 |

0.71875 |

0.46875 |1321.01419067 |

0.00001305 | meleleh |

|

776 |

0.71875 |

0.50000 |1324.88003985 |

0.00002668 | meleleh |

|

777 |

0.71875 |

0.53125 |1321.01419059 |

0.00001832 | meleleh |

|

778 |

0.71875 |

0.56250 |1309.37713473 |

0.00001124 | meleleh |

|

779 |

0.71875 |

0.59375 |1289.84984229 |

0.00000177 | meleleh |

|

780 |

0.71875 |

0.62500 |1262.23224572 | -0.00000206 | meleleh |

|

781 |

0.71875 |

0.65625 |1226.24071938 | -0.00000980 | meleleh |

|

782 |

0.71875 |

0.68750 |1181.50457736 |

|

805 |

0.75000 |

0.37500 |1188.12987471 | -0.00000761 | ----- |

| |

806 | 807 |

0.75000 | 0.75000 |

0.40625 |1213.57027420 | 0.43750 |1231.54936794 |

0.00001870 | meleleh | 0.00003070 | meleleh |

|

808 |

0.75000 |

0.46875 |1242.26040812 |

0.00003675 | meleleh |

|

809 |

0.75000 |

0.50000 |1245.81808166 |

0.00000521 | meleleh |

|

810 |

0.75000 |

0.53125 |1242.26040807 |

0.00000146 | meleleh |

|

811 |

0.75000 |

0.56250 |1231.54936783 |

0.00000913 | meleleh |

|

812 |

0.75000 |

0.59375 |1213.57027402 | -0.00000255 | meleleh |

|

813 |

0.75000 |

0.62500 |1188.12987446 | -0.00001469 | ----- |

| 1082 |

1.00000 |

0.78125 |100.00000000 |

0.00000000 | ----- |

| 1083 |

1.00000 |

0.81250 |100.00000000 |

0.00000000 | ----- |

| 1084 |

1.00000 |

0.84375 |100.00000000 |

0.00000000 | ----- |

... 0.00000671 | ----- |

0.00001269 | meleleh | 0.00000491 | ----- |

...

0.00000925 | meleleh |

...

0.00001370 | ----- |

...

...

124

| 1085 | 1.00000 | 0.87500 |100.00000000 | 0.00000000 | ----- | | 1086 | 1.00000 | 0.90625 |100.00000000 | 0.00000000 | ----- | | 1087 | 1.00000 | 0.93750 |100.00000000 | 0.00000000 | ----- | | 1088 | 1.00000 | 0.96875 |100.00000000 | 0.00000000 | ----- | | 1089 | 1.00000 | 1.00000 |100.00000000 | 0.00000000 | ----- | --------------------------------------------------------------------Keterangan ----- : tidak meleleh maksimum absolut Error relatif : maks|R(x,y)|= terjadi pada koordinat x= 0.18750 y= 0.18750 dengan besar suhu dan error relatif

u(x,y)= R(x,y)=

0.0013

737.1813 -0.0013

10B) Mengaproksimasi bentuk lelehan logam menggunakan grid IxJ = 1025 x 1025 » Persamaan Panas steady state: Uxx + Uyy = f ------------------------------------------Jenis: 1. Laplace (tidak ada sumber panas internal) 2. Poisson (ada sumber panas internal) * pilih 1-2 : ? 2 f = ? -20000 UKURAN PLAT LOGAM -----------------lebar = ? 1 tinggi = ? 1 Level grid. L =10 Banyak Titik Grid : 1025 X 1025

= 1050625

Syarat batas sisi ke 3 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 Suhu : u = 100 Syarat batas sisi ke 4 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 Suhu : u = 100

* Suhu lebur logam : ? 1200 --------------

w=1.9311 Nilai toleransi =0.000001 ----------------------------------...

SYARAT BATAS ------------

Hasil dalam bentuk grafik ---------------------------...

Syarat batas sisi ke 1 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 Suhu : u = 100

>>

Syarat batas sisi ke 2 1. Dirichlet 2. Neumann 3. Campuran (Robyn) * pilih 1-3 : 1 Suhu : u = 100

125

LAMPIRAN 11

Program untuk Contoh 4.1 dengan Nilai Awal Iterasi

u

( 0)

=0

close all; clear all; clc; % PROGRAM_1A.m % Program ini digunakan untuk menyelesaikan PD POISSON 2D % dengan fungsi sumber f dan ukuran domain serta nilai % batasnya dapat diubah. Banyaknya titik grid : I x J dapat diinput % % Program ini menggunakan iterasi SOR. % Tebakan awal iterasi U=0 % % Contoh: % j % y U(x,1)=x^2 + 1 % J 1 ---------------- -----# % . | | % . | | % . | | % . | | % . u(0,y)=y^2 | Uxx + Uyy = 4 |u(1,y)=y^2 + 1 % . | | % . | | % . | | % 3 | | % 2 | | % 1 #-------------- -------+--- x % 0 1 % U(x,0)=x^2 % 1 2 3 ............... .I i %----------------------------------------------------------------------------------close all; clear all; %Bentuk PDP: Uxx + Uyy = f '); f=4; disp(' ');

disp('Dimensi Domain Segi Empat'); disp('--------------------------'); lebar= input(' lebar = ?'); tinggi= input(' tinggi = ?'); disp('--------------------------'); I=input('I:inde ks maksimum searah sb-x. i=1,2,...,I. J=input('J:inde ks maksimum searah sb-y. j=1,2,...,J.

I = ?'); J = ?');

h=lebar/(I-1); % ukuran kisi searah sb-x k=tinggi/(J-1); % ukuran kisi searah sb-y %------------Set syarat batas Dirichlet-------------------%----(Bagian ini dapat diubah sesuai nilai batasnya)------for i=1:I x=(i-1)*h; % Nilai U sisi bawah y=0; U(i,1)=x.^2; % Nilai U sisi atas

% fungsi ini dapat diganti

126

y=tinggi; U(i,J)=x.^2+y.^2; % fungsi ini dapat diganti end A1=U(1,1); B1=U(I,1); C1=U(1,J); D1=U(I,J);

% % % %

sudut sudut sudut sudut

kiri bawah kanan bawah kiri atas kanan atas

for j=1:J y=(j-1)*k; % Nilai U sisi kiri x=0; U(1,j)= y.^2; % fungsi ini dapat diganti % Nilai U sisi kanan x=lebar; U(I,j)= x.^2 + y.^2; % fungsi ini dapat diganti end A2=U(1,1); B2=U(I,1); C2=U(1,J); D2=U(I,J);

% % % %



U(1,1)= (A1+A2)/2; % Nilai rata-rata U sudut kiri bawah U(I,1)= (B1+B2)/2; % Nilai rata-rata U sudut kanan bawah U(1,J)= (C1+C2)/2; % Nilai rata-rata U sudut kiri atas U(I,J)= (D1+D2)/2; % Nilai rata-rata U sudut kanan atas %---------------------end nilai batas-----------------------

%Set tebakan awal titik interior for i=2:I-1 for j=2:J-1 U(i,j)=0; end % for j end % for i Ukonv=0; banyakU=I*J; % banyaknya U yang belum diketahui nilainya epsilon=0.01 % nilai toleransi perbedaan nilai U %Mj=(cos(pi/(I+1)) + cos(pi/(J+1)))/2 w=1.2 % 2/(1+sqrt(1-Mj.^2)) dh=1/(h*h); dk=1/(k*k); dhk= 2*(dh+dk); iter=0; %----------------------- ITERASI SOR---------------------------------while Ukonv1&i1&j1&i
127

if bedaU < epsilon, Ukonv=Ukonv+1 ; end; end % for j end% for i fprintf('\nIterasi ke'); fprintf('%5d ',iter); fprintf('banyak U konvergen');fprintf('%5d ',Ukonv); fprintf(' //');fprintf('%5d',banyakU);

end % while %----------------------END ITERASI SOR------------------------------

% Warning ! % KOORDINAT GAMBAR INI MASIH TERBALIK DAN i,j BELUM DI UBAH X, Y % % ---------------------------------------------------------------

%plot solusi numerik pada domain pcolor(U(1:I,1:J)); colorbar vert; shading interp; colormap(jet); title('Grafik Solusi Numerik'); xlabel('i'); ylabel('j') drawnow; %Hitung error globalnya for i=1:I x=(i-1)*h; for j=1:J y=(j-1)*k; u_eksak(i,j)=x. ^2 + y.^2; end; end; error= U(1:I,1:J)-u_eksak(1:I,1:J);

% dpt diganti sesuai solusi analitisnya

%Plot errornya figure; mesh(error); shading interp; colormap(jet); title('Grafik error global '); xlabel('i'); ylabel('j');zlabel('error global'); colorbar vert; fprintf('\nSelesai\n')

128

LAMPIRAN 12 Program untuk Contoh 4.1 dengan Nilai Awal Iterasi % % % % % % % %

u

(0)

Dimodifikasi

PROGRAM_1B.m Program ini digunakan untuk menyelesaikan PD POISSON 2D dengan fungsi sumber f dan ukuran domain serta nilai batasnya dapat diubah. Level grid L berhubungan dengan banyaknya titik grid: 2^L+1 x 2^L+1 Program ini menggunakan iterasi SOR dengan nilai awal iterasi diseting otomatis mengikuti syarat batasnya. Contoh: ( seperti pada lampiran 11.A )

close all; clear all; % PDP: Uxx + Uyy = f ; f=4; % nilai ini dapat diganti disp('Domain segi empat berukuran :'); disp('----------------------------'); lebar= input(' lebar = '); tinggi=input(' tinggi = '); L=input('Level grid. L ='); disp('-------------------------'); I=2.^L+1; J=2.^L+1; fprintf('Banyak Titik Grid :'); fprintf('%5d ',I);fprintf('X');fprintf('%5d',J); disp(' '); for level=1: L if level>1, %konversi titik interior dari level sebelumnya U(3:2:2*I-3,3:2:2*J-3)=U(2:1:I-1,2:1:J-1); end; %if I=(2.^level)+1; J=(2.^level)+1; h=lebar/(I-1); k=tinggi/(J-1); R=h*h+k*k; dh=1/(h*h);dk=1/(k*k); dhk=2*(dh+dk); %------------Set syarat batas Dirichlet----------------------for i=1:I x=(i-1)*h; % Nilai U sisi bawah y=0; U(i,1)=x.^2; % fungsi ini dapat diganti % Nilai U sisi atas y=tinggi; U(i,J)=x.^2+y.^2; % fungsi ini dapat diganti end A1=U(1,1); % sudut kiri bawah B1=U(I,1); % sudut kanan bawah C1=U(1,J); % sudut kiri atas D1=U(I,J); % sudut kanan atas for j=1:J y=(j-1)*k; % Nilai U sisi kiri x=0; U(1,j)= y.^2; % fungsi ini dapat diganti sesuia nilai batasnya % Nilai U sisi kanan x=lebar; U(I,j)= x.^2 + y.^2; % fungsi ini dapat diganti sesuai nilai batasnya end A2=U(1,1);

% sudut kiri bawah

129

B2=U(I,1); C2=U(1,J); D2=U(I,J); U(1,1)= U(I,1)= U(1,J)= U(I,J)=

% sudut kanan bawah % sudut kiri atas % sudut kanan atas

(A1+A2)/2; (B1+B2)/2; (C1+C2)/2; (D1+D2)/2;

% % % %

Nilai Nilai Nilai Nilai

rata-rata rata-rata rata-rata rata-rata

U U U U



%--------------------end syarat batas-------------------------%---------------set tebakan awal titik interior---------------switch level case 1, U(2,2)=(U(1,1)+U(I,1)+U(1,J)+U(I,J))/4 -R*f/4; otherwise, %hitung U(i,j) titik interior menggunakan skema x for i=2:2:(I-1) for j=2:2:(J-1) U(i,j)=( U(i-1,j-1)+U(i+1,j-1)+U(i-1,j+1)+U(i+1,j+1) )/4 ... -R*f/4; end % for i end% for j %hitung U(i,j) titik interior menggunakan skema + for j=2:1:J-1 switch mod(j,2) case 0, % j genap for i=3:2:I-2 U(i,j)=(dh*(U(i-1,j)+U(i+1,j))+ dk*(U(i,j-1)+U(i,j+1)))/dhk ... -f/dhk; end; % for i otherwise % ganjil for i=2:2:I-1 U(i,j)=(dh*(U(i-1,j)+U(i+1,j))+ dk*(U(i,j-1)+U(i,j+1)))/dhk ... -f/dhk; end; % for i end;%switch mod end; % for j=2 end; % switch level end; % for level %------------end tebakan awal titik interior----------------------%-------------------- ITERASI SOR---------------------------------Ukonv=0; banyakU=I*J; % banyaknya U yang belum diketahui nilainya epsilon=0.01 % nilai toleransi perbedaan nilai U %Mj=(cos(pi/(I+1)) + cos(pi/(J+1)))/2; %w= 2/(1+sqrt(1-Mj.^2); % parameter SOR w=1.2 dh=1/(h*h); dk=1/(k*k); dhk= 2*(dh+dk); iter=0; pause; while Ukonv1&i1&j
130

- f)/dhk; end;%if h==k end;%if i>1&i
% Warning ! % KOORDINAT GAMBAR INI MASIH TERBALIK DAN i,j BELUM DI UBAH X, Y % % --------------------------------------------------------------% plot solusi numerik pada domain solusinumerik= U(1:I,1:J); pcolor(solusinumerik); colorbar vert shading interp colormap(jet) title('Solusi Numerik'); xlabel('nilai i'); ylabel('nilai j') drawnow % end % while. end ini untuk melihat gambar tiap iterasi(tidak diaktifkan) %Diketahui solusi eksaknya: u(x,y)=x.^2 + y.^2 %Hitung error globalnya :error(i,j)=U(i,j)= x.^2 + y.^2 - U(i,j); % variabel error(i,j) menimpa U(i,j) untuk menghemat memori for i=1:I x=(i-1)*h; for j=1:J y=(j-1)*k; U(i,j)=U(i,j)- (x.^2 + y.^2 ); end; end; %Plot errornya figure; mesh(U(1:I,1:J)) shading interp colormap(jet) title('Grafik Nilai error '); xlabel('nilai i'); ylabel('nilai colorbar vert fprintf('\nSelesai\n')

j');zlabel('err or');

131

LAMPIRAN 13 Program untuk Menyelesaikan PD Poisson dan Laplace 2D dengan Tipe Syarat Batas Dirichlet, Neumann, Robin. clc; close all; clear all; %program6.m % Menyelesaikan PD Poisson & laplace % dengan domain segi empat dengan 3 tipe syarat bataa dapat dipilih beripe Dirichlet, Neumman, % Robbin % % j % y Syarat batas 2 % J 1 ---------------- -----# % . | | % . | | % . | | % . | | % .Syarat batas 3 | Uxx + Uyy = 4 | Sayarat batas 4 % . | | % . | | % 3 | | % 2 | | % 1 #-------------- -------+--- x % 0 1 % Syarat batas 1 % 1 2 3 ............... .I i % % %-----------------------------------------------------------------------------------

disp(' BENTUK UMUM : Uxx + Uyy = f '); disp(' -------------------------------'); disp(' Jenis Persamaan Differensial: '); disp(' 1. Laplace ') disp(' 2. Poisson ') PD=input(' * pilih 1-2 : ? '); if PD==1, f=0; else f=input(' end; %if

f = ? ');

disp('UKURAN DOMAIN') disp('--------------') lebar= input(' lebar = ? '); tinggi= input(' tinggi = ? '); I=input('I:inde ks maks.pilahan pd arah sb-x. i=1,2,.,I. J=input('J:inde ks maks.pilahan pd arah sb-y. j=1,2,.,J.

I = ? '); J = ? ');

h=lebar/(I-1); % lebar kisi k=tinggi/(J-1);% tinggi kisi

disp(' ');

132

% s : indeks sisi s=1,2,,3,4 disp('SYARAT BATAS'); disp('------------'); for s=1:4

fprintf('\n Syarat batas sisi ke');fprintf('% 2d ',(s)); fprintf('\n 1. Dirichlet '); fprintf('\n 2. Neumann '); fprintf('\n 3. Campuran (Robyn)\n'); tipe = input(' * pilih 1-3 : '); switch tipe case 1, TipeBC(s)=1; % syarat batas sisi ke s bertipe Dirichlet BC(s)=input(' u = '); case 2, TipeBC(s)=3; % syarat batas bertipe Neummann dianggap alfa(s)=0; % syarat batas bertipe Robyn dgn alfa=0 switch s case 1, fprintf(' du/dy '); case 2, fprintf(' du/dy '); case 3 fprintf(' du/dx '); case 4, fprintf(' du/dx '); otherwise % tidak ada lagi sisi end %switch s BC(s)=input(' = '); %syarat batas sisi ke s bertipe Robyn otherwise, TipeBC(s)=3; switch s case 1, disp(' dU/dy=alfa.u+betta'); case 2, disp(' dU/dy=alfa.u+betta'); case 3 disp(' dU/dx=alfa.u+betta'); case 4, disp(' dU/dx=alfa.u+betta'); otherwise % tidak ada lagi sisi Robyn end %switch s alfa(s)=input(' BC(s)=input('

alfa ='); betta = ');

end % switch tipe end; %for s % inisialisasi U pada batas for i=1:I U(i,1)=0; U(i,J)=0; end; for j=1:J U(1,j)=0; U(I,j)=0; end;

133

% set nilai U pada syarat batas bertipe Dirichlet if TipeBC(1)==1, for i=1:I U(i,1)=BC(1); end; A1=U(1,1); B1=U(I,1); end; % if if TipeBC(2)==1, for i=1:I U(i,J)=BC(2); end; C2=U(1,J); D2=U(I,J); end; % if if TipeBC(3)==1, for j=1:J U(1,j)=BC(3); end; A3=U(1,1); C3=U(1,J); end; % if if TipeBC(4)==1, for j=1:J U(I,j)=BC(4); end; B4=U(I,1); D4=U(I,J); end; % if

%----- nilai U pada titik pojok mengikuti sisi Dirichlet----% pojok A if TipeBC(1)==1 & TipeBC(3)==1 , U(1,1)=(A1+A3)/2; end; if TipeBC(1)==1 & TipeBC(3)==3, U(1,1)=A1; end; if TipeBC(1)==3 & TipeBC(3)==1, U(1,1)=A3; end; % pojok B if TipeBC(1)==1 & TipeBC(4)==1, U(I,1)=(B1+B4)/2; end; if TipeBC(1)==1 & TipeBC(4)==3, U(I,1)=B1; end; if TipeBC(1)==3 & TipeBC(4)==1, U(I,1)=B4; end; % pojok C if TipeBC(2)==1 & TipeBC(3)==1, U(1,J)=(C2+C3)/2; end; if TipeBC(2)==1 & TipeBC(3)==3, U(1,J)=C2; end; if TipeBC(2)==3 & TipeBC(3)==1,

134

U(1,J)=C3; end;

% pojok D if TipeBC(2)==1 & TipeBC(4)==1, U(I,J)=(D2+D4)/2; end; if TipeBC(2)==1 & TipeBC(4)==3, U(I,J)=D2; end; if TipeBC(2)==3 & TipeBC(4)==1, U(I,J)=D4; end; % set tebakan awal U pada sisi dgn tipe Neumann atau Robyn % sisi ke 1 if TipeBC(1)==3 & TipeBC(2)==1, for i=2:I-1 U(i,1)=(BC(2)-tinggi*BC(1)-tinggi)/(1+tinggi*alfa(1)); end; end; % if % sisi ke 2 if TipeBC(1)==1 & TipeBC(2)==3, for i=2:I-1 U(i,J)=BC(1)+tinggi*BC(2)/(1-tinggi*alfa(2)); end; end; % if % sisi ke 3 if TipeBC(3)==3 & TipeBC(4)==1, for j=2:J-1 U(1,j)=BC(4)-lebar*BC(3)/(1+lebar*alfa(3)); end; end; % if % sisi ke 4 if TipeBC(3)==1 & TipeBC(4)==3, for j=2:J-1 U(I,j)=BC(3)+lebar*BC(4)/(1-lebar*alfa(4)); end; end; % if % Tebakan awal pojok yang dibentuk 2 sisi Neumann atu Robyn % pojok A if TipeBC(1)==3 & TipeBC(3)==3, A1=(U(1,2)-k*BC(1))/(1+k*alfa(1)); A3=U(2,1)-h*BC(3)/(1+h*alfa(3)); U(1,1)=(A1+A3)/2; end; % tebakan awal pojok B if TipeBC(1)==3 & TipeBC(4)==3, B1=U(I,2)-k*BC(1)/(1+k*alfa(1)); B4=U(I-1,1)-h*BC(4)/(1+h*alfa(4)); U(I,1)=(B1+B4)/2; end;

% tebakan awal pojok C if TipeBC(2)==3 & TipeBC(3)==3, C2=U(1,J-1)+k*BC(2)/(1-k*alfa(2)); C3=U(2,J)-h*BC(3)/(1+h*alfa(3));

135

U(1,J)=(C2+C3)/2; end; % tebakan awal pojok D if TipeBC(2)==3 & TipeBC(4)==3, D2=U(I,J-1)+k*BC(2)/(1-k*alfa(2)); D4=U(I-1,J)+h*BC(4)/(1-h*alfa(4)); U(I,j)=(D2+D4)/2; end; H=lebar/2; K=tinggi/2; r=sqrt(H*H+K*K); aver = (U(1,1)+U(I,1)+U(1,J)+U(I,J)-r*r*f)/4; %tebakan awal titik interior for i=2:I-1 for j=2:J-1 U(i,j)=0;%aver end end

Ukonv=0; % inisialisasi banyaknya titik konvergen banyakU=I*J; % banyaknya U seluruhnya epsilon=10.^(-6); % nilai toleransi perbedaan nilai U Mj=(cos(pi/(I+1)) + cos(pi/(J+1)))/2; w= 1.5;%2/(1+sqrt(1-Mj.^2)) pause; dh=1/(h*h); dk=1/(k*k); dhk= 2*(dh+dk); iter=0; %-------------------Bagian iterasi-------------------------while Ukonv1&i1&j1&i1&i1&j1&j
%-------skema pojok yang dibentuk 2 sisi Robyn-----------------

136

elseif i==1&j==1 & TipeBC(1)==3 & TipeBC(3)==3, % pojok kiri bawah (pojok A) U(1,1)=(dh*(2*U(2,1)-2*h*BC(3))+dk*(2*U(1,2)-2*k*BC(1))- f)/... (dhk+dh*2*h*alfa(3)+dk*2*k*alfa(1));

elseif i==I&j==1& TipeBC(1)==3 & TipeBC(4)==3, % pojok kanan bawah (pojok B) U(I,1)=(dh*(2*U(I-1,1)+2*h*BC(4))+dk*(2*U(I,2)-2*k*BC(1))- f)/... (dhk-dh*2*h*alfa(4)+dk*2*k*alfa(1)); elseif i==1&j==J & TipeBC(2)==3 & TipeBC(3)==3, % pojok kiri atas (pojok C) U(1,J)=(dh*(2*U(2,J)-2*h*BC(3))+dk*(2*U(1,J-1)+2*k*BC(2))(dhk+dh*2*h*alfa(3)-dk*2*k*alfa(2));

f)/...

elseif i==I&j==J & TipeBC(2)==3 & TipeBC(4)==3, % pojok kanan atas (pojok D) U(I,J)=(dh*(2*U(I-1,1)+2*h*BC(4))+dk*(2*U(I,J-1)+2*k*BC(2))- f)/... (dhk-dh*2*h*alfa(4)-dk*2*k*alfa(2)); else %--sisi Dirichlet tidak punya skema iterasi %--karena nilai U nya telah diketahui end % if

% skema SOR U(i,j)=w*U(i,j)+(1-w)*temp; % kriteria konvergensi beda=abs(U(i,j)-temp); if beda < epsilon, Ukonv=Ukonv+1 ; end; end % for j end% for i fprintf('\nIterasi ke'); fprintf('%5d ',iter); fprintf('banyak U konvergen');fprintf('%5d ',Ukonv); fprintf(' //');fprintf('%5d',banyakU);

% Warning ! % KOORDINAT GAMBAR INI MASIH TERBALIK DAN i,j BELUM DI UBAH X, Y % % ---------------------------------------------------------------

% %

% plot nilai U pada domain arsiran= U(1:I,1:J); pcolor(arsiran) % colorbar vert % shading interp % title('Grafik Nilai U pada domain'); % xlabel('nilai i'); ylabel('nilai j'); %drawnow;

%end % while %figure;

137

% plot 3D solusi numerik % surf(U(1:I,1:J ));shading interp;colormap( jet); % title('Solusi Numeik'); %xlabel('nilai i'); ylabel('nilai j');zlabel('nilai U'); %colorbar vert; %-------menghitung nilai error melalui fungsi Residu(i,j)-----------for i=1:I for j=1:J if

i>1&i1&j
%--- skema Residu sisi Neumann atau Robyn (tak termasuk ujungnya)-----% sisi ke 1 (sisi bawah) elseif i>1&i1&i1&j1&j
elseif i==I&j==1& TipeBC(1)==3 & TipeBC(4)==3, % pojok kanan bawah (pojok B) Residu(I,1)= dh*(2*U(I-1,1)+2*h*BC(4))+dk*(2*U(I,2)-2*k*BC(1))- f-... (dhk-dh*2*h*alfa(4)+dk*2*k*alfa(1))*U(I,1); elseif i==1&j==J & TipeBC(2)==3 & TipeBC(3)==3, % pojok kiri atas (pojok C) Residu(1,J)= dh*(2*U(2,J)-2*h*BC(3))+dk*(2*U(1,J-1)+2*k*BC(2))- f-... (dhk+dh*2*h*alfa(3)-dk*2*k*alfa(2))*U(1,J); elseif i==I&j==J & TipeBC(2)==3 & TipeBC(4)==3, % pojok kanan atas (pojok D) Residu(I,J)= dh*(2*U(I-1,1)+2*h*BC(4))+dk*(2*U(I,J-1)+2*k*BC(2))- f-... (dhk-dh*2*h*alfa(4)-dk*2*k*alfa(2))*U(I,J); else %--sisi Dirichlet nilainya telah ada jadi tak punya error jadi Residunya=0 Residu(i,j)=0; end % if i>1&i1&j
138

Rbesar=max(abs(Residu(i,j))); end % for j end% for i Rbesarmutlak(iter)=Rbesar;

% Plot 3D Nilai fungsi Residu(i,j) %figure; %mesh(Residu(1:I,1:J));shading interp;colormap(jet); %title('Grafik error'); %xlabel(' i '); ylabel(' j ');zlabel(' R(i,j)'); %colorbar vert;

% plot nomor iterasi VS maks|R(i,j,iter)| GBR=plot(Rbesarmutlak(1:iter)); set(GBR,'Color','blue') title('Grafik maks|R| tiap iterasi'); xlabel('Nomor iterasi (n)|'); ylabel('Nilai Maks|R| '); drawnow; pause; if iter==100000, break; end %if end % while--------------------End iterasi----------------------fprintf('\nSelesai\n');

139

Lampiran 14 Program untuk Menghitung Spectral Radius Matrik Iterasi untuk Contoh 4.3 >> % Menghitung spectral radius matrik iterasi untuk Contoh 4.3 % Untuk contoh lainnya menyesuiakan. close all; clear all; clc; %Grid Io X Jo = 33 x 33 Io = 33; Jo = 33;

n=( Io -1); % orde matrik T adalah nxn % orde matrik A adalah mxm dengan m=n^2 m= n*n; C=[4 -1 zeros(1,n-2)]; T=toeplitz(C); I=eye(size(T)); O=zeros(size(T));

A=[T -I O O O O O O O O O O O O O O O O O O O O O O O O O O O O O O;... -I T -I O O O O O O O O O O O O O O O O O O O O O O O O O O O O O;... O -I T -I O O O O O O O O O O O O O O O O O O O O O O O O O O O O;... O O -I T -I O O O O O O O O O O O O O O O O O O O O O O O O O O O;... O O O -I T -I O O O O O O O O O O O O O O O O O O O O O O O O O O;... O O O O -I T -I O O O O O O O O O O O O O O O O O O O O O O O O O;... O O O O O -I T -I O O O O O O O O O O O O O O O O O O O O O O O O;... O O O O O O -I T -I O O O O O O O O O O O O O O O O O O O O O O O;... O O O O O O O -I T -I O O O O O O O O O O O O O O O O O O O O O O;... O O O O O O O O -I T -I O O O O O O O O O O O O O O O O O O O O O;... O O O O O O O O O -I T -I O O O O O O O O O O O O O O O O O O O O;... O O O O O O O O O O -I T -I O O O O O O O O O O O O O O O O O O O;... O O O O O O O O O O O -I T -I O O O O O O O O O O O O O O O O O O;... O O O O O O O O O O O O -I T -I O O O O O O O O O O O O O O O O O;... O O O O O O O O O O O O O -I T -I O O O O O O O O O O O O O O O O;... O O O O O O O O O O O O O O -I T -I O O O O O O O O O O O O O O O;... O O O O O O O O O O O O O O O -I T -I O O O O O O O O O O O O O O;... O O O O O O O O O O O O O O O O -I T -I O O O O O O O O O O O O O;... O O O O O O O O O O O O O O O O O -I T -I O O O O O O O O O O O O;... O O O O O O O O O O O O O O O O O O -I T -I O O O O O O O O O O O;... O O O O O O O O O O O O O O O O O O O -I T -I O O O O O O O O O O;... O O O O O O O O O O O O O O O O O O O O -I T -I O O O O O O O O O;... O O O O O O O O O O O O O O O O O O O O O -I T -I O O O O O O O O;... O O O O O O O O O O O O O O O O O O O O O O -I T -I O O O O O O O;... O O O O O O O O O O O O O O O O O O O O O O O -I T -I O O O O O O;... O O O O O O O O O O O O O O O O O O O O O O O O -I T -I O O O O O;... O O O O O O O O O O O O O O O O O O O O O O O O O -I T -I O O O O;... O O O O O O O O O O O O O O O O O O O O O O O O O O -I T -I O O O;... O O O O O O O O O O O O O O O O O O O O O O O O O O O -I T -I O O;... O O O O O O O O O O O O O O O O O O O O O O O O O O O O -I T -I O;... O O O O O O O O O O O O O O O O O O O O O O O O O O O O O -I T -I;... O O O O O O O O O O O O O O O O O O O O O O O O O O O O O O -2*I T];

140

for s=1:n-1 t=s*n; A(t,t-1)=-2; A(t,t+1)=0; A(t+1,t)=0; end; %for s A(m,m-1)=-2;

clear C T

O;

disp(' Ukuran Matrik A'); disp(size(A)) disp(' A = '); disp(' disp(A(1:33,1:33)) disp(' '); disp(' disp(' '); disp('

33 baris dan 33 kolom pertama');

...') 33 baris dan 33 kolom terakhir');

disp(A(n*n-32:n*n,n*n-32:n*n))

% A = D-E-F D=diag(diag(A)); disp(' Ukuran matrik D '); disp(size(D)) disp('D = '); disp(' 7 baris dan 7 kolom pertama'); disp(D(1:7,1:7)) disp(' ...') disp(' 7 baris dan 7 kolom terakhir'); disp(D(n*n-6:n*n,n*n-6:n*n)) E=-tril(A,-1); disp(' Matrik E '); disp(size(E)) disp('E = '); disp(' 7 baris dan 7 kolom pertama'); disp(E(1:7,1:7)) disp(' ...') disp(' 7 baris dan 7 kolom terakhir'); disp(E(n*n-6:n*n,n*n-6:n*n)) F=-triu(A,+1); disp(' Ukuran matrik F '); disp(size(F)) disp('F = '); disp(' 7 baris dan 7 kolom pertama'); disp(F(1:7,1:7)) disp(' ...') disp(' 7 baris dan 7 kolom terakhir'); disp(F(n*n-6:n*n,n*n-6:n*n))

disp(' Parameter SOR : w=1.2 ');w=1.2; M=inv(D-w*E)*((1-w)*D+w*F); disp(' ukuran matrik iterasi :M '); disp(' Matrik iterasi M=inv(D-w*E)*((1-w)*D+w*F)'); disp('M = ');

141

disp(' 7 baris dan 7 kolom pertama'); disp(M(1:7,1:7)) disp(' ...') disp(' 7 baris dan 7 kolom terakhir'); disp(M(n*n-6:n*n,n*n-6:n*n)) pause; M=M(1:n*n,1:n*n); disp(' '); disp(' KOMPONEN MATRIK ITERASI M '); disp('-------------------------------- ') disp(' Baris Kolom '); disp(' no i j M(i,j) '); disp('--------------------------------'); nomor=0; for i=1:n for j=1:n; nomor=nomor+1; fprintf('\n |'); fprintf('%5d ',nomor);fprintf('|'); fprintf('%5d ',i);fprintf('|'); fprintf('%5d ',j);fprintf('|'); fprintf('%8.5f ',M(i,j));fprintf('|'); end; end; disp('... ');

nomor=(n*n-7).^2; for i=n*n-7:n*n for j=n*n-7:n*n nomor=nomor+1; if mod(nomor,500)==0,pause;end; fprintf('\n |'); fprintf('%5d ',nomor);fprintf('|'); fprintf('%5d ',i);fprintf('|'); fprintf('%5d ',j);fprintf('|'); fprintf('%8.5f ',M(i,j));fprintf('|'); end; end; disp(' '); disp('---------------------------------------'); disp(' Spectral radius matrik M=max(abs(eig(M))') spectral_radius=max(abs(eig(M)))

142

DAFTAR PUSTAKA

[1]

2

Abayomi, K., 2003, Iterative Methods for the Laplace Equation in R , http://www.columbia.eduww.columbia.edu/~gb2030/courses/E6302/ Presentations/numanalysisproject.pdf

[2]

Bailey,W., 2003, The SOR algorithm & its Application to Numerical Solution of Eliptic Partial Differential Equation , Dublin Institute of

Technology, Ireland.

[3]

Bassaruddin, T., 1994, Metode Beda Hingga untuk Persamaan Differensial, Elex Media Komputindo, Jakarta.

[4]

Constantinides,A, 1987, Applied Numerical Methods with Personal Computer , Graw-Hill,Mc-Inc, New York.

[5]

Luenberger,D.G., 1969, Optimization By Vector Space Methods , John Wiley & Sons, inc., New York.

[6]

Marek,I., Szyld, D,B., 2005, Algebraic Schwarz Methods For The Numerical Solution Of Markov Chains , Departement of mathematics,

Temple University, Philadelphia, Pennysilvania 19122-6094 USA. http://www.math.temple.edu/~szyld/reports/Schw_mark_rep.pdf .

Penyelesaian numerik Persamaan differensial Poisson 2D, dengan syarat batas Dirichlet, Neumann, Robin dilengkapi Analisis Real dan Program Matlab

Recommend Documents