METODE NUMERIS UNTUK MENYELESAIKAN MODEL PERGERAKAN LAPISAN FLUIDA YANG MELIBATKAN MINYAK DAN AIR

PLAGIAT MERUPAKAN TINDAKAN TIDAK TERPUJI

METODE NUMERIS UNTUK MENYELESAIKAN MODEL PERGERAKAN LAPISAN FLUIDA YANG MELIBATKAN MINYAK DAN AIR

Skripsi Diajukan untuk Memenuhi Salah Satu Syarat Memperoleh Gelar Sarjana Sains Program Studi Matematika

Oleh : Friska Dwi Mesra Mangadil NIM: 133114007

PROGRAM STUDI MATEMATIKA, JURUSAN MATEMATIKA FAKULTAS SAINS DAN TEKNOLOGI UNIVERSITAS SANATA DHARMA YOGYAKARTA 2017 i


NUMERICAL METHOD FOR SOLVING A MOTION MODEL OF FLUID LAYERS INVOLVING OIL AND WATER Thesis Presented as a Partial Fulfillment of the Requirements to Obtain the Degree ofSarjana Sains in Mathematics

By : Friska Dwi Mesra Mangadil Student Number: 133114007

MATHEMATICS STUDY PROGRAM, DEPARTMENT OF MATHEMATICS FACULTY OF SCIENCE AND TECHNOLOGY SANATA DHARMA UNIVERSITY YOGYAKARTA 2017

ii


iii


iv


v


MOTTO

“Segala perkara dapat kutanggung didalam Dia yang memberi kekuatan kepadaku” (Filipi 4:13) “Karena masa depan sungguh ada, dan harapanmu tidak akan hilang” (Amsal 23:18) “Untuk segala sesuatu ada waktunya” (Pengkhotbah 3:1-15)

vi


HALAMAN PERSEMBAHAN

Karya ini kupersembahkan untuk: Tuhan Yesus Kristus yang senantiasa menyertai, mendengarkan, dan selalu memberi perlindungan. Bapak Mesak Luas, Ibu A. Irawati Pongkapadang, Laurance Feien Eka Prakasa Mangadil, S.E., Falerio Ishak Renfeika Mangadil, dan Akhmalia Fiabel Hawari Mangadil yang selalu mendukung, dan memberi keceriaan ketika rindu menghampiri.

vii


ABSTRAK

Pergerakan lapisan fluida merupakan salah satu masalah yang biasanya muncul pada bidang perminyakan. Pergerakan lapisan fluida dapat diselesaikan dan disimulasikan dengan menggunakan banyak metode dan aplikasi komputer. Dalam skripsi ini, dibahas mengenai penyelesaian masalah pergerakan lapisan fluida yang melibatkan minyak dan air. Masalah pergerakan lapisan fluida diselesaikan dengan menggunakan dua metode yaitu metode volume hingga LaxFriedrichs dan metode beda hingga. Metode volume hingga bekerja dengan cara membagi domain ruang menjadi beberapa bagian kemudian dihitung rata-rata kuantitas untuk masing-masing bagian. Metode beda hingga bekerja dengan menghampiri solusi masalah secara titik demi titik. Pengujian dilakukan menggunakan simulasi numeris. Analisis hasil simulasi dilakukan dengan membandingkan hasil solusi numeris dengan solusi eksak, untuk kasus yang mempunyai solusi eksak.

Kata kunci: lapisan fluida, metode beda hingga, metode volume hingga.

viii


ABSTRACT

The motion of fluid layers is one of problems that usually happened in petroleum engineering. The motion of fluid layer can be solved and simulated using many methods and computer application. In this undergraduate thesis, the solution to the problem of motion of fluid layersinvolving oil and water will be discussed. The problem of motion of fluid layers can be solved using two methods: Lax-Friedrichs finite volume method and finite difference method. The finite volume method works by dividing the spatial domain into a finite number of cells, then calculating the average quantity of each cell. The finite difference method works by approaching the solution to the problem point by point. Test cases were done using numerical simulations. Simulation result analysis was conducted by comparing numerical solutions with the analytical ones, for cases having analytical solutions.

Keywords: fluid layers, finite difference method, finite volume method

ix


KATA PENGANTAR

Puji dan syukur penulis panjatkan kepada Tuhan Yesus Kristus yang telah mencurahkan rahmat dan Roh KudusNya sehingga penulis dapat mengerjakan dan menyelesaikan skripsi ini dengan baik. Skripsi ini dibuat dengan tujuan memenuhi salah satu syarat untuk memperoleh gelar Sarjana Sains pada Program Studi Matematika, Fakultas Sains dan Teknologi, Universitas Sanata Dharma. Penulis menyadari bahwa penulis melibatkan banyak pihak untuk membantu dalam menghadapi berbagai macam tantangan, kesulitan, dan hambatan. Oleh karena itu, pada kesempatan ini penulis mengucapkan terima kasih kepada: 1. Bapak Sudi Mungkasi, S.Si., M.Math.Sc., Ph.D., selakuDekan Fakultas Sains dan Teknologi sekaligus sebagaiDosen Pembimbing Skripsi. 2. Bapak Hartono, S.Si., M.Sc., Ph.D., selaku Kaprodi Matematika. 3. Ibu M. V. Any Herawati, S.Si., M.Si., selaku Dosen Pembimbing Akademik. 4. Romo Prof. Dr. Frans Susilo, S.J., Bapak Ir. Ig. Aris Dwiatmoko, M.Sc., Bapak Dr. rer. nat. Herry P. Suryawan, S.Si., M.Si., dan Ibu Lusia Krismiyati Budiasih, S.Si., M.Si. selaku dosen-dosen Prodi Matematika yang telah memberikan banyak pengetahuan kepada penulis selama proses perkuliahan. 5. Bapak/Ibu dosen/tenaga kependidikan Fakultas Sains dan Teknologi yang telah berdinamika bersama selama penulis berkuliah. x


6. Kedua orang tua, kakak, dan dua adik yang telah membantu dan mendukung penulis selama proses pengerjaan skripsi. 7. Teman-teman Matematika 2013: Inge, Yui, Sorta, Melisa, Agung, Laras, Ambar, Yuni, Rey, Dion, Wahyu, Indra, Bintang, Tia, Lya, Andre, Sisca, Natali, Yola, Sari, Dita, dan Kristo yang selalu memotivasi, memberi masukan, dan masih banyak yang tidak bisa disebutkan satu persatu. Terima kasih atas kebersamaan dan kekompakan ini. 8. Kakak, teman-teman dan adik-adik: Mbak Tiwi, Ignatia, Gege, Nando dan Edo, terimakasih untuk semangat dan dukungannya selama penulis berkuliah dan menulis skripsi ini. 9. Kak Mike yang memberi dukungan dengan membantu memperbaiki penulisan bahasa Inggris penulis yang tidak beraturan dan memberi tantangan kepada penulis untuk segera mungkin menyelesaikan skripsi ini. 10. Kak Sri dan Cleo yang memberi dukungan dengan hampir setiap bulan mengingatkan penulis agar mengerjakan skripsi dengan semangat dari kota dan negeri seberang. 11. Mas Susilo yang memberi dukungan dengan mempersilakan penulis mengerjakan skripsi di laboratorium berhari-hari. 12. Pralana Anggi yang selalu siap membantu apabila laptop penulis mengalami gangguan. 13. Pemuda-pemudi GKN Gloria yang memberi dukungan semangat dan doa kepada penulis dalam menyelesaikan skripsi ini.

xi


xii


xiii


DAFTAR ISI

HALAMAN JUDUL...........................................................................................i HALAMAN PERSETUJUAN PEMBIMBING .................................................iii HALAMAN PENGESAHAN .............................................................................iv PERNYATAAN KEASLIAN KARYA .............................................................v HALAMAN PERSEMBAHAN .........................................................................vii ABSTRAK ..........................................................................................................viii ABSTRACT ........................................................................................................ix KATA PENGANTAR ........................................................................................x PERNYATAAN PERSETUJUAN PUBLIKASI ...............................................xiii DAFTAR ISI .......................................................................................................xiv DAFTAR GAMBAR ..........................................................................................xvi BAB I PENDAHULUAN ...................................................................................1 A. Latar Belakang....................................................................................1 B. Rumusan Masalah ..............................................................................4 C. Batasan Masalah .................................................................................5 D. Tujuan Penulisan ................................................................................5 E. Manfaat Penelitian ..............................................................................6 F. Metode Penelitian ...............................................................................6 G. Sistematika Penulisan .........................................................................6

xiv


BAB II PERSAMAAN DIFERENSIAL ............................................................9 A. Integral ................................................................................................9 B. Klasifikasi Persamaan Diferensial ......................................................12 C. Nilai Eigen dan Vektor Eigen .............................................................15 D. Klasifikasi Persamaan Diferensial Parsial Orde Dua .........................16 E. Penurunan Numeris ............................................................................18 F. Matriks Tridiagonal ............................................................................25 BAB III MODEL PERGERAKAN LAPISAN FLUIDA ...................................27 A. Penurunan Persamaan Gerak Lapisan Fluida Satu Dimensi ..............27 B. Masalah Pergerakan Fluida ................................................................31 C. Solusi Analitis Masalah Pergerakan Lapisan Fluida ..........................33 D. Metode Volume Hingga Lax-Friedrichs.............................................35 E. Metode Beda Hingga untuk Model Pergerakan Lapisan Fluida.........42 BAB IV ANALISIS HASIL SIMULASI ...........................................................50 A. Metode Volume Hingga Lax-Friedrichs.............................................50 B. Metode Beda Hingga ..........................................................................52 BAB V PENUTUP ..............................................................................................55 A. Kesimpulan .........................................................................................55 B. Saran ...................................................................................................55 DAFTAR PUSTAKA .........................................................................................57 LAMPIRAN ........................................................................................................58

xv


DAFTAR GAMBAR

Gambar 1.1 Dua plat rata dengan jarak 10 cm berisi lapisan fluida. Plat atas ditarik ke kanan dengan kecepatan konstan ..................................... 2 Gambar 2.1 Ilustrasi fungsi satu variabel ............................................................ 10 Gambar 2.2 Ilustrasi geometri vektor eigen ........................................................ 16 Gambar 2.3 a. Hampiran beda majub. Hampiran beda mundurc. Hampiran beda pusat ................................................................................................. 22 Gambar 3.1 Kawat satu dimensi dengan energi panas yang mengalir masuk dan keluar ............................................................................................... 28 Gambar 3.2 Ilustrasi diskretisasi domain ruang .................................................. 36 Gambar 3.3 Hasil simulasi masalah pergerakan lapisan fluida dengan metode volume hingga saat 𝑡 = 50 .............................................................. 42 Gambar 3.4Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 2 dengan metode beda hingga ............................................................. 46 Gambar 3.5Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 1 dengan metode beda hingga ................................................ 47 Gambar 3.6Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 0.5 dengan metode beda hingga............................................. 47 Gambar 3.7Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 0.25 dengan metode beda hingga .......................................... 48 Gambar 3.8Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 0.125dengan metode beda hingga ........................................... 48

xvi


Gambar 3.9Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 0.0625dengan metode beda hingga ......................................... 49

Gambar 4.1 Ilustrasi geometri galat metode volume hingga Lax-Friedrichs ...... 52 Gambar 4.2 Ilustrasi geometri galat metode beda hingga ................................... 54

xvii


BAB I PENDAHULUAN

A. Latar Belakang Masalah Persamaan diferensial adalah persamaan yang menyatakan hubungan suatu fungsi dengan turunan-turunannya. Persamaan diferensial biasanya diklasifikasikan menjadi dua, yaitu pesamaan diferensial biasa dan persamaan diferensial parsial. Persamaan diferensial biasa merupakan persamaan diferensial yang memuat satu variabel bebas. Persamaan diferensial parsial sebenarnya hampir sama dengan persamaan diferensial biasa, perbedaannya terletak pada banyaknya variabel bebas. Pada persamaan diferensial parsial terdapat lebih dari satu variabel bebas, sehingga terdapat turunan parsial. Fluida adalah zat yang dapat mengalir atau biasa disebut zat alir. Pada prinsipnya, fluida adalah semua jenis zat cair dan zat gas. Fluida biasanya banyak digunakan dalam kehidupan sehari-hari, contohnya adalah minyak dan air. Minyak adalah zat cair yang mengandung lemak dan memiliki suatu kekentalan. Jika dilihat dari asalnya, minyak dapat dikelompokkan menjadi minyak nabati dan hewani. Minyak yang telah diolah banyak digunakan oleh masyarakat, seperti untuk memasak.

1

PLAGIAT MERUPAKAN TINDAKAN TIDAK TERPUJI 2

Minyak dan air memiliki massa jenis yang berbeda.Jika keduanya dimasukkan ke dalam suatu wadah, keduanya tidak akan tercampur menjadi satu cairan, melainkan akan terpisah dengan air berada di bawah minyak. Hal itu disebabkan karena massa jenis air lebih besar daripada massa jenis minyak.

𝑥 = 10

9 8 7

6

5 4 3

2 1 𝑥=0

Gambar 1.1. Dua plat rata dengan jarak 10 cm berisi lapisan fluida. Plat atas ditarik ke kanan dengan kecepatan konstan. Dalam skripsi ini akan dibahas tentang pengaruh pergerakan minyak terhadap pergerakan air. Minyak dan air dalam kasus ini diletakkan di antara dua plat rata sehingga air akan berada di bawah minyak seperti terlihat pada Gambar 1.1. Pengaruh pergerakan minyak terhadap pergerakan air dilihat dengan menggunakan kajian. Kajian dilakukan terhadap minyak dan air yang berada di antara dua plat horizontal. Ketika plat yang berada di atas minyak ditarik dengan kecepatan konstan, maka akan terbentuk gelombang di permukaan minyak. Gelombang di permukaan minyak inilah yang akan mempengaruhi pergerakan air. Pertanyaan yang timbul adalah seberapa besarkah pengaruhnya?


Pengaruh

pergerakan

minyak

terhadap

pergerakan

air

akan

diprediksikan dengan mempertimbangkan jarak kedua plat dan waktu yang berbeda. Maksudnya ketika plat atas ditarik dengan kecepatan konstan, akan dilihat seberapa besar pengaruh yang muncul terhadap kedua cairan ini pada waktu tertentu. Misalnya plat atas ditarik dengan kecepatan konstan dan dalam waktu 1 detik, akan dilihat berapa besar pengaruh pergerakan minyak terhadap air dalam waktu tersebut. Pengaruh pergerakan minyak terhadap pergerakan air di sini akan dicari dengan menggunakan metode volume hingga dan metode beda hingga. Metode volume hingga adalah salah satu metode penyelesaian persamaan diferensial parsial. Metode volume hingga bekerja dengan mendiskretkan domain ruang ke dalam interval, kemudian dihitung rata-rata kuantitas untuk masing-masing interval. Perhitungan rata-rata ini menghasilkan fluks, maka dalam metode volume hingga selain mendiskretkan ruang ke dalam interval dan menghitung rata-rata tiap interval, harus dihitung pula fluks agar hasil yang didapat stabil. Metode beda hingga merupakan suatu metode yang menghampiri penyelesaian model matematika titik demi titik. Metode ini menggunakan pendekatan ekspansi Taylor di suatu titik acuan. Metode beda hingga unggul dalam kemudahan komputasi. Pada skripsi ini akan dibahas mengenai seberapa besar pengaruh pergerakan minyak terhadap pergerakan air dengan metode volume hingga dan beda hingga.


Pergerakan minyak dan air yang akan diselesaikan di sini menggunakan persamaan gerak fluida, yaitu: 𝜕𝑣minyak 𝜕 2 𝑣minyak = 𝜇minyak 𝜕𝑡 𝜕𝑥 2

( (1)

dan 𝜕𝑣air 𝜕 2 𝑣air = 𝜇air 𝜕𝑡 𝜕𝑥 2

( (2)

dengan 𝑣minyak (𝑥, 𝑡) adalah kecepatan minyak, 𝑣air (𝑥, 𝑡) adalah kecepatan air, 𝜇minyak (𝑥, 𝑡) adalah kekentalan minyak, dan 𝜇air (𝑥, 𝑡) adalah kekentalan air. Di sini, variabel bebas x mewakili domain ruang dan variabel t melambangkan nilai waktu.Hubungan di titik perbatasan antara minyak dan air diberikan oleh: 𝑣minyak = 𝑣air

( (3)

dan

𝜇minyak

𝜕𝑣minyak 𝜕𝑣air = 𝜇air , 𝜕𝑥 𝜕𝑥

( (4)


B. Rumusan Masalah Berdasarkan latar belakang yang telah dipaparkan, penulis mengadakan penelitian terhadap masalah-masalah berikut: 1.

Bagaimana memperoleh persamaan pergerakan lapisan fluida yang dipengaruhi plat yang ditarik dengan kecepatan konstan?

2.

Bagaimana menyelesaikan persamaan pergerakan lapisan fluida dengan metode volume hingga dan beda hingga?

C. Pembatasan Masalah Penulis akan membatasi penulisan agar menjadi lebih terarah dan tidak menyimpang dari masalah yang akan dibahas, yaitu: 1.

Persamaan gelombang yang akan dibahas adalah gelombang lapisan minyak dan air.

2.

Persamaan diselesaikan dengan metode analitis,metode volume hingga dan metode beda hingga.

3.

Masalah pergerakan lapisan fluida yang diselesaikan merupakan masalah pergerakan lapisan fluida satu dimensi.

D. Tujuan Penulisan Tujuan yang ingin dicapai penulis, selain untuk memenuhi syarat tugas akhir dalam Program Studi Matematika Universitas Sanata Dharma, juga untuk:


1.

Mencari seberapa besar pengaruh pergerakan minyak terhadap pergerakan air dengan menggunakan metode volume hingga dan beda hingga.

2.

Memperluas wawasan pembaca tentang aplikasi matematika dalam pengaruh pergerakan suatu cairan terhadap cairan lain yang memiliki massa jenis yang berbeda.

E. Manfaat Penulisan Manfaat dari penulisan skripsi ini adalah sebagai berikut: 1. Penulis memperoleh pengetahuan baru selama mengerjakan tulisan ini. 2. Pembaca mendapat gambaran tentang aplikasi matematika dalam pengaruh pergerakan suatu cairan terhadap cairan lain yang memiliki massa jenis berbeda. 3. Skripsi ini dapat dijadikan referensi bagi peneliti lain.

F. Metode Penelitian Metode penelitian yang digunakan oleh penulis dalam menulis skripsi adalah studi pustaka, yaitu dengan mempelajari buku maupun jurnal yang berkaitan dengan metode volume hingga dan beda hingga khususnya dalam mencari seberapa besar pengaruh pergerakan minyak terhadap pergerakan air; komputasi yang digunakan pada penelitian ini adalah komputasi numeris


khususnya untuk metode volume hingga dan metode beda hingga; selain itu juga akan dilakukan simulasi dengan komputer.

G. Sistematika Penulisan Skripsi ini ditulis menggunakan sistematika berikut: BAB I PENDAHULUAN A. Latar Belakang B. Rumusan Masalah C. Batasan Masalah D. Tujuan Penulisan E. Manfaat Penulisan F. Metode Penulisan G. Sistematika Penulisan BAB II PERSAMAAN DIFERENSIAL A. Integral B. Klasifikasi Persamaan Diferensial C. Nilai Eigen dan Vektor Eigen D. Klasifikasi Persamaan Diferensial Parsial Orde Dua E. Penurunan Numeris F. Matriks Tridiagonal BAB III MODEL PERGERAKAN LAPISAN FLUIDA


A. Penurunan Persamaan Gerak Lapisan Fluida Satu Dimensi B. Masalah Pergerakan Lapisan Fluida C. Solusi Analitis Masalah Pergerakan Lapisan Fluida D. Metode Volume Hingga Lax-Friedrichs E. Metode Beda Hingga BAB IV ANALISIS HASIL SIMULASI A. Metode Volume Hingga Lax-Friedrichs B. Metode Beda Hingga BAB V PENUTUP A. Kesimpulan B. Saran DAFTAR PUSTAKA


BAB II PERSAMAAN DIFERENSIAL Dalam bab ini akan dibahas mengenai dasar teori dari skripsi ini. Dasar teori dari skripsi ini meliputi integral, klasifikasi persamaan diferensial, nilai eigen dan vektor eigen, matriks tridiagonal, klasifikasi persamaan diferensial parsial orde dua, dan penurunan numeris. A.

Integral Pada bagian ini akan dibahas mengenai integral yang meliputi definisi dan

contoh dari integral tentu dan teorema fundamental kalkulus. Definisi 2.1 Suatu fungsi 𝐹 disebut anti turunan dari 𝑓 pada interval 𝐼, jika 𝐹 ′ (𝑥) = 𝑓(𝑥) untuk setiap 𝑥 dalam interval 𝐼. Contoh 2.1 Carilah suatu anti turunan dari 𝑓(𝑥) = 𝑥 2 . Penyelesaian: Fungsi 𝐹(𝑥) = 𝑥 3 bukanlah anti turunannya, karena turunan dari 𝑥 3 adalah 3𝑥 2 . 1

1

Akan tetapi hal ini menyarankan 𝐹(𝑥) = 3 𝑥 3 , yang memenuhi 𝐹 ′ (𝑥) = 3 3𝑥 2 = 1

𝑥 2 . Dengan demikian, suatu anti turunan dari 𝑓 adalah 3 𝑥 3 . 9


Anti turunan dinotasikan dengan ∫ … 𝑑𝑥. Notasi tersebut menunjukkan anti turunan terhadap 𝑥. Anti turunan biasanya disebut integral tak tentu. 1. Integral Tentu

𝑦

𝑦 = 𝑓(𝑥)

𝑎

𝑏

𝑥

Perhatikan Gambar 2.1 berikut ini: Gambar 2.1. Ilustrasi fungsi satu variabel.


Untuk menghitung luas di bawah kurva 𝑦 = 𝑓(𝑥), dapat dilakukan dengan aproksimasi, yaitu dengan membagi interval [𝑎, 𝑏] oleh partisi 𝑝 = {𝑥0 , 𝑥1 , … , 𝑥𝑛 } ke dalam n subinterval yaitu [𝑥0 , 𝑥1 ], [𝑥1 , 𝑥2 ], …[𝑥𝑛−1 , 𝑥𝑛 ]. Panjang subinterval ke-𝑖 ditulis dengan ∆𝑥𝑖 = 𝑥𝑖 − 𝑥𝑖−1. Selanjutnya dipilih sebarang 𝑥𝑖∗ dari [𝑥0 , 𝑥1 ], [𝑥1 , 𝑥2 ], …[𝑥𝑛−1 , 𝑥𝑛 ] dengan 𝑖 = 1, 2, … , 𝑛. Total luas di bawah kurva dapat dihitung dengan 𝑓(𝑥1∗ )∆𝑥1 + 𝑓(𝑥2∗ )∆𝑥2 + ⋯ + 𝑓(𝑥𝑛∗ )∆𝑥𝑛 = ∑𝑛𝑖=1 𝑓(𝑥𝑖∗ )∆𝑥𝑖 yang disebut jumlahan Riemann fungsi 𝑓 pada interval [𝑎, 𝑏], sebagai pendekatan luas daerah di bawah kurva 𝑦 = 𝑓(𝑥) dan diatas sumbu 𝑥. Semakin banyak subinterval seragam yang digunakan artinya ∆𝑥𝑖 → 0, maka semakin baik pula aproksimasi luasan tersebut dan semakin dekat dengan luasan yang sebenarnya. Dengan demikian, luas daerah = lim ∑𝑛𝑖=1 𝑓(𝑥𝑖∗ )∆𝑥𝑖 . ∆𝑥𝑖 →0

Definisi 2.2 Misalkan𝑓 suatu fungsi yang didefinisikan pada selang tertutup [𝑎, 𝑏]. Jika 𝑛

lim ∑ 𝑓(𝑥𝑖∗ )∆𝑥𝑖

∆𝑥𝑖 →0

𝑖=1

ada, maka nilai limit tersebut dinamakan integral tentu 𝑓 dari 𝑎 ke 𝑏 dan ditulis 𝑏

sebagai ∫𝑎 𝑓(𝑥)𝑑𝑥 = lim ∑𝑛𝑖=1 𝑓(𝑥𝑖∗ )∆𝑥𝑖 . ∆𝑥𝑖 →0


2. Teorema Fundamental Kalkulus Pada bagian ini hanya akan diberikan teorema fundamental kalkulus, tidak dibahas mengenai pembuktiannya. Teorema 2.1 (Teorema Nilai Rata-Rata) Jika 𝑓 fungsi kontinu pada [𝑎, 𝑏], maka terdapat 𝑐 ∈ [𝑎, 𝑏], sehingga berlaku 𝑏 1 𝑓(𝑐) = ∫ 𝑓(𝑥)𝑑𝑥. 𝑏−𝑎 𝑎

Teorema 2.2 (Teorema Fundamental Kalkulus I) 𝑥

Jika 𝑓 fungsi kontinu pada [𝑎, 𝑏], maka 𝐹(𝑥) = ∫𝑎 𝑓(𝑡)𝑑𝑡 kontinu pada [𝑎, 𝑏] dan terdiferensial pada (𝑎, 𝑏) dan berlaku 𝑑

𝑥

𝐹 ′ (𝑥) = 𝑑𝑥 ∫𝑎 𝑓(𝑡)𝑑𝑡 = 𝑓(𝑥). Teorema 2.3 (Teorema Fundamental Kalkulus II) Jika 𝑓 fungsi kontinu pada setiap titik dalam [𝑎, 𝑏] dan 𝐹 adalah antiturunan dari 𝑓 pada [𝑎, 𝑏], maka 𝑏

∫𝑎 𝑓(𝑥)𝑑𝑥 = 𝐹(𝑏) − 𝐹(𝑎).


Bukti dari ketiga teorema yang disebut di atas dapat dilihat pada buku karangan Thomas (2010). B. Klasifikasi Persamaan Diferensial Berikut ini akan dibahas mengenai klasifikasi persamaan diferensial. Klasifikasi tersebut meliputi definisi dan contoh persamaan diferensial, persamaan diferensial biasa, persamaan diferensial parsial, dan orde persamaan diferensial. Definisi 2.3 Persamaan diferensial adalah persamaan yang melibatkan variabelvariabel tak bebas dari fungsi yang tidak diketahui dan turunan terhadap variabelvariabel bebas dari fungsi tersebut. Contoh 2.3 Persamaan-persamaan di bawah ini merupakan contoh persamaan diferensial: 𝑑𝑦 = −𝑦 2 , 𝑑𝑡

(2.1)

𝑑𝑦 𝑑2𝑦 + 2 2 = 0, 𝑑𝑥 𝑑𝑥

(2.2)

𝜕𝑢 𝜕 2 𝑢 𝜕 2 𝑢 = + , 𝜕𝑡 𝜕𝑥 2 𝜕𝑦 2

(2.3)

𝜕 2𝑢 𝜕 2𝑢 𝜕 2𝑢 + + = 0. 𝜕𝑥 2 𝜕𝑦 2 𝜕𝑧 2

(2.4)


Definisi 2.4 Persamaan diferensial biasa adalah persamaan diferensial yang melibatkan turunan biasa atas satu atau lebih variabel tak bebas terhadap satu variabel bebas. Contoh 2.4 Persamaan (2.1) dan (2.2) adalah contoh persamaan diferensial biasa. Pada persamaan (2.1) variabel 𝑦 adalah variabel terikat atau tak bebas dan variabel 𝑡 adalah variabel bebas. Pada persamaan (2.2) variabel 𝑦 adalah varabel tak bebas dan variabel 𝑥 adalah variabel bebas.

Definisi 2.5 Persamaan diferensial parsial adalah persamaan diferensial yang melibatkan turunan parsial atas satu atau lebih variabel tak bebas terhadapvariabel bebas, dengan catatan bahwa banyaknya variabel bebas dalam persamaan tersebut adalah lebih dari satu. Contoh 2.5 Persamaan (2.3) dan (2.4) merupakan contoh persamaan diferensial parsial. Pada persamaan (2.3) variabel 𝑥, 𝑦, dan 𝑡 merupakan variabel bebas dan variabel 𝑢 merupakan variabel tak bebas. Pada persamaan (2.4) variabel 𝑥, 𝑦, dan 𝑧 merupakan variabel bebas dan variabel 𝑢 merupakan variabel tak bebas.


Definisi 2.6 Orde persamaan diferensial adalah tingkat tertinggi dari turunan yang muncul dalam persamaan diferensial tersebut. Contoh 2.6 Persamaan (2.1) adalah persamaan diferensial biasa orde pertama karena tingkat tertinggi yang muncul adalah tingkat satu. Persamaan (2.2) adalah contoh persamaan diferensial biasa orde dua karena tingkat turunan yang muncul adalah tingkat dua. Persamaan (2.3) dan (2.4) adalah persamaan diferensial parsial orde dua karena tingkat tertinggi dari turunanparsial yang muncul adalah tingkat dua. C.

Nilai Eigen dan Vektor Eigen Berikut akan dibahas mengenai nilai eigen dan vektor eigen beserta

dengan contohnya. Definisi 2.8 Misalkan 𝑨 adalah suatu matriks 𝑛 × 𝑛. Skalar 𝜆 disebut sebagai suatu nilai eigen atau nilai karakteristik dari 𝑨 jika terdapat suatu vektor taknol 𝑥̅ , sehingga 𝑨𝑥̅ = 𝜆𝑥̅ . Vektor 𝑥̅ disebut vektor eigen atau vektor karakteristik yang bersesuaian dengan nilai eigen λ dari 𝑨. Contoh 2.8 4 −2 2 ) dan 𝑥̅ = ( ). Karena 1 1 1

Misalkan 𝑨 = (


4 𝑨𝑥̅ = ( 1

6 −2 2 2 ) ( ) = ( ) = 3 ( ) = 3𝑥̅ 3 1 1 1

maka dari persamaan ini dapat dilihat bahwa 𝜆 = 3 adalah nilai eigen dari 𝐴 dan 2 𝑥̅ = ( ) merupakan vektor eigen yang bersesuaian dengan𝜆 = 3 tersebut, seperti 1 yang dijelaskan oleh Leon (2001).

𝑦

3 3𝑥̅

1 𝑥̅

2

6

𝑥

Secara geometris, perkalian matriks 𝑨 dengan vektor 𝑥̅ memiliki kelipatan 3 terhadap vektor 𝑥̅ . Ilustrasi secara geometris ditunjukkan dalam Gambar 2.2. Gambar 2.2. Ilustrasi geometri vektor eigen. D.

Klasifikasi Persamaan Diferensial Parsial Orde Dua Pada bagian ini akan dibahas tentang menentukan jenis suatu persamaan

diferensial parsial orde dua.


Persamaan diferensial parsial orde dua, yang linear homogen, dan memiliki koefisien konstan berbentuk 𝑎𝑢𝑥𝑥 + 2𝑏𝑢𝑥𝑦 + 𝑐𝑢𝑦𝑦 + 𝑑𝑢𝑥 + 𝑒𝑢𝑦 + 𝑓𝑢 = 0 dengan𝑢 = 𝑢(𝑥, 𝑦) dan 𝑎, 𝑏, 𝑐, 𝑑, 𝑒, 𝑓 adalah konstanta. Tiga suku pertama bentuk persamaan diferensial parsial linear homogen orde dua di atas disebut bagian utama persamaan diferensial parsial dan digunakan untuk menentukan jenis persamaan diferensial parsial. Dipandang bagian utama persamaan diferensial parsial:

2

𝑎𝑢𝑥𝑥 + 2𝑏𝑢𝑥𝑦 + 𝑐𝑢𝑦𝑦 = 𝑎

𝑎 𝑏

Matriks koefisien (

2

2

𝜕 𝑢 𝜕 𝑢 𝜕 𝑢 𝜕 + 2𝑏 +𝑐 2 =( 2 𝜕𝑥 𝜕𝑥𝜕𝑦 𝜕𝑦 𝜕𝑥

𝜕 𝑎 )( 𝜕𝑦 𝑏

𝜕 𝑏 𝜕𝑥 ) 𝜕 . 𝑐 (𝜕𝑦)

𝑏 ) merupakan matriks simetri yang mempunyai nilai eigen 𝑐

berupa bilangan real, 𝑎 𝑏

det ((

1 0 𝑏 𝑎−𝜆 )−𝜆( )) = det ( 0 1 𝑐 𝑏

𝑏 ) = (𝑎 − 𝜆)(𝑐 − 𝜆)−𝑏 2 𝑐−𝜆 (2.9)

= 𝜆2 − (𝑎 + 𝑐)𝜆 + 𝑎𝑐 − 𝑏 2

Jika 𝜆1 dan 𝜆2 adalah nilai eigen dari matriks 𝑨 = (

𝑎 𝑏

𝑏 ) maka 𝑐

persamaan karakteristiknya adalah (𝜆 − 𝜆1 )(𝜆 − 𝜆2 ) = 0

(2.10)


↔ 𝜆2 − (𝜆1 + 𝜆2 )𝜆 + 𝜆1 𝜆2 = 0 dari (2.9) dan (2.10) didapat: a. 𝜆1 + 𝜆2 = 𝑎 + 𝑐 = trace(𝐴) b. 𝜆1 𝜆2 = 𝑎𝑐 − 𝑏 2 = det(𝐴) Persamaan diferensial parsial disebut parabolik jika 𝑎𝑐 − 𝑏 2 = 0,yang artinya𝜆1 𝜆2 = 0; dengan kata lain, salah satu nilai eigennya bernilai 0. Persamaan diferensial parsial disebut eliptik apabila 𝑎𝑐 − 𝑏 2 > 0,yang artinya𝜆1 𝜆2 > 0; dengan kata lain, kedua nilai eigennya positif atau kedua nilai eigennya negatif. Persamaan

diferensial

parsial

disebut

hiperbolik

jika 𝑎𝑐 − 𝑏 2 < 0,yang

artinya𝜆1 𝜆2 < 0; dengan kata lain, salah satu nilai eigennya positif dan salah satu nilai eigennya negatif. E.

Penurunan Numeris Pada bagian ini akan dibahas mengenai penurunan numeris dan contohnya,

serta penjelasan tentang tiga pendekatan dalam menghitung turunan numeris yaitu pendekatan beda maju, beda pusat dan beda mundur. Definisi 2.9 Turunan suatu fungsi didefinisikan sebagai 𝑓 ′ (𝑥) = lim

ℎ→0

𝑓(𝑥+ℎ)−𝑓(𝑥) ℎ

.


Bila fungsi 𝑓(𝑥) diberikan secara eksplisit, maka kita dapat menentukan fungsi turunannya, 𝑓′(𝑥), 𝑓′′(𝑥), …, 𝑓 (𝑛+1) (𝑥), lalu menggunakannya untuk menghitung nilai turunan fungsi di 𝑥 = 𝑡. Namun demikian, seringkali fungsi 𝑓(𝑥) tidak diketahui secara eksplisit, tetapi hanya memiliki beberapa titik data saja. Pada kasus seperti ini, nilai turunan fungsi secara analitis susah untuk dicari. Seringkali 𝑓(𝑥) diketahui secara eksplisit, namun karena bentuk yang sulit maka untuk mencari hasil turunan fungsinya juga sulit, misalnya pada fungsi-fungsi berikut: √cos(2𝑥 2 )+𝑥tan(3𝑥)

(a). 𝑓(𝑥) = sin(𝑥)+𝑒 𝑥 −2𝑥/ cos(𝑥) (b). 𝑓(𝑥) = 𝑥𝑒 2𝑥+2 ln(4𝑥 2 ) Perhitungan nilai turunan pada fungsi (a) dan (b) dapat dikerjakan secara numeris. Nilai turunan yang diperoleh merupakan nilai hampiran dan diharapkan menghasilkan nilai galat yang kecil. 1. Tiga Pendekatan dalam Menghitung Turunan Numeris Turunan adalah limit dari hasil bagi selisih: yaitu pengurangan dua buah nilai yang besar (𝑓(𝑥 + ℎ) − 𝑓(𝑥)) dan membaginya dengan bilangan yang kecil (ℎ). Misal diberikan nilai-nilai 𝑥 di 𝑥0 − ℎ, 𝑥0 , dan 𝑥0 + ℎ, serta nilai fungsi untuk nilai-nilai 𝑥 tersebut. Titik-titik yang diperoleh adalah (𝑥−1 , 𝑓−1 ), (𝑥0 , 𝑓0 ), dan (𝑥1 , 𝑓1 ), yang dalam hal ini 𝑥−1 = 𝑥0 − ℎ dan 𝑥1 = 𝑥0 + ℎ. Terdapat tiga pendekatan dalam menghitung 𝑓′(𝑥0 ):


a. Hampiran Beda Maju Diketahui fungsi 𝑦 = 𝑓(𝑥0 ). Akan ditunjukkan 𝑓′(𝑥0 ) dengan hampiran beda maju 𝑓(𝑥0 + ℎ) − 𝑓(𝑥0 ) ℎ→0 ℎ

𝑓 ′ (𝑥0 ) = lim

≈

𝑓(𝑥0 + ℎ) − 𝑓(𝑥0 ) ℎ

=

𝑓1 − 𝑓0 ℎ

b. Hampiran Beda Mundur Diketahui fungsi 𝑦 = 𝑓(𝑥0 ). Akan ditunjukkan 𝑓′(𝑥0 ) dengan hampiran beda mundur 𝑓(𝑥0 ) − 𝑓(𝑥0 − ℎ) ℎ→0 ℎ

𝑓 ′ (𝑥0 ) = lim

≈

𝑓(𝑥0 ) − 𝑓(𝑥0 − ℎ) ℎ

=

𝑓0 − 𝑓−1 ℎ

c. Hampiran Beda Pusat Diketahui fungsi 𝑦 = 𝑓(𝑥0 ). Akan ditunjukkan 𝑓′(𝑥0 ) dengan hampiran beda pusat


𝑓(𝑥0 + ℎ) − 𝑓(𝑥0 − ℎ) ℎ→0 2ℎ

𝑓 ′ (𝑥0 ) = lim

≈

𝑓(𝑥0 + ℎ) − 𝑓(𝑥0 − ℎ) 2ℎ

=

𝑓1 − 𝑓−1 2ℎ

Tafsiran geometris dari ketiga pendekatan di atas diperlihatkan pada Gambar 2.3. 𝑦

𝑦

𝑦1

𝑦0

𝑦−1

ℎ

𝑥−1

𝑥0

𝑦 = 𝑓(𝑥)

𝑦0

𝑦 = 𝑓(𝑥)

𝑥1

ℎ

𝑥

𝑥−1

(a)

𝑥0

𝑥1

𝑥

(b) 𝑦 𝑦0

𝑦 = 𝑓(𝑥)

𝑦−1

2ℎ 𝑥−1

𝑥0

𝑥1

𝑥

(c) Gambar 2.3. (a) Hampiran beda maju. (b) Hampiran beda mundur. (c) Hampiran beda pusat.


2. Penurunan Rumus Turunan dengan Deret Taylor Misal diberikan titik-titik (𝑥𝑖 , 𝑓𝑖 ) dengan 𝑖 = 0,1,2, … , 𝑛, yang dalam hal ini 𝑥𝑖 = 𝑥0 + 𝑖ℎ dan 𝑓𝑖 = 𝑓(𝑥𝑖 ). Selanjutnya akan dihitung 𝑓′(𝑥), yang dalam hal ini 𝑥 = 𝑥0 + 𝑠ℎ, 𝑠 ∈ ℝ dengan ketiga pendekatan sebelumnya (beda maju, beda mundur, beda pusat). a. Hampiran Beda Maju Uraikan 𝑓(𝑥𝑖+1 ) di sekitar 𝑥𝑖 :

𝑓(𝑥𝑖+1 ) = 𝑓(𝑥𝑖 ) +

(𝑥𝑖+1 − 𝑥𝑖 )2 ′′ (𝑥𝑖+1 − 𝑥𝑖 ) ′ 𝑓 (𝑥𝑖 ) + 𝑓 (𝑥𝑖 ) + ⋯ 1! 2!

𝑓𝑖+1

= 𝑓𝑖 + ℎ𝑓𝑖′ + ℎ2 ⁄2𝑓′′(𝑥𝑖 ) + ⋯

ℎ𝑓𝑖 ′

= 𝑓𝑖+1 − 𝑓𝑖 − ℎ2 ⁄2𝑓 ′′ (𝑥𝑖 ) − ⋯

𝑓𝑖 ′

=

𝑓𝑖 ′

𝑓𝑖+1 − 𝑓𝑖 − ℎ⁄2𝑓 ′′ (𝑥𝑖 ) − ⋯ ℎ

=

𝑓𝑖+1 − 𝑓𝑖 − ℎ⁄2𝑓′′(𝑡) ℎ

(2.11)


𝑓𝑖 ′

=

𝑓𝑖+1 − 𝑓𝑖 + 𝑂(ℎ) ℎ

yang dalam hal ini, 𝑂(ℎ) = − ℎ⁄2𝑓′′(𝑡), untuk suatu 𝑡 dengan𝑥𝑖 < 𝑡 < 𝑥𝑖+1 . Untuk nilai-nilai 𝑓 di 𝑥0 dan 𝑥1 rumusnya menjadi:

𝑓0′ =

𝑓1 − 𝑓0 + 𝑂(ℎ) ℎ

(2.12)

yang dalam hal ini 𝑂(ℎ) = − ℎ⁄2𝑓′′(𝑡), untuk suatu 𝑡 dengan𝑥0 < 𝑡 < 𝑥1. b. Hampiran Beda Mundur Uraikan 𝑓(𝑥𝑖−1 ) di sekitar 𝑥𝑖 :

𝑓(𝑥𝑖−1 )

= 𝑓(𝑥𝑖 ) +

(𝑥𝑖−1 − 𝑥𝑖 ) ′ (𝑥𝑖−1 − 𝑥𝑖 )2 ′′ 𝑓 (𝑥𝑖 ) + 𝑓 (𝑥𝑖 ) + ⋯ 1! 2!

𝑓𝑖−1

= 𝑓𝑖 − ℎ𝑓𝑖′ + ℎ2 ⁄2𝑓′′(𝑥𝑖 ) + ⋯

ℎ𝑓𝑖 ′

= 𝑓𝑖 − 𝑓𝑖−1 + ℎ2 ⁄2𝑓′′(𝑥𝑖 ) + ⋯

𝑓𝑖

𝑓𝑖 ′

𝑓𝑖 ′

=

(2.13)

𝑓𝑖 − 𝑓𝑖−1 + ℎ⁄2𝑓 ′′ (𝑥𝑖 ) + ⋯ ℎ

=

𝑓𝑖 − 𝑓𝑖−1 + ℎ⁄2𝑓′′(𝑡) ℎ

=

𝑓𝑖 − 𝑓𝑖−1 + 𝑂(ℎ) ℎ

yang dalam hal ini, 𝑂(ℎ) = ℎ⁄2𝑓′′(𝑡), untuk suatu 𝑡 dengan𝑥𝑖−1 < 𝑡 < 𝑥𝑖 .


Untuk nilai-nilai 𝑓 di 𝑥0 dan 𝑥1 persamaan rumusnya menjadi: 𝑓1 − 𝑓0 + 𝑂(ℎ) ℎ

𝑓1′ =

(2.14)

yang dalam hal ini 𝑂(ℎ) = ℎ⁄2𝑓′′(𝑡), untuk suatu 𝑡 dengan 𝑥0 < 𝑡 < 𝑥1 . c. Hampiran Beda Pusat Kurangkan persamaan (2.13) dari persamaan (2.11): 𝑓𝑖+1 − 𝑓𝑖−1

= 2ℎ𝑓𝑖′ + ℎ3 ⁄3𝑓𝑖 ′′′ + ⋯

2ℎ𝑓𝑖′

= 𝑓𝑖+1 − 𝑓𝑖−1 − ℎ3 ⁄3𝑓𝑖′′′ − ⋯

𝑓𝑖 ′

=

𝑓𝑖 ′

𝑓𝑖+1 − 𝑓𝑖−1 − ℎ2 ⁄6𝑓𝑖′′′ − ⋯ 2ℎ

=

𝑓𝑖+1 − 𝑓𝑖−1 + 𝑂(ℎ2 ) 2ℎ

yang dalam hal ini, 𝑂(ℎ2 ) = −ℎ2 ⁄6𝑓′′′(𝑡), untuk suatu 𝑡 dengan 𝑥𝑖−1 < 𝑡 < 𝑥𝑖+1. Untuk nilai-nilai 𝑓 di 𝑥−1 dan 𝑥1 persamaan rumusnya menjadi:

𝑓1′ =

𝑓2 − 𝑓0 + 𝑂(ℎ2 ) 2ℎ

(2.12)

yang dalam hal ini 𝑂(ℎ2 ) = −ℎ2 ⁄6𝑓′′′(𝑡), untuk suatu 𝑡 dengan 𝑥0 < 𝑡 < 𝑥2 . Hampiran beda pusat ini memiliki hasil yang lebih baik dari dua hampiran lainnya karena orde galatnya adalah 𝑂(ℎ2 ).


F. Matriks Tridiagonal Pada bagian ini akan dibahas mengenai definisi matriks tridiagonal dan contohnya. Definisi 2.10 Misalkan 𝑛 ≥ 3. Matriks 𝑻 = (𝑡𝑖𝑗 ) ∈ ℝ𝑛×𝑛 disebut matriks tridiagonal jika elemen-elemen yang berada pada selain diagonal utama dan dua diagonal berdekatan bernilai nol, yaitu 𝑡𝑖𝑗 = 0 jika |𝑖 − 𝑗| > 1, 𝑖, 𝑗 ∈ {1,2, … , 𝑛} matriks tersebut juga sering disebut tiga diagonal. Untuk penjelasan lebih jelasnya dapat dilihat pada buku karangan Süli dan Mayers (2007). Contoh 2.9 Berikut ditunjukan beberapa matriks 1 𝑨 = (7 7

2 3 1 2 5 3), 𝑩 = (2 3 8 9 0 1

0 1 0 2), 𝑪 = (0 1 2 0 0

0 0). 1

Dari ketiga matriks di atas, matriks tridiagonal ditunjukan oleh matriks 𝑩 dan 𝑪. Matriks 𝑨 bukan matriks tridiagonal karena𝑎13 ≠ 0 dan 𝑎31 ≠ 0. Matriks 𝑪adalah matriks identitas. Matriks 𝑪memenuhi definisi matriks tridiagonal karena 𝑎13 = 𝑎31 = 0.


BAB III MODEL PERGERAKAN LAPISAN FLUIDA Pada bab ini akan dibahas tentang pemodelan pergerakan lapisan fluida, penurunan gerak fluida satu dimensi, serta metode volume hingga dan metode beda hingga untuk model pergerakan lapisan fluida. A. Penurunan Persamaan Gerak Lapisan Fluida Persamaan gerak fluida pada kasus ini dideskripsikan dengan persamaan panas seperti yang dijelaskan oleh Caldwell dan Ng Douglas (2004). Hal ini dikarenakan gerakan fluida seperti menjalar dari sumber gerakan. Plat atas yang ditarik secara konstan adalah sumber gerakan awal, kecepatan fluida yang bersentuhan langsung dengan plat sama dengan kecepatan plat yang ditarik secara konstan tersebut, sedangkan kecepatan fluida yang berada jauh dari plat atas tersebut memiliki kecepatan yang lebih kecil dari pada kecepatan fluida yang bersentuhan langsung dengan plat yang ditarik. Persamaan panas dapat juga disebut sebagai persamaan difusi. Persamaan panas dapat diformulasikan dengan merumuskan persamaan aliran panas (Haberman, 2004). Misalkan kawat penampang 𝐴 berorientasi terhadap arah 𝑥 seperti yang diilustrasikan pada Gambar 3.1. Jumlah energi panas per satuan volume sebagai variabel yang tidak diketahui disebut kepadatan energi panas.

26


𝐴

∅(𝑥 + ∆𝑥, 𝑡)

∅(𝑥, 𝑡)

𝑥=0

𝑥

𝑥 + ∆𝑥

𝑥=𝐿

Gambar 3.1 Kawat satu dimensi dengan energi panas yang mengalir masuk dan keluar. Di sini,𝐴 adalah luas penampang kawat, dan ∅(𝑥, 𝑡) adalah energi panas yang lewat di penampang kawat pada posisi 𝑥 dan waktu 𝑡. Asumsikan pada setiap waktu 𝑡, suhu di dalam kawat pada posisi 𝑥seragam yaitu 𝑢(𝑥, 𝑡), tetapi berbeda bila dibandingkan suhu penampang kawat posisi yang lain. Akan dicari distribusi suhu penampang kawat pada setiap posisi 𝑥 dan pada setiap waktu 𝑡, yaitu 𝑢(𝑥, 𝑡), ∀(𝑥, 𝑡). Misalkan 𝑐 konstanta yang menyatakan berapa banyak energi yang dibutuhkan oleh satu unit massa suatu benda untuk menaikkan suhu sebesar 1 derajat. Segmen kawat dari 𝑥 ke 𝑥 + ∆𝑥 mempunyai massa: 𝜌=

𝑚 , 𝑉

𝑚 = 𝜌 ∗ 𝑉 = 𝜌 ∗ 𝐴 ∗ ∆𝑥, dengan 𝜌 adalah kepadatan kawat, 𝑚 adalah massa, dan 𝑉 merupakan volume kawat. Sehingga untuk menaikan suhu segmen kawat sebesar 1 derajat dibutuhkan energi sebanyak 𝜌 ∗ 𝑐 ∗ 𝐴 ∗ ∆𝑥. Apabila suhunya naik dari 0 ke 𝑢(𝑥, 𝑡) maka energi yang dibutuhkan sebesar 𝜌 ∗ 𝑐 ∗ 𝐴 ∗ ∆𝑥 ∗ 𝑢(𝑥, 𝑡). Jadi, total energi panas pada segmen tersebut untuk 𝑡 > 0 adalah 𝐸 = lim ∑ 𝜌 ∗ 𝑐 ∗ 𝐴 ∗ ∆𝑥 ∗ 𝑢(𝑥, 𝑡), ∆𝑥→∞

atau

𝑥


𝑥+∆𝑥

𝐸 = ∫ 𝜌 ∗ 𝑐 ∗ 𝐴 ∗ 𝑢(𝑠, 𝑡) 𝑑𝑠. 𝑥

Fluks Panas Fluks panas adalah laju perubahan energi panas yang melewati suatu penampang. Fluks dapat dihitung dengan cara: 𝑥+∆𝑥

𝜕𝐸 𝜕 Fluks = = ∫ 𝜌 ∗ 𝑐 ∗ 𝐴 ∗ 𝑢(𝑠, 𝑡) 𝑑𝑠 𝜕𝑡 𝜕𝑡 𝑥

𝑥+∆𝑥

𝜕𝑢(𝑠, 𝑡) = ∫ 𝜌∗𝑐∗𝐴∗ 𝑑𝑠 𝜕𝑡

(3.1)

𝑥

atau dengan cara: Fluks = 𝐴∅(𝑥, 𝑡) − 𝐴∅(𝑥 + ∆𝑥, 𝑡) = −𝐴[∅(𝑥 + ∆𝑥, 𝑡) − ∅(𝑥, 𝑡)].

(3.2)

Karena panas menjalar dari benda bersuhu tinggi ke rendah dan banyak energi berbanding dengan perbedaan suhu di antara 2 titik (Hukum Newton Pendinginan) maka: ∅(𝑥, 𝑡) = −𝑘

𝜕𝑢(𝑥, 𝑡) . 𝜕𝑥

(3.3)

Substitusi persamaan (3.3) ke persamaan (3.2) didapat Fluks

= −𝐴 (−𝑘

𝜕𝑢(𝑥 + ∆𝑥, 𝑡) 𝜕𝑢(𝑥, 𝑡) +𝑘 ) 𝜕𝑥 𝜕𝑥

𝜕𝑢(𝑥 + ∆𝑥, 𝑡) 𝜕𝑢(𝑥, 𝑡) = 𝐴𝑘 ( − ) 𝜕𝑥 𝜕𝑥 𝑥+∆𝑥

= 𝐴𝑘 ∫ 𝑥

𝜕 𝜕𝑢(𝑠, 𝑡) ( ) 𝑑𝑠 𝜕𝑥 𝜕𝑥

(3.4)


Dari persamaan (3.1) dan (3.4) didapat 𝑥+∆𝑥

∫ 𝜌∗𝑐∗𝐴∗ 𝑥

𝑥+∆𝑥 𝜕𝑢(𝑠, 𝑡) 𝜕 𝜕𝑢(𝑠, 𝑡) 𝑑𝑠 = 𝐴𝑘 ∫ ( ) 𝑑𝑠 𝜕𝑡 𝜕𝑥 𝜕𝑥 𝑥

𝑥+∆𝑥

𝐴∫

(𝜌 ∗ 𝑐 ∗

𝑥

(𝜌 ∗ 𝑐 ∗

𝜕𝑢 𝜕 2𝑢 − 𝑘 ∗ 2 ) 𝑑𝑠 = 0 𝜕𝑡 𝜕𝑥

𝜕𝑢 𝜕 2𝑢 − 𝑘 ∗ 2) = 0 𝜕𝑡 𝜕𝑥

atau 𝜕𝑢 𝑘 𝜕 2𝑢 = 𝜕𝑡 𝜌𝑐 𝜕𝑥 2 𝜕𝑢 𝜕 2𝑢 =𝐾 2 𝜕𝑡 𝜕𝑥

(3.5)

dengan 𝐾 = 𝑘⁄𝜌𝑐 adalah koefisien difusi. Pada persamaan gerak fluida, koefisien difusi 𝐾 diganti dengan 𝜇 yang berarti kekentalan fluida. Pada kasus pergerakan lapisan fluida, persamaan (3.5) diberikan subskrip minyak dan air guna membedakan antara persamaan gerak untuk minyak dan persamaan gerak untuk air, seperti pada persamaan (1) dan persamaan (2). Persamaan (3.5) merupakan persamaan diferensial parsial parabolik. Hal ini dikarenakan bagian utama persamaan diferensialnya berbentuk: 𝐾

𝜕 2𝑢 = 0, 𝜕𝑥 2

sehingga det(𝑨) = 0 atau dengan kata lain salah satu nilai eigen dari persamaan tersebut bernilai 0. Di sini,


𝑨=(

𝑎 𝑏

𝐾 𝑏 )=( 0 𝑐

0 ) 0

seperti yang dijelaskan pada subbab “D. Klasifikasi Persamaan Diferensial Parsial Orde Dua”. B. Masalah Pergerakan Fluida Diketahui persamaan gerak dari lapisan fluida untuk air, yaitu: 𝜕𝑣air 𝜕 2 𝑣air = 𝜇air 𝜕𝑡 𝜕𝑥 2

(3.6)

𝜕𝑣minyak 𝜕 2 𝑣minyak = 𝜇minyak 𝜕𝑡 𝜕𝑥 2

(3.7)

dan untuk minyak yaitu:

dengan x adalah variabel ruang, t adalah variabel waktu, 𝑣air adalah kecepatan air, 𝑣minyak adalah kecepatan minyak, 𝜇air adalah kekentalan air dan 𝜇minyak adalah kekentalan minyak. Hubungan di titik perbatasan antara minyak dan air ditunjukkan dengan: 𝑣minyak = 𝑣air

(3.8)

𝜕𝑣minyak 𝜕𝑣air = 𝜇air . 𝜕𝑥 𝜕𝑥

(3.9)

dan 𝜇minyak


Persamaan (3.9) merupakan definisi dari tegangan gesek yaitu:

𝜏=𝜇

dengan

𝜇

merupakan

viskositas,

𝑑𝑉 𝑑𝑦

𝑉

merupakan

kecepatan

fluida,

merepresentasikan jarak dua plat rata yang disusun secara horisontal, dan

𝑦 𝑑𝑉 𝑑𝑦

merepresentasikan gradien dari kecepatan fluida. Pada kasus dalam skripsi ini, jarak dua plat rata horisontal direpresentasikan dengan 𝑥, dan viskositas diberikan subscrip minyak dan air sebagai pembeda koefisien viskositas untuk minyak dan air. Dalam skripsi ini tidak akan dibahas lebih lanjut tentang bagaimana mendapatkan definisi tegangan gesek. Materi tentang vikositas dan tegangan gesek dapat dilihat pada buku-buku atau jurnal tentang mekanika fluida seperti yang ditulis oleh Crowe. C. T., Elger D. F., Williams B. C., dan Roberson. J. A. pada buku berjudul Engineering Fluid Mechanics (2010). Akan disimulasikan pergerakan lapisan fluida dalam kasus ini adalah antara minyak dan air yang berada diantara dua plat rata dengan jarak 10 cm menggunakan metode volume hingga Lax-Friedrichs, dan metode beda hingga dengan menggunakan MATLAB. Pada kasus ini plat atas ditarik sehingga mempunyai kecepatan konstan sebesar 7cm/s seperti pada Gambar 1.1. Pada kasus ini, terdapat dua masalah nilai awal yang melibatkan kecepatan dua fluida, sebagai berikut:


𝜕𝑣air 𝜕 2 𝑣air = 𝜇air ,0 < 𝑥 < 6 𝜕𝑡 𝜕𝑥 2 𝑣air (𝑥, 0) = 0 𝑣air (0, 𝑡) = 0 𝑣air (6, 𝑡) = 𝑣minyak (6, 𝑡) 𝜕𝑣minyak 𝜕𝑣air 𝜇minyak | = 𝜇air | 𝜕𝑥 𝜕𝑥 𝑥=6 { 𝑥=6

(3.10)

𝜕𝑣minyak 𝜕 2 𝑣minyak = 𝜇minyak , 6 < 𝑥 < 10 𝜕𝑡 𝜕𝑥 2 𝑣minyak (𝑥, 0) = 0 𝑣minyak (10, 𝑡) = 7 𝑣minyak (6, 𝑡) = 𝑣air (6, 𝑡) 𝜕𝑣minyak 𝜕𝑣air 𝜇minyak | = 𝜇air | 𝜕𝑥 𝜕𝑥 𝑥=6 { 𝑥=6

(3.11)

dengan 𝑣air (𝑥, 𝑡)adalah kecepatan air, 𝑣minyak (𝑥, 𝑡) adalah kecepatan minyak, 𝜇air sebagai kekentalan air, and 𝜇minyak kekentalan minyak. Untuk model pergerakan lapisan fluida, dibuat beberapa asumsi sebagai berikut: 1. Pada plat atas tidak terdapat kekentalan. 2. Plat atas bergerak secara konstan yaitu 7 cm/s. 3. Kekentalan fluida minyak dan air diberikan oleh Caldwel dan Ng Douglas, K. S. (2004). 4. Kekentalan fluida minyak dan air diasumsikan tetap, tidak berubah terhadap suhu. 5. Aliran fluida hanya dalamsatu arah, yaitu arah yang tegak lurus sumbu x. 6. Aliran fluidanya bersifat laminer.


C. Solusi Analitis Masalah Pergerakan Lapisan Fluida Masalah pergerakan lapisan fluida sangat sulit diselesaikan secara analitis untuk kasus aliran tak tunak. Akan tetapi, penyelesaian numeris dapat dibandingkan dengan solusi analitis untuk kasus aliran tunak. Aliran tak tunak adalah kondisi dimana komponen aliran berubah terhadap waktu, dan aliran tunak adalah kondisi dimana komponen aliran tidak berubah terhadap waktu. Untuk kasus aliran tunak, solusi analitis tidak bergantung terhadap waktu. Dengan demikian, untuk kasus aliran tunak, solusi analitis 𝑣air (𝑥, 𝑡) = 𝑣air (𝑥) dan 𝑣minyak (𝑥, 𝑡) = 𝑣minyak (𝑥). Dalam kasus aliran tunak persamaan (3.10) menjadi 𝜕 2 𝑣air = 0, 𝜕𝑥 2

0 < 𝑥 < 6,

𝑣air (0) = 0.

(3.12) (3.13)

Persamaan (3.12) memiliki penyelesaian𝑣air (𝑥) = 𝐴𝑥 + 𝐵.Karena 𝑣air (0) = 0 maka 𝐵 = 0, sehingga penyelesaian untuk persamaan (3.12) adalah 𝑣air (𝑥) = 𝐴𝑥,

0 ≤ 𝑥 ≤ 6.

(3.14)

𝜕𝑣air =𝐴 𝜕𝑥 Persamaan (3.11) untuk kasus aliran tunak dapat ditulis menjadi 𝜕 2 𝑣minyak = 0, 𝜕𝑥 2

6 ≤ 𝑥 ≤ 10,

(3.15)

𝑣minyak (10) = 7,

(3.16)

𝑣minyak (6) = 𝑣air (6),

(3.17)


𝜇minyak

𝜕𝑣minyak 𝜕𝑣air | = 𝜇air | . 𝜕𝑥 𝜕𝑥 𝑥=6 𝑥=6

(3.18)

Persamaan (3.15) menghasilkanpenyelesaian: 𝑣minyak (𝑥) = 𝐶𝑥 + 𝐷,

6 ≤ 𝑥 ≤ 10,

(3.19)

𝜕𝑣minyak =𝐶 𝜕𝑥 dan pada titik batas yakni persamaan (3.16) dan (3.17) ditulis menjadi 10𝐶 + 𝐷 = 7

(3.20)

6𝐶 + 𝐷 = 6𝐴.

(3.21)

dan

Selanjutnya karena 𝑣air𝑥 = 𝐴 dan 𝑣minyak𝑥 = 𝐶 maka persamaan (3.18) berlaku 𝜇minyak (𝐶) = 𝜇air (𝐴).

(3.22)

Eliminasi persamaan (3.20) dan (3.21) sehingga mendapat 𝐴=

7 − 4𝐶 . 6

(3.23)

Substitusi persamaan (3.23) ke persamaan (3.22) akan menghasilkan 𝐶=

7𝜇air , 6𝜇minyak + 4𝜇air

(3.24)

substitusikan pula persamaan (3.24) ke persamaan (3.22) sehingga didapat 𝐴=

7𝜇minyak , 6𝜇minyak + 4𝜇air

(3.25)

substitusikan persamaan (3.24) dan (3.25) ke persamaan (3.21) didapat 𝐷=

21(𝜇minyak − 𝜇air ) . 6𝜇minyak + 4𝜇air

(3.26)


Berikut adalah solusi aliran tunak yang dihasilkan dengan mensubstitusikan persamaan (3.24), (3.25), dan (3.26) ke persamaan (3.14) dan (3.19):

𝑣air (𝑥) =

𝑣minyak (𝑥) =

7𝜇minyak 𝑥, 6𝜇minyak + 4𝜇air

0 ≤ 𝑥 ≤ 6,

21(𝜇minyak − 𝜇air ) 7𝜇air 𝑥+ , 6𝜇minyak + 4𝜇air 3𝜇minyak + 2𝜇air

(3.27)

(3.28)

6 ≤ 𝑥 ≤ 10,

dengan 𝑣air (𝑥) adalah kecepatan air, 𝑣minyak (𝑥) adalah kecepatan minyak, 𝜇air menyatakan kekentalan air dan 𝜇minyak menyatakan kekentalan minyak. Solusi di atas akan digunakan dalam perhitungan simulasi numeris dengan MATLAB.

D. Metode Volume Hingga Lax-Friedrichs Pada bagian ini dibahas mengenai skema metode volume hingga, perhitungan flux secara numeris dalam metode volume hingga dan solusi numeris metode volume hingga Lax-Friedrichs. 1. Skema Metode Volume Hingga Persamaan diferensial parsial hukum kekekalan berbentuk 𝑞𝑡 + 𝑓(𝑞)𝑥 = 0 atau ditulis 𝜕 𝜕 𝑞(𝑥, 𝑡) + 𝑓(𝑞(𝑥, 𝑡)) = 0. 𝜕𝑡 𝜕𝑥


Skema metode volume hingga berdasar pada pendiskretan domain pada ruang ke dalam interval, seperti ditunjukkan dalam Gambar 3.2.

𝑥𝑖+1

𝑥𝑖−1

𝑥𝑖−3

𝑥𝑖+3

2

2

2

𝑥𝑖−1

2

𝑥𝑖+1

𝑥𝑖

Gambar 3.2.Ilustrasi diskretisasi domain ruang.

Di sini ∆𝑥 = 𝑥𝑖 − 𝑥𝑖−1 atau ∆𝑥 = 𝑥𝑖+1 − 𝑥𝑖−1 . Domain waktu didiskretkan 2

2

menjadi 𝑡 𝑛 = 𝑛 ∙ ∆𝑡 dengan 𝑛 = 0,1,2,3, …. Misalkan 𝑄𝑖𝑛 adalah nilai pendekatan rata-rata volume kuantitas 𝑞(𝑥, 𝑡) dalam interval ke-i pada waktu 𝑡 𝑛 , yaitu: 𝑥 1

1

𝑄𝑖𝑛 ≈ ∆𝑥 ∫𝑥 𝑖+12 𝑞(𝑥, 𝑡 𝑛 )𝑑𝑥. 𝑖−2

𝑛 Misalkan pula 𝐹𝑖− 1 adalah pendekatan dari rata-rata fluks(debit material) 2

𝑓(𝑞(𝑥, 𝑡)) di titik 𝑥𝑖−1 , yaitu 2

𝑛+1

𝑛 𝐹𝑖− 1 2

1 𝑡 ≈ ∫ ∆𝑡 𝑡 𝑛

𝑓 (𝑞 (𝑥𝑖−1 , 𝑡)) 𝑑𝑡. 2

Bentuk integral dari hukum kekekalan diberikan oleh: 𝜕

𝑥 1

∫ 𝑖+2 𝑞(𝑥, 𝑡 𝑛 )𝑑𝑥 = − [𝑓 (𝑞 (𝑥𝑖+1 , 𝑡)) − 𝑓 (𝑞 (𝑥𝑖−1 , 𝑡))], 𝜕𝑡 𝑥 1 𝑖− 2

2

2

dengan nilai-nilai pendekatan diperoleh untuk 𝑞𝑡 + 𝑓(𝑞)𝑥 = 0, yaitu


𝑄𝑖𝑛+1 − 𝑄𝑖𝑛 =− ∆𝑡

𝑛 𝑛 𝐹𝑖+ 1 −𝐹 1 𝑖− 2

2

∆𝑥

atau dapat ditulis menjadi ∆𝑡

𝑛 𝑛 𝑄𝑖𝑛+1 = 𝑄𝑖𝑛 − ∆𝑥 (𝐹𝑖+ 1 − 𝐹 1 ). 𝑖− 2

2

Persamaan di atas merupakan skema volume hingga bagi𝑞𝑡 + 𝑓(𝑞)𝑥 = 0. Skema metode volume hingga tersebut konsisten dengan skema metode beda hingga karena 𝑄𝑖𝑛+1 − 𝑄𝑖𝑛 =− ∆𝑡

𝑛 𝑛 𝐹𝑖+ 1 −𝐹 1 𝑖− 2

2

∆𝑥

dapat ditulis menjadi 𝑄𝑖𝑛+1 − 𝑄𝑖𝑛 + ∆𝑡

𝑛 𝑛 𝐹𝑖+ 1 −𝐹 1 𝑖− 2

2

∆𝑥

=0

yang merupakan suatu bentuk diskret dari 𝑞𝑡 + 𝑓(𝑞)𝑥 = 0. 2. Perhitungan FluksSecara Numeris dalam Metode Volume Hingga Diberikan persamaan diferensial parsial dengan bentuk hukum kekekalan 𝑞𝑡 + 𝑓(𝑞)𝑥 = 0. 𝑛 𝑛 Misal 𝑄𝑖𝑛 ≈ 𝑞(𝑥𝑖 , 𝑡 𝑛 ) dan 𝐹𝑖− 1 ≈ 𝑓 (𝑞 (𝑥 1 , 𝑡 )), seperti telah dijelaskan pada 𝑖− 2

2

bagian Skema Metode Volume Hingga di muka.Skema metode volume hingga untuk persamaan di atas adalah 𝑄𝑖𝑛+1 = 𝑄𝑖𝑛 −

∆𝑡 𝑛 𝑛 (𝐹 1 − 𝐹𝑖− 1 ). ∆𝑥 𝑖+2 2

Diketahui 𝑄𝑖𝑛 merupakan nilai kuantitas numeris di titik 𝑥𝑖 dan pada waktu 𝑡 𝑛 .Oleh karena itu,fluks di titik 𝑥𝑖 pada waktu 𝑡 𝑛 diketahui, yaitu


𝐹𝑖𝑛 ≈ 𝑓(𝑞(𝑥𝑖 , 𝑡 𝑛 )) ≈ 𝑓(𝑄𝑖𝑛 ).

Metode Stabil dan Tidak Stabil Metode numeris dikatakan stabil apabila galat atau error yang muncul disetiap iterasi tidak membesar terlalu cepat pada iterasi-iterasi berikutnya. Jika galat yang muncul pada suatu iterasi membesar menuju tak hingga maka metode tersebut dikatakan tidak stabil. Teori tentang kestabilan tidak akan dibahas pada skripsi ini. Teori kestabilan dapat dilihat dalam buku-buku referensi misalnya LeVeque (1992,2002).

1. Flukstak stabil 𝑛 Akan didefinisikan rata-rata fluks pada titik 𝑥𝑖−1 berdasarkan pada 𝑄𝑖−1 2

dan 𝑄𝑖𝑛 , sebagai berikut: 1

𝑛 𝑛 𝑛 𝑛 𝑛 𝐹𝑖− 1 = 𝐹(𝑄𝑖−1 , 𝑄𝑖 ) = [𝑓(𝑄𝑖−1 ) + 𝑓(𝑄𝑖 )]. 2 2

Dengan demikian, skema metode volume hingga menjadi 𝑄𝑖𝑛+1 = 𝑄𝑖𝑛 −

∆𝑡 𝑛 𝑛 (𝐹 1 − 𝐹𝑖− 1) ∆𝑥 𝑖+2 2

menjadi 𝑄𝑖𝑛+1 = 𝑄𝑖𝑛 −

∆𝑡 1 1 𝑛 𝑛 ) + 𝑓(𝑄𝑖𝑛 )] − [𝑓(𝑄𝑖−1 ) + 𝑓(𝑄𝑖𝑛 )]) ( [𝑓(𝑄𝑖+1 ∆𝑥 2 2 ∆𝑡

𝑛 𝑛 𝑄𝑖𝑛+1 = 𝑄𝑖𝑛 − 2∆𝑥 [𝑓(𝑄𝑖+1 ) − 𝑓(𝑄𝑖−1 )] .


Akan tetapi, skema metode volume hingga ini tidak stabil. 2. Fluks Lax-Friedrichs Skema Lax-Friedrichs adalah skema yang memodifikasi skema metode volume hingga di atas, dengan 1 𝑛 𝑛 ) 𝑄𝑖𝑛 = (𝑄𝑖+1 + 𝑄𝑖−1 2 sehingga skema Lax-Friedrichsmenjadi 1

∆𝑡

𝑛 𝑛 𝑛 𝑛 )− 𝑄𝑖𝑛+1 = 2 (𝑄𝑖+1 + 𝑄𝑖−1 [𝑓(𝑄𝑖+1 ) − 𝑓(𝑄𝑖−1 ]. 2∆𝑥

Skema Lax-Friedrichs ini stabil untuk ∆𝑡 yang cukup kecil. 3.

Solusi Numeris Metode Volume Hingga Lax-Friedrichs Masalah pergerakan lapisan fluida dapat diselesaikan dengan menggunakan

metode volume hingga Lax-Friedrichs. Diberikan persamaan lapisan fluida (3.10) dan (3.11) yaitu ∂𝑣air 𝜕 2 𝑣air = 𝜇air , 𝜕𝑡 𝜕𝑥 2

0≤𝑥≤6

atau dapat ditulis 𝑣air 𝑡 + (−𝜇air 𝑣air 𝑥 )𝑥 = 0

(3.29)

dan 𝜕𝑣minyak 𝜕 2 𝑣minyak = 𝜇minyak , 𝜕𝑡 𝜕𝑥 2

6 ≤ 𝑥 ≤ 10

atau 𝑣minyak + (−𝜇minyak 𝑣minyak )𝑥 = 0. 𝑡

𝑥

Persamaan (3.29) mempunyai skema metode volume hingga

(3.30)


∆𝑡

𝑛 𝑛 𝑄𝑖𝑛+1 = 𝑄𝑖𝑛 − ∆𝑥 (𝐹𝑖+ 1 − 𝐹 1 ). 𝑖− 2

2

Jadi, jika diketahui persamaan (3.29) maka didapat 𝑞 = 𝑣air dan 𝑓(𝑞) = 𝑛 𝑛 −𝜇air 𝑣air𝑥 .Sekarang akan dicari fluks 𝐹𝑖+ 1 dan 𝐹 1 dari persamaan (3.29), yaitu: 𝑖− 2

𝑛 𝐹𝑖+ 1 2

𝑛 𝐹𝑖− 1 2

2

=

1 ∆𝑥 𝑛 [𝑓(𝑄𝑖+1 ) + 𝑓(𝑄𝑖𝑛 )] − (𝑄 𝑛 − 𝑄𝑖𝑛 ) 2 2∆𝑡 𝑖+1

=

1 ∆𝑥 ((𝑣air )𝑛𝑖+1 − (𝑣air )𝑛𝑖 ) [(−𝜇air 𝑣air𝑥 )𝑛𝑖+1 + (−𝜇air 𝑣air𝑥 )𝑛𝑖 ] − 2 2∆𝑡

=

1 ∆𝑥 𝑛 𝑛 (−𝜇air ) [(𝑣air𝑥 ) + (𝑣air𝑥 ) ] − ((𝑣air )𝑛𝑖+1 − (𝑣air )𝑛𝑖 ), 𝑖+1 𝑖 2 2∆𝑡

=

1 ∆𝑥 𝑛 [𝑓(𝑄𝑖𝑛 ) + 𝑓(𝑄𝑖−1 )] − (𝑄 𝑛 − 𝑄𝑖𝑛 ) 2 2∆𝑡 𝑖+1

=

1 ∆𝑥 ((𝑣air )𝑛𝑖 − (𝑣air )𝑛𝑖−1 ) [(−𝜇air 𝑣air𝑥 )𝑛𝑖 + (−𝜇air 𝑣air𝑥 )𝑛𝑖−1 ] − 2 2∆𝑡

=

1 ∆𝑥 𝑛 𝑛 (−𝜇air ) [(𝑣air𝑥 ) + (𝑣air𝑥 ) ] − ((𝑣air )𝑛𝑖 − (𝑣air )𝑛𝑖−1 ). 𝑖 𝑖−1 2 2∆𝑡

Persamaan (3.30) juga mempunyai skema metode volume hingga ∆𝑡

𝑛 𝑛 𝑄𝑖𝑛+1 = 𝑄𝑖𝑛 − ∆𝑥 (𝐹𝑖+ 1 − 𝐹 1 ). 𝑖− 2

2

Jadi, jika diketahui persamaan (3.30) maka didapat 𝑞 = 𝑣minyak dan 𝑓(𝑞) = 𝑛 𝑛 −𝜇minyak 𝑣minyak𝑥 .Sekarang akan dicari fluks𝐹𝑖+ 1 dan 𝐹 1 dari persamaan (3.30), 𝑖− 2

yaitu: 𝑛 𝐹𝑖+ 1 2

=

1 ∆𝑥 𝑛 𝑛 [𝑓(𝑄𝑖+1 ) + 𝑓(𝑄𝑖𝑛 )] − (𝑄𝑖+1 − 𝑄𝑖𝑛 ) 2 2∆𝑡

2


=

1 [(−𝜇minyak 𝑣minyak𝑥 )𝑛𝑖+1 + (−𝜇minyak 𝑣minyak𝑥 )𝑛𝑖 ] 2 −

=

1 𝑛 𝑛 (−𝜇minyak ) [(𝑣minyak𝑥 ) + (𝑣minyak𝑥 ) ] 𝑖+1 𝑖 2 −

𝑛 𝐹𝑖− 1 2

∆𝑥 𝑛 𝑛 ((𝑣minyak ) − (𝑣minyak ) ) 𝑖+1 𝑖 2∆𝑡

∆𝑥 𝑛 𝑛 ((𝑣minyak ) − (𝑣minyak ) ), 𝑖+1 𝑖 2∆𝑡

=

1 ∆𝑥 𝑛 𝑛 [𝑓(𝑄𝑖𝑛 ) + 𝑓(𝑄𝑖−1 )] − (𝑄𝑖+1 − 𝑄𝑖𝑛 ) 2 2∆𝑡

=

1 [(−𝜇minyak 𝑣minyak𝑥 )𝑛𝑖 + (−𝜇minyak 𝑣minyak𝑥 )𝑛𝑖−1 ] 2 −

=

∆𝑥 𝑛 𝑛 ((𝑣minyak ) − (𝑣minyak ) ) 𝑖 𝑖−1 2∆𝑡

1 𝑛 𝑛 (−𝜇minyak ) [(𝑣minyak𝑥 ) + (𝑣minyak𝑥 ) ] 𝑖 𝑖−1 2 −

∆𝑥 𝑛 𝑛 ((𝑣minyak ) − (𝑣minyak ) ). 𝑖 𝑖−1 2∆𝑡

Hasil simulasi penyelesaian masalah pergerakan lapisan fluida dengan metode volume hingga Lax-Friedrichs dengan menggunakan perangkat lunak MATLAB ditunjukkan oleh Gambar 3.3. Pada hasil simulasi pergerakan lapisan fluida diberikan nilai 𝜇air = 1 dan 𝜇minyak = 3, program dijalankan dengan ∆𝑥 = 0.1, ∆𝑡 = 0.001 ∗ ∆𝑥, dan 𝑡𝐹𝑖𝑛𝑎𝑙 = 50.


Gambar 3.3. Hasil simulasi penyelesaian masalah pergerakan lapisan fluida dengan metode volume hingga saat 𝑡 = 50. Terlihat pada gambar bahwa terjadi patahan pada saat 𝑥 = 6. Hal ini terjadi karena diketahui hubungan di titik perbatasan antara minyak dan air yang sudah dijelaskan sebelumnya.

E. Metode Beda Hingga untuk Model Pergerakan Lapisan Fluida Pada bagian ini dibahas mengenai skema metode beda hingga untuk model pergerakan lapisan fluida, dan solusi numeris metode beda hingga untuk model pergerakan lapisan fluida.


1. Skema Metode Beda Hingga untuk Model Pergerakan Lapisan Fluida Persamaan (3.10) dan (3.11) tidak dapat diselesaikan secara terpisah, karena terdapat beberapa kondisi yang saling berkaitan satu dengan yang lainnya. Dengan menggunakan skema implisit, persamaan gerak fluida untuk air pada persamaan (3.10) dan persamaan gerak fluida untuk minyak pada persamaan (3.11) dapat ditulis menjadi: 𝑣air 𝑙+1 − 𝑣air 𝑙𝑖 𝑣air 𝑙+1 − 2𝑣air 𝑙+1 + 𝑣air 𝑙+1 𝑖 𝑖+1 𝑖 𝑖−1 = 𝜇air , (∆𝑥)2 ∆𝑡

(3.31)

untuk 𝑖 = 1,2,3, … , 𝑚, 𝑣minyak 𝑙+1 − 𝑣minyak 𝑙 𝑗

𝑗

∆𝑡 = 𝜇minyak

𝑣minyak 𝑙+1 − 2𝑣minyak 𝑙+1 + 𝑣minyak 𝑙+1 𝑗+1

𝑗

(∆𝑥)2

𝑗−1

,

(3.32)

untuk 𝑗 = 𝑚, 𝑚 + 1, 𝑚 + 2, … , 𝑛 − 1, Persamaan (3.31) dan (3.32) dapat ditulis ulang menjadi: 𝜇air 2𝜇air 1 𝜇air 1 ( ) 𝑣air 𝑙+1 −( + ) 𝑣air 𝑙+1 +( ) 𝑣air 𝑙+1 = − 𝑣air 𝑙𝑖 , 2 𝑖−1 2 𝑖 2 𝑖+1 (∆𝑥) (∆𝑥) (∆𝑥) ∆𝑡 ∆𝑡

(3.33)

untuk𝑖 = 1, 2, 3, … , 𝑚, (

𝜇minyak 2𝜇minyak 1 𝑙+1 ) 𝑣 − ( + ) 𝑣minyak 𝑙+1 minyak 𝑗−1 𝑗 (∆𝑥)2 (∆𝑥)2 ∆𝑡 +(

𝜇minyak 1 ) 𝑣minyak 𝑙+1 = − 𝑣minyak 𝑙 , 2 𝑗+1 𝑗 (∆𝑥) ∆𝑡

untuk𝑗 = 𝑚, 𝑚 + 1, 𝑚 + 2, … , 𝑛 − 1.

(3.34)


Untuk posisi di𝑖 = 1 dan syarat awal 𝑣air (0, 𝑡) = 𝑣air 𝑙+1 = 0 , persamaan 𝑖 (3.33) dapat ditulis menjadi: 2𝜇air 1 𝜇air 1 −( + ) 𝑣air1𝑙+1 + ( ) 𝑣air 𝑙+1 = − 𝑣air1𝑙 2 2 2 (∆𝑥) (∆𝑥) ∆𝑡 ∆𝑡

(3.35)

Pada posisi batas antara minyak dan air (𝑖 = 𝑚dan𝑗 = 𝑚), persamaan (3.33) dan (3.34) menjadi (

𝜇air 2𝜇air 1 𝜇air ) 𝑣air 𝑙+1 −( + ) 𝑣air 𝑙+1 +( ) 𝑣 𝑙+1 2 𝑚−1 2 𝑚 (∆𝑥) (∆𝑥) (∆𝑥)2 air 𝑚+1 ∆𝑡 1 = − 𝑣air 𝑙𝑚 , ∆𝑡

(3.36)

𝜇minyak 2𝜇minyak 1 ( ) 𝑣minyak 𝑙+1 − ( + ) 𝑣minyak 𝑙+1 2 𝑚−1 𝑚 (∆𝑥) (∆𝑥)2 ∆𝑡 𝜇minyak 1 𝑙+1 𝑙 +( ) 𝑣 = − 𝑣 . minyak 𝑚+1 (∆𝑥)2 ∆𝑡 minyak 𝑚

(3.37)

Kondisi pada posisi batas dapat dijabarkan menjadi: 𝜇minyak

𝜕𝑣minyak 𝜕𝑣air | = 𝜇air | , 𝜕𝑥 𝜕𝑥 𝑥=6 𝑥=6

𝑣minyak 𝑙+1 − 𝑣minyak 𝑙+1 𝑣air 𝑙+1 − 𝑣air 𝑙+1 𝑚+1 𝑚−1 𝑚+1 𝑚−1 𝜇air = μminyak . 2∆𝑥 2∆𝑥

(3.38)

(3.39)

Persamaan (3.38) dapat ditulis sebagai: 𝑣air 𝑙+1 = 𝑚+1

𝜇minyak (𝑣minyak 𝑙+1 − 𝑣minyak 𝑙+1 ) + 𝑣air 𝑙+1 . 𝑚−1 𝑚+1 𝑚−1 𝜇air

(3.40)

Substitusi persamaan (3.40) ke persamaan (3.36), didapat: 2𝜇air 2𝜇air 1 ( ) 𝑣air 𝑙+1 −( + ) 𝑣air 𝑙+1 2 𝑚−1 2 𝑚 (∆𝑥) (∆𝑥) ∆𝑡 𝜇minyak 1 +( ) (𝑣minyak 𝑙+1 − 𝑣minyak 𝑙+1 ) = − 𝑣air 𝑙𝑚 . 2 𝑚+1 𝑚−1 (∆𝑥) ∆𝑡

(3.41)


Jumlahkan persamaan (3.37) dengan persamaan (3.41) sehingga didapat: (

2𝜇minyak 1 2𝜇air 2𝜇air 1 𝑙+1 𝑙+1 ) v − ( + ) 𝑣 − ( + ) 𝑣minyak 𝑙+1 air air 𝑚−1 𝑚 𝑚 (∆𝑥)2 (∆𝑥)2 ∆𝑡 (∆𝑥)2 ∆𝑡 (3.42)

2𝜇minyak 1 1 𝑙+1 𝑙 +( ) 𝑣 = − 𝑣 − 𝑣minyak 𝑙 . minyak air 𝑚 𝑚+1 𝑚 (∆𝑥)2 ∆𝑡 ∆𝑡 Karena 𝑣air 𝑙+1 = 𝑣minyak 𝑙+1 , maka persamaan (3.42) dapat ditulis menjadi: 𝑚 𝑚

2𝜇air + 2𝜇minyak 2 2𝜇air 𝑙+1 ( ) 𝑣 − ( + ) 𝑣minyak 𝑙+1 air 𝑚−1 𝑚 (∆𝑥)2 (∆𝑥)2 ∆𝑡 2𝜇minyak 1 1 𝑙+1 𝑙 𝑙 +( ) 𝑣 = − 𝑣 − 𝑣 . minyak air 𝑚 𝑚+1 (∆𝑥)2 ∆𝑡 ∆𝑡 minyak 𝑚

(3.43)

Saat 𝑗 = 𝑚 + 1 persamaan (3.34) menjadi: (

𝜇minyak 2𝜇minyak 1 ) 𝑣minyak 𝑙+1 − ( + ) 𝑣minyak 𝑙+1 2 𝑚 𝑚+1 (∆𝑥) (∆𝑥)2 ∆𝑡 𝜇minyak 1 +( ) 𝑣minyak 𝑙+1 = − 𝑣minyak 𝑙 . 2 𝑚+2 𝑚+1 (∆𝑥) ∆𝑡

(3.44)

Pada plat atas (𝑗 = 𝑛 − 1), persamaan (3.34) dapat ditulis menjadi: 𝜇minyak 2𝜇minyak 1 ( ) 𝑣minyak 𝑙+1 − ( + ) 𝑣minyak 𝑙+1 2 𝑛−2 𝑛−1 (∆𝑥) (∆𝑥)2 ∆𝑡 𝜇minyak 1 = − 𝑣minyak 𝑙 −( ) 𝑣minyak 𝑙+1 . 𝑚+1 𝑛 (∆𝑥)2 ∆𝑡

(3.45)

Karena 𝑣minyak (10, 𝑡) = 𝑣minyak 𝑙+1 = 7, maka: 𝑛

𝜇minyak 2𝜇minyak 1 𝑙+1 ( ) 𝑣 − ( + ) 𝑣minyak 𝑙+1 minyak 𝑛−2 𝑛−1 (∆𝑥)2 (∆𝑥)2 ∆𝑡 7𝜇minyak 1 = − 𝑣minyak 𝑙 −( ). 𝑚+1 (∆𝑥)2 ∆𝑡

(3.46)


Persamaan (3.33), (3.34), (3.35), (3.43), (3.44) dan (3.46) adalah persamaan yang mewakili semua titik diantara 0 sampai 10. Keenam persamaan merupakan sistem tridiagonal yang dapat diselesaikan dengan menggunakan perintah 𝐴\𝑥 pada MATLAB. Misalkan: 𝜇air b=( ), (∆𝑥)2

c=(

2𝜇air 1 + ), 2 (∆𝑥) ∆𝑡

𝜇minyak d=( ), (∆𝑥)2

2𝜇minyak 1 2𝜇air + 2𝜇minyak 2 e=( + ) , dan f = ( + ) 2 (∆𝑥) (∆𝑥)2 ∆𝑡 ∆𝑡 contoh membentuk sistem tridiagonal dengan ∆𝑥 = 2 adalah sebagai berikut: −𝑐 𝑏 𝐴= 0 0 0 [ 0

𝑏 0 −𝑐 𝑏 𝑏 −𝑐 0 2𝑏 0 0 0 0

0 0 0 0 0 0 𝑏 0 0 , −𝑓 2𝑑 0 𝑑 −𝑒 𝑑 0 𝑑 −𝑒]

𝑣air1𝑙

𝑣air 𝑙2 𝑣air 𝑙3

𝑥̅ = 𝑣 𝑙 , minyak 4

𝑣minyak 𝑙

5

𝑙 [𝑣minyak ] 6

1 𝑣 𝑙 ∆𝑡 air1 1 − 𝑣air 𝑙2 ∆𝑡 1 − 𝑣air 𝑙3 ∆𝑡 𝑏̅ = . 1 1 𝑙 𝑙 − 𝑣air 4 − 𝑣minyak 4 ∆𝑡 ∆𝑡 1 − 𝑣minyak 𝑙 5 ∆𝑡 1 − 𝑣minyak 𝑙 [ ] 6 ∆𝑡 −

Sistem tridiagonal di atas merupakan penyelesaian pada metode beda hingga. Dengan variasi ∆𝑥 yang berbeda akan terbentuk sistem tridiagonal yang besarnya berbeda pula. Semakin kecil ∆𝑥 maka semakin besar sistem tridiagonalnya.


Sistem tridiagonal akan diselesaikan dengan menggunakan perintah 𝐴\𝑥 pada MATLAB.

2. Solusi Numeris Metode Beda Hingga untuk Masalah Pergerakan Lapisan Fluida. Hasil simulasi pergerakan lapisan fluida dengan metode beda hingga dengan menggunakan program MATLAB ditunjukkan dalam Gambar 3.4 sampai dengan Gambar

3.9.

Simulasi

ini

dilakukan

untuk

beberapa

nilai∆𝑥 =

2, 1, 0.5, 0.25, 0.125, 0.0625 dan ∆𝑡 = 0.5 ∗ ∆𝑥. Untuk jarak kedua plat adalah 10 cm dan waktu 50 detik.

Gambar 3.4. Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 2 dengan metode beda hingga.


Gambar 3.5. Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 1 dengan metode beda hingga.

Gambar 3.6. Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 0.5 dengan metode beda hingga.


Gambar 3.7. Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 0.25 dengan metode beda hingga.

Gambar 3.8. Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 0.125 dengan motode beda hingga.


Gambar 3.9. Hasil simulasi masalah pergerakan lapisan fluida untuk ∆𝑥 = 0.0625 dengan metode beda hingga. Terlihat pada gambar-gambar hasil simulasi untuk metode beda hingga bahwa terjadi patahan pada saat 𝑥 = 6. Hal ini terjadi karena diketahui hubungan di titik perbatasan antara minyak dan air yang sudah dijelaskan sebelumnya.


BAB IV ANALISIS HASIL SIMULASI Pada bab ini akan dibahas mengenai hasil simulasi numeris untuk metode beda hingga dan metode volume hingga Lax-Friedrichs. Simulasi numeris dilakukan dengan menggunakan MATLAB dengan jarak antara plat bawah dan plat atas adalah 10 cm, dan plat atas ditarik dengan kecepatan konstan 7 cm/s. Galat atau error dihitung dengan menggunakan rumus 𝑛

1 Galat = ∑|𝑣eksak (𝑖) − 𝑣numeris (𝑖)| 𝑛 𝑖=1

dengan 𝑣eksak adalah nilai eksak di titik 𝑥𝑖 , 𝑣numeris adalah nilai numeris di titik 𝑥𝑖 , dan 𝑛 adalah banyaknya data yang ada di domain ruang. Menghitung galat saja masih belum cukup, seberapa cepat suatu metode konvergen juga harus diperhatikan. Untuk mengetahui seberapa cepat konvergen dari simulasi ini, dihitung dengan menggunakan rumus:

Perbandingan Galat =

𝑒𝑟𝑟𝑜𝑟𝑖+1 𝑒𝑟𝑟𝑜𝑟𝑖

Dengan 𝑒𝑟𝑟𝑜𝑟𝑖+1 merupakan galat pada titik 𝑥𝑖+1 dan merupakan galat pada titik 𝑥𝑖+1.

51


A.

Metode Volume Hingga Lax-Friedrichs Pada Bab sebelumnya telah dibahas tentang solusi numeris untuk masalah

pergerakan lapisan fluida dengan metode volume hingga Lax-Friedrichs. Solusi yang didapat adalah hasil siulasi dengan program MATLAB. Simulasi ini menggunakan ∆𝑥 = 0.1, ∆𝑡 = 0.01 ∗ ∆𝑥dan 𝑡𝐹𝑖𝑛𝑎𝑙 = 50 dengan nilai 𝜇air = 1dan 𝜇minyak = 3. Grafik simulasi ditunjukkan pada Gambar 3.3. Terlihat pada gambar bahwa solusi numeris mendekati solusi eksak. Namun masih terdapat ruang antara solusi numeris dan solusi eksak yang didapat, dengan kata lain masih terdapat galat. Berikut merupakan galat solusi numeris untuk masalah pergerakan lapisan fluida menggunakan metode volume hingga LaxFriedrichs, seperti dirangkum pada Tabel 4.1, untuk beberapa variansi nilai ∆𝑥. Tabel 4.1. Galat hasil simulasi metode volume hingga Lax-Friedrichs ∆𝒙

Galat atau error

1

0.7463

0.5

0.7197

0.96435

0.25

0.6604

0.91760

0.125

0.3752

0.56814

Perbandingan Galat


Dari Tabel 4.1 terlihat bahwa semakin kecil ∆𝑥semakin kecil juga galat yang dihasilkan. Ketika diambil∆𝑥yang sangat kecil maka galat yang dihasilkan akan semakin kecil. Hal ini dikarenakan banyaknya langkah pada ruang di sumbu 𝑥. Simulasi pada metode volume hingga Lax-Friedrichs dalam skripsi ini berhenti pada saat ∆𝑥 = 0.125 karena perhitungan yang sangat lama. Walaupun berhenti pada saat ∆𝑥 = 0.125, galat pada metode ini sudah mendekati 0 dan perbandingan galatnya adalah 0.56814 yang berarti kecepatan konvergensi metode ini adalah 0.5 dan sudah cukup baik. Ilustrasi galat secara grafik ditunjukkan pada Gambar 4.1.

Gambar 4.1. Ilustrasi geometris galat metode volume hingga Lax-Friedrichs.


Terlihat pada Gambar 4.1 bahwa semakin besar ∆𝑥 maka semakin besar juga galat yang akan muncul. Kelebihan dari metode ini adalah komputasi yang cukup mudah. Kekurangan dari metode ini adalah untuk ∆𝑥 yang semakin kecil,diperlukan ∆𝑡 yang juga semakin kecil agar perhitungan stabil (metodenya stabil), sehingga waktu perhitungan menjadi sangat lama. B.

Metode Beda Hingga Pada Bab III telah dibahas tentang solusi numeris untuk masalah

pergerakan lapisan fluida dengan metode beda hingga. Solusi yang diperoleh merupakan hasil simulasi dengan menggunakan MATLAB. Grafik simulasi ditunjukkan pada Gambar 3.4 sampai dengan Gambar 3.9.Semakin kecil ∆𝑥gambar yang ditunjukkan semakin mulus dan mendekati solusi eksak. Simulasi ini menggunakan nilai 𝜇air = 1 dan 𝜇minyak = 3, program dijalankan menggunakan∆𝑥 bervariasi dengan∆𝑡 = 0.001 ∗ ∆𝑥, dan 𝑡𝐹𝑖𝑛𝑎𝑙 = 50. Terlihat pada gambar, semakin kecil∆𝑥dan semakin besar𝑡, selisih antara solusi numeris dan solusi eksaknya semakin kecil. Dengan kata lain, solusi numerisnya mendekati solusi eksak. Berikut merupakan hasil galat simulasi numeris untuk masalah pergerakan lapisan fluida menggunakan metode beda hingga.


Tabel 4.2. Galat hasil simulasi metode beda hingga ∆𝒙

Galat atau error

Perbandingan Galat

2

0.0898

1

0.0964

1.07350

0.5

0.07888

0.81743

0.25

0.0547

0.69416

0.125

0.0334

0.61060

0.0625

0.0187

0.55988

Tabel 4.2 menunjukkan semakin kecil ∆𝑥 semakin kecil juga galat yang dihasilkan.Hal ini dikarenakan banyaknya langkah pada ruang di sumbu 𝑥. Simulasi pada metode beda hingga dalam skripsi ini berhenti pada saat ∆𝑥 = 0.0625 karena galat yang didapat dari metode ini sudah sangat kecil dan perbandingan galatnya adalah 0.55988 yang berarti kecepatan konvergensi metode ini adalah 0.5 dan sudah cukup baik. Ilustrasi galat secara grafik ditunjukkan pada Gambar 4.2.


Gambar 4.2. Ilustrasi geometris galat metode beda hingga. Gambar 4.2 menunjukkan bahwa semakin besar ∆𝑥 maka semakin besar juga galat yang akan muncul.Komputasi yang mudah dan tidak memerlukan waktu lama merupakan kelebihan metode ini.Kekurangan metode ini adalah semakin kecil∆𝑥mengakibatkan semakin besar matriks tridiagonal yang dibutuhkan, sehingga memerlukan memori komputer yang juga semakin besar.


BAB V KESIMPULAN DAN SARAN

A. Kesimpulan Persamaan diferensial mampu memodelkan pergerakan lapisan fluida yang melibatkan minyak dan air. Dalam kasus ini, model pergerakan lapisan fluida antara minyak dan air diselesaikan dengan menggunakan metode volume hingga Lax-Friedrichs dan metode beda hingga. Hasil numeris yang diperoleh menunjukkan kecepatan dua fluida yang berada dalam dua plat rata horisontal dengan plat atas ditarik secara konstan 7 cm/s mendekati hasil analitisnya. Lebih lanjut lagi, dari hasil simulasi dengan menggunakan program MATLAB, solusi metode beda hingga lebih baik daripadametode volume hingga Lax-Friedrichs, hal ini dilihat dari hasil galat yang diperoleh. Metode beda hingga memiliki galat yang mendekati 0 lebih cepat daripada metode volume hingga LaxFriedrichs seiring dengan ∆𝑥 yang semakin kecil, dan pehitungan menggunakan metode beda hingga lebih cepat bila dibandingkan dengan metode volume hingga Lax-Friedrichs. B. Saran Penulis sadar bahwa dalam penulisan skripsi ini masih banyak kekurangan. Oleh sebab itu, penulis mengharapkan kelak ada yang melanjutkan penelitian ini. Tulisan ini hanya membahas penyelesaian model pergerakan

57


lapisan fluida dengan metode volume hingga Lax-Friedrichs dan metode beda hingga. Penulis berharap di waktu yang akan datang, ada yang melanjutkan penulisan ini dengan metode yang lain dan lebih baik.


DAFTAR PUSTAKA Blanchard. P., Devaney, R.L. & Hall, G.R. (2012). Differential Equations,Fourth edition. Boston: Brooks/Cole. Bleecker, D. & Csordas, G. (1992). Basic Partial Differential Equations. New York: Van Nostrand Reinhold. Caldwell, J. & Ng Douglas, K. S. (2004). Mathematical Modelling. New York: Kluwer Academic Publisher. Crowe, C. T., Elger, D. F., Williams, B. C. & Roberson, J. A. (2010). Engineering Fluid Mechanics, Ninth Edition. New Jersey: John Wiley & Sons (Asia) Pte Ltd. Haberman, R.(2004). Applied Partial Differential Equations, Fourth edition. New Jersey: Pearson Prentice Hall. LeVeque R.J. (2002). Finite Volume Methods for Hyperbolic Problems. Cambridge: Cambridge University Press. LeVeque R J 1992 Numerical Methods for Conservation Laws (Basel: Springer)

Leon, S. J. (2001). Aljabar Linear dan Aplikasinya. Jakarta: Erlangga. Munir, R. (2008). Metode Numerik. Bandung: Informatika Bandung. Purcell, E J. (1981). Kalkulus dan Geometri Analitis edisi ketiga. Jakarta: Erlangga. Süli, E. & Mayers, D F. (2007). An Introduction to Numerical Analysis. Cambridge: Cambridge University Press. Thomas, G. B. (2010). Thomas’ Calculus Early Transcendentals. Boston: Pearson Education.

59


LAMPIRAN Berikut ini merupakan code pada program MATLAB untuk solusi analitik, solusi numeris metode volume hingga Lax-Friedrichs dan metode beda hingga beserta code untuk menggambar galatnya. A. Metode Volume Hngga Lax-Friedrichs tic clc clear close all dx=0.1; %delta x dt=0.001*dx; %delta t x=0:dx:10; %jarak x k = find(x==6); tFinal=50; nx=length(x); mum=3; mua=1; va=zeros(1,nx); %matriks untuk v air vm=zeros(1,nx); %matriks untuk v minyak %kondisi yang diberikan va(1:2)=0; vm(nx-1:nx)=7; vaNew=va; vmNew=vm; v=zeros(nx,1); for j=1:nx if j<=k v(j)=(7*mum*x(j))/(6*mum+4*mua); else v(j)=((7*mua*x(j))/(6*mum+4*mua))+((21*mum21*mua)/(3*mum+2*m ua)); end end for n=3:ceil(tFinal/dt) for i=3:nx-2 if i==k-1 %Untuk air Fim1a=-mua*(va(i)-va(i-2))/(2*dx); Fia=-mua*(va(i+1)-va(i-1))/(2*dx); Fip1a=(-mum*vm(i+2)+mua*va(i))/(2*dx); %Untuk minyak Fim1m=-mua*(va(i)-va(i-2))/(2*dx); Fim=-mua*(va(i+1)-va(i-1))/(2*dx);

60


Fip1m=(-mum*vm(i+2)+mua*va(i))/(2*dx); %F untuk mencari Va FUa=0.5*(Fip1a+Fia) - dx/(2*dt)*(va(i+1)-va(i)); FLa=0.5*(Fia+Fim1a) - dx/(2*dt)*(va(i)-va(i-1)); %F untuk mencari Vm FUm=0.5*(Fip1m+Fim) - dx/(2*dt)*(va(i+1)-va(i)); FLm=0.5*(Fim+Fim1m) - dx/(2*dt)*(va(i)-va(i-1)); %Untuk mencari nilai Va dan Vm vaNew(i)=va(i) - dt/dx*(FUa-FLa); vmNew(i)=vm(i) - dt/dx*(FUm-FLm); elseif i==k %Untuk air Fim1a=-mua*(va(i) - va(i-2))/(2*dx); Fia=(-mum*vm(i+1) + mua*va(i-1))/(2*dx); Fip1a=(-mum*vm(i+2) + mua*va(i))/(2*dx); %Untuk minyak Fim1m=-mua*(va(i) - va(i-2))/(2*dx); Fim=(-mum*vm(i+1) + mua*va(i-1))/(2*dx); Fip1m=(-mum*vm(i+2) + mua*va(i))/(2*dx); %F untuk mencari Va FUa=0.5*(Fip1a+Fia) - dx/(2*dt)*(vm(i+1)-va(i)); FLa=0.5*(Fia+Fim1a) - dx/(2*dt)*(va(i)-va(i-1)); %F untuk mencari Vm FUm=0.5*(Fip1m+Fim) - dx/(2*dt)*(vm(i+1)-va(i)); FLm=0.5*(Fim+Fim1m) - dx/(2*dt)*(va(i)-va(i-1)); %Untuk mencari nilai Va dan Vm vaNew(i)=va(i) - dt/dx*(FUa-FLa); vmNew(i)=vm(i) - dt/dx*(FUm-FLm); elseif i==k+1 %Untuk air Fim1a=(-mum*vm(i)+mua*va(i-2))/(2*dx); Fia=(-mum*vm(i+1)+mua*va(i-1))/(2*dx); Fip1a=(-mum*vm(i+2)+mum*vm(i))/(2*dx); %Untuk minyak Fim1m=(-mum*vm(i)+mua*va(i-2))/(2*dx); Fim=(-mum*vm(i+1)+mua*va(i-1))/(2*dx); Fip1m=-mum*(vm(i+2)-vm(i))/(2*dx); %F untuk mencari Va FUa=0.5*(Fip1a+Fia) - dx/(2*dt)*(vm(i+1)-vm(i)); FLa=0.5*(Fia+Fim1a) - dx/(2*dt)*(vm(i)-va(i-1)); %F untuk mencari Vm FUm=0.5*(Fip1m+Fim) - dx/(2*dt)*(vm(i+1)-vm(i)); FLm=0.5*(Fim+Fim1m) - dx/(2*dt)*(vm(i)-va(i-1)); %Untuk mencari nilai Va dan Vm vaNew(i)=va(i) - dt/dx*(FUa-FLa); vmNew(i)=vm(i) - dt/dx*(FUm-FLm); else %Untuk air Fim1a=-mua*(va(i) - va(i-2))/(2*dx); Fia=-mua*(va(i+1) - va(i-1))/(2*dx); Fip1a=-mua*(va(i+2) - va(i))/(2*dx); %Untuk minyak


Fim1m=-mum*(vm(i) - vm(i-2))/(2*dx); Fim=-mum*(vm(i+1) - vm(i-1))/(2*dx); Fip1m=-mum*(vm(i+2) - vm(i))/(2*dx); %F untuk mencari Va FUa=0.5*(Fip1a+Fia) - dx/(2*dt)*(va(i+1)-va(i)); FLa=0.5*(Fia+Fim1a) - dx/(2*dt)*(va(i)-va(i-1)); %F untuk mencari Vm FUm=0.5*(Fip1m+Fim) - dx/(2*dt)*(vm(i+1)-vm(i)); FLm=0.5*(Fim+Fim1m) - dx/(2*dt)*(vm(i)-vm(i-1)); %Untuk mencari nilai Va dan Vm vaNew(i)=va(i) - dt/dx*(FUa-FLa); vmNew(i)=vm(i) - dt/dx*(FUm-FLm); end end va = vaNew; vm = vmNew; %kondisi yang diberikan va(1:2)=0; vm(nx-1:nx)=7; V=[va(1:k) vm(k+1:nx)]; plot(V,x,v,x) xlabel('Velocity of fluid layer') ylabel('Position of fluid layer') legend('Numerical solution','Exact solution') pause(0.00000000001) end Error = sum(abs(V-v'))/nx toc

B. Gambar Galat Metode Volume Hingga Lax-Friedrichs clc clear all x=[0.125,0.25,0.5,1]; Galat=[0.3754,0.6604,0.7197,0.7463]; plot(x,Galat,'r*-') legend('Galat Kecepatan Fluida','Location','southeast') xlabel('Steps in space domain') ylabel('Error') title('Galat Perhitungan Kecepatan Lapisan Fluida')

C. Metode Beda Hingga ∆𝑥 = 2 close all clear clc mua=1;%input('Kekentalan air: '); mum=3;%input('Kekentalan minyak: '); dx=2;%input('delta x: ') dt=0.5*dx; %input('delta t: ')


k=6; t=0:dt:50; x=0:dx:10; xA=0:dx:6; nxA=length(xA); xM=6:dx:10; nxM=length(xM); n=length(x);%input('masukkan besar matriks: '); %perhitungan Analitik v=zeros(n,1); for j=1:n if x(j)<=k v(j)=(7*mum*x(j))/(6*mum+4*mua); else v(j)=((7*mua*x(j))/(6*mum+4*mua))+((21*mum21*mua)/(3*mum+2*mua)); end end %perhitungan numerik V=zeros(n,1); g=1/dt; for i=1:length(t) V0=V; %nilai awal V0(n,1)=7; b=mua/(dx^2); c=(2*mua/(dx^2))+g; d=mum/(dx^2); e=(2*mum/(dx^2))+g; f=((2*mua+2*mum)/(dx^2))+(2*g); dimin=ones(nxA-1,1); %diagonal utama u/ 0:6 diminb=ones(nxA-2,1); %diagonal atas u/ 0:6 di=ones(1,1); %untuk 6 diplus=ones(nxM-1,1); %diagonal utama dan atas u/7:10 diplusb=ones(nxM-2,1); %diagonal bwh u/7:10 for j=0:dx:6 q3=b*diminb; w3=-c*dimin; r3=b*dimin; k1=-g*V0(1:nxA-1,1); end for j=6 q2=2*b*di; w2=-f*di; r2=2*d*di; k2=-g*V0(nxA,1)-g*V0(nxA,1); end for j=6:dx:10 q1=d*diplus; w1=-e*diplus; r1=d*diplusb; k3=-g*V0(nxA+1:n-1,1);


k4=-g*V0(n,1)-7*d; end q=[q3;q2;q1]; w=[w3;w2;w1]; r=[r3;r2;r1]; k=[k1;k2;k3;k4]; A=diag(q,-1)+diag(w)+diag(r,1); V=A\k; V(1,1)=0; V(n,1)=7; V; if t(i)==1 V05=V; elseif t(i)==5 V5=V; elseif t(i)==50 V50=V; end end plot(V05,x,'r+-', V5,x,'g.-',V50,x,'b--', v,x,'k-') legend('Numerical velocity at t=1', 'Numerical velocity at t=5', 'Numerical velocity at t=50', 'Exact steady velocity', 'Location', 'southeast') title('Kecepatan Fluida saat delta x = 2') xlabel('v (velocity)') ylabel('x (position)') error=sum(abs(V-v))/n

D. Metode Beda Hingga ∆𝑥 = 1 close all clear clc mua=1;%input('Kekentalan air: '); mum=3;%input('Kekentalan minyak: '); dx=1;%input('delta x: ') dt=0.5*dx; %input('delta t: ') k=6; t=0:dt:50; x=0:dx:10; xA=0:dx:6; nxA=length(xA); xM=6:dx:10; nxM=length(xM); n=length(x);%input('masukkan besar matriks: '); %perhitungan Analitik v=zeros(n,1); for j=1:n if x(j)<=k v(j)=(7*mum*x(j))/(6*mum+4*mua); else


v(j)=((7*mua*x(j))/(6*mum+4*mua))+((21*mum21*mua)/(3*mum+2*mua)); end end %perhitungan numerik V=zeros(n,1); g=1/dt; for i=1:length(t) V0=V; %nilai awal V0(n,1)=7; b=mua/(dx^2); c=(2*mua/(dx^2))+g; d=mum/(dx^2); e=(2*mum/(dx^2))+g; f=((2*mua+2*mum)/(dx^2))+(2*g); dimin=ones(nxA-1,1); %diagonal utama u/ 0:6 diminb=ones(nxA-2,1); %diagonal atas u/ 0:6 di=ones(1,1); %untuk 6 diplus=ones(nxM-1,1); %diagonal utama dan atas u/7:10 diplusb=ones(nxM-2,1); %diagonal bwh u/7:10 for j=0:dx:6 q3=b*diminb; w3=-c*dimin; r3=b*dimin; k1=-g*V0(1:nxA-1,1); end for j=6 q2=2*b*di; w2=-f*di; r2=2*d*di; k2=-g*V0(nxA,1)-g*V0(nxA,1); end for j=6:dx:10 q1=d*diplus; w1=-e*diplus; r1=d*diplusb; k3=-g*V0(nxA+1:n-1,1); k4=-g*V0(n,1)-7*d; end q=[q3;q2;q1]; w=[w3;w2;w1]; r=[r3;r2;r1]; k=[k1;k2;k3;k4]; A=diag(q,-1)+diag(w)+diag(r,1); V=A\k; V(1,1)=0; V(n,1)=7; V if t(i)==1 V05=V; elseif t(i)==5 V5=V;


elseif t(i)==50 V50=V; end end plot(V05,x,'r+-', V5,x,'g.-',V50,x,'b--', v,x,'k-') legend('Numerical velocity at t=1', 'Numerical velocity at t=5', 'Numerical velocity at t=50', 'Exact steady velocity', 'Location', 'southeast') title('Kecepatan Fluida saat delta x = 1') xlabel('v (velocity)') ylabel('x (position)') error=sum(abs(V-v))/n

E. Metode Beda Hingga ∆𝑥 = 0.5 close all clear clc mua=1;%input('Kekentalan air: '); mum=3;%input('Kekentalan minyak: '); dx=0.5;%input('delta x: ') dt=0.5*dx; %input('delta t: ') k=6; t=0:dt:50; x=0:dx:10; xA=0:dx:6; nxA=length(xA); xM=6:dx:10; nxM=length(xM); n=length(x);%input('masukkan besar matriks: '); %perhitungan Analitik v=zeros(n,1); for j=1:n if x(j)<=k v(j)=(7*mum*x(j))/(6*mum+4*mua); else v(j)=((7*mua*x(j))/(6*mum+4*mua))+((21*mum21*mua)/(3*mum+2*mua)); end end %perhitungan numerik V=zeros(n,1); g=1/dt; for i=1:length(t) V0=V; %nilai awal V0(n,1)=7; b=mua/(dx^2); c=(2*mua/(dx^2))+g;


d=mum/(dx^2); e=(2*mum/(dx^2))+g; f=((2*mua+2*mum)/(dx^2))+(2*g); dimin=ones(nxA-1,1); %diagonal utama u/ 0:6 diminb=ones(nxA-2,1); %diagonal atas u/ 0:6 di=ones(1,1); %untuk 6 diplus=ones(nxM-1,1); %diagonal utama dan atas u/7:10 diplusb=ones(nxM-2,1); %diagonal bwh u/7:10 for j=0:dx:6 q3=b*diminb; w3=-c*dimin; r3=b*dimin; k1=-g*V0(1:nxA-1,1); end for j=6 q2=2*b*di; w2=-f*di; r2=2*d*di; k2=-g*V0(nxA,1)-g*V0(nxA,1); end for j=6:dx:10 q1=d*diplus; w1=-e*diplus; r1=d*diplusb; k3=-g*V0(nxA+1:n-1,1); k4=-g*V0(n,1)-7*d; end q=[q3;q2;q1]; w=[w3;w2;w1]; r=[r3;r2;r1]; k=[k1;k2;k3;k4]; A=diag(q,-1)+diag(w)+diag(r,1); V=A\k; V(1,1)=0; V(n,1)=7; V; if t(i)==1 V05=V; elseif t(i)==5 V5=V; elseif t(i)==50 V50=V; end end plot(V05,x,'r+-', V5,x,'g.-',V50,x,'b--', v,x,'k-') legend('Numerical velocity at t=1', 'Numerical velocity at t=5', 'Numerical velocity at t=50', 'Exact steady velocity', 'Location', 'southeast') title('Kecepatan Fluida saat delta x = 0.5') xlabel('v (velocity)') ylabel('x (position)') error=sum(abs(V-v))/n


F. Metode Beda Hingga ∆𝑥 = 0.25 clear clc close all mua=1;%input('Kekentalan mum=3;%input('Kekentalan dx=0.25;%input('delta x: dt=0.5*dx; %input('delta

air: '); minyak: '); ') t: ')

k=6; t=0:dt:50; x=0:dx:10; xA=0:dx:6; nxA=length(xA); xM=6:dx:10; nxM=length(xM); n=length(x);%input('masukkan besar matriks: '); %perhitungan Analitik v=zeros(n,1); for j=1:n if x(j)<=k v(j)=(7*mum*x(j))/(6*mum+4*mua); else v(j)=((7*mua*x(j))/(6*mum+4*mua))+((21*mum21*mua)/(3*mum+2*mua)); end end %perhitungan numerik V=zeros(n,1); g=1/dt; for i=1:length(t) V0=V; %nilai awal V0(n,1)=7; b=mua/(dx^2); c=(2*mua/(dx^2))+g; d=mum/(dx^2); e=(2*mum/(dx^2))+g; f=((2*mua+2*mum)/(dx^2))+(2*g); dimin=ones(nxA-1,1); %diagonal utama u/ 0:6 diminb=ones(nxA-2,1); %diagonal atas u/ 0:6 di=ones(1,1); %untuk 6 diplus=ones(nxM-1,1); %diagonal utama dan atas u/7:10 diplusb=ones(nxM-2,1); %diagonal bwh u/7:10 for j=0:dx:6 q3=b*diminb; w3=-c*dimin; r3=b*dimin; k1=-g*V0(1:nxA-1,1); end for j=6


q2=2*b*di; w2=-f*di; r2=2*d*di; k2=-g*V0(nxA,1)-g*V0(nxA,1); end for j=6:dx:10 q1=d*diplus; w1=-e*diplus; r1=d*diplusb; k3=-g*V0(nxA+1:n-1,1); k4=-g*V0(n,1)-7*d; end q=[q3;q2;q1]; w=[w3;w2;w1]; r=[r3;r2;r1]; k=[k1;k2;k3;k4]; A=diag(q,-1)+diag(w)+diag(r,1); V=A\k; V(1,1)=0; V(n,1)=7; V; if t(i)==1 V05=V; elseif t(i)==5 V5=V; elseif t(i)==50 V50=V; end end plot(V05,x,'r+-', V5,x,'g.-',V50,x,'b--', v,x,'k-') legend('Numerical velocity at t=1', 'Numerical velocity at t=5', 'Numerical velocity at t=50', 'Exact steady velocity', 'Location', 'southeast') title('Kecepatan Fluida saat delta x = 0.25') xlabel('v (velocity)') ylabel('x (position)') error=sum(abs(V-v))/n

G. Metode Beda Hingga ∆𝑥 = 0.125 close all clear clc mua=1;%input('Kekentalan air: '); mum=3;%input('Kekentalan minyak: '); dx=0.125;%input('delta x: ') dt=0.5*dx; %input('delta t: ') k=6; t=0:dt:50; x=0:dx:10; xA=0:dx:6; nxA=length(xA); xM=6:dx:10;


nxM=length(xM); n=length(x);%input('masukkan besar matriks: '); %perhitungan Analitik v=zeros(n,1); for j=1:n if x(j)<=k v(j)=(7*mum*x(j))/(6*mum+4*mua); else v(j)=((7*mua*x(j))/(6*mum+4*mua))+((21*mum21*mua)/(3*mum+2*mua)); end end %perhitungan numerik V=zeros(n,1); g=1/dt; for i=1:length(t) V0=V; %nilai awal V0(n,1)=7; b=mua/(dx^2); c=(2*mua/(dx^2))+g; d=mum/(dx^2); e=(2*mum/(dx^2))+g; f=((2*mua+2*mum)/(dx^2))+(2*g); dimin=ones(nxA-1,1); %diagonal utama u/ 0:6 diminb=ones(nxA-2,1); %diagonal atas u/ 0:6 di=ones(1,1); %untuk 6 diplus=ones(nxM-1,1); %diagonal utama dan atas u/7:10 diplusb=ones(nxM-2,1); %diagonal bwh u/7:10 for j=0:dx:6 q3=b*diminb; w3=-c*dimin; r3=b*dimin; k1=-g*V0(1:nxA-1,1); end for j=6 q2=2*b*di; w2=-f*di; r2=2*d*di; k2=-g*V0(nxA,1)-g*V0(nxA,1); end for j=6:dx:10 q1=d*diplus; w1=-e*diplus; r1=d*diplusb; k3=-g*V0(nxA+1:n-1,1); k4=-g*V0(n,1)-7*d; end q=[q3;q2;q1]; w=[w3;w2;w1]; r=[r3;r2;r1]; k=[k1;k2;k3;k4];


A=diag(q,-1)+diag(w)+diag(r,1); V=A\k; V(1,1)=0; V(n,1)=7; V; if t(i)==1 V05=V; elseif t(i)==5 V5=V; elseif t(i)==50 V50=V; end end plot(V05,x,'r+-', V5,x,'g.-',V50,x,'b--', v,x,'k-') legend('Numerical velocity at t=1', 'Numerical velocity at t=5', 'Numerical velocity at t=50', 'Exact steady velocity', 'Location', 'southeast') title('Kecepatan Fluida saat delta x = 0.125') xlabel('v (velocity)') ylabel('x (position)') error=sum(abs(V-v))/n

H. Metode Beda Hingga ∆𝑥 = 0.0625 close all clear clc mua=1;%input('Kekentalan air: '); mum=3;%input('Kekentalan minyak: '); dx=0.0625;%input('delta x: ') dt=0.5*dx; %input('delta t: ') k=6; t=0:dt:50; x=0:dx:10; xA=0:dx:6; nxA=length(xA); xM=6:dx:10; nxM=length(xM); n=length(x);%input('masukkan besar matriks: '); %perhitungan Analitik v=zeros(n,1); for j=1:n if x(j)<=k v(j)=(7*mum*x(j))/(6*mum+4*mua); else v(j)=((7*mua*x(j))/(6*mum+4*mua))+((21*mum21*mua)/(3*mum+2*mua)); end end %perhitungan numerik V=zeros(n,1);


g=1/dt; for i=1:length(t) V0=V; %nilai awal V0(n,1)=7; b=mua/(dx^2); c=(2*mua/(dx^2))+g; d=mum/(dx^2); e=(2*mum/(dx^2))+g; f=((2*mua+2*mum)/(dx^2))+(2*g); dimin=ones(nxA-1,1); %diagonal utama u/ 0:6 diminb=ones(nxA-2,1); %diagonal atas u/ 0:6 di=ones(1,1); %untuk 6 diplus=ones(nxM-1,1); %diagonal utama dan atas u/7:10 diplusb=ones(nxM-2,1); %diagonal bwh u/7:10 for j=0:dx:6 q3=b*diminb; w3=-c*dimin; r3=b*dimin; k1=-g*V0(1:nxA-1,1); end for j=6 q2=2*b*di; w2=-f*di; r2=2*d*di; k2=-g*V0(nxA,1)-g*V0(nxA,1); end for j=6:dx:10 q1=d*diplus; w1=-e*diplus; r1=d*diplusb; k3=-g*V0(nxA+1:n-1,1); k4=-g*V0(n,1)-7*d; end q=[q3;q2;q1]; w=[w3;w2;w1]; r=[r3;r2;r1]; k=[k1;k2;k3;k4]; A=diag(q,-1)+diag(w)+diag(r,1); V=A\k; V(1,1)=0; V(n,1)=7; V; if t(i)==1 V05=V; elseif t(i)==5 V5=V; elseif t(i)==50 V50=V; end end plot(V05,x,'r+-', V5,x,'g.-',V50,x,'b--', v,x,'k-') legend('Numerical velocity at t=1', 'Numerical velocity at t=5', 'Numerical velocity at t=50', 'Exact steady velocity', 'Location', 'southeast')


title('Kecepatan Fluida saat delta x = 0.0625') xlabel('v (velocity)') ylabel('x (position)') error=sum(abs(V-v))/n

I. Gambar Galat Metode Beda Hingga clc clear all x=[0.0625,0.125,0.25,0.5,1,2]; Galat=[0.0187,0.0334,0.0547,0.0788,0.0964,0.0898]; plot(x,Galat,'r*-') legend('Galat Kecepatan Fluida','Location','southeast') xlabel('Steps in space domain') ylabel('Error') title('Galat Perhitungan Kecepatan Lapisan Fluida')

METODE NUMERIS UNTUK MENYELESAIKAN MODEL PERGERAKAN LAPISAN FLUIDA YANG MELIBATKAN MINYAK DAN AIR

Recommend Documents