PDP linear orde 2 Agus Yodi Gunawan Pada bagian ini akan dipelajari tiga jenis persamaan diferensial parsial (PDP) linear orde dua yang biasa dijumpai pada masalah-masalah dunia nyata, yaitu persamaan panas (persamaan difusi), persamaan Laplace, dan persamaan gelombang. Beberapa teknik penyelesaian akan disajikan untuk masing-masing persamaan tersebut, antara lain teknik pemisahan variabel, teknik transformasi integral, dan pengantar metode numerik beda hingga. Prasyarat: Masalah Sturm-Liouville Regular, Deret Fourier.
1
Klasifikasi
Bentuk umum PDP linear orde dua diberikan oleh ∂ 2u ∂2u ∂u ∂u ∂ 2u + B + C +D +E + F u = G, (1.1) 2 2 ∂x ∂x∂t ∂t ∂x ∂t dimana A, B, C, · · · , G merupakan fungsi konstan atau fungsi yang bergantung pada x dan t, u(x, t) disebut variabel terikat yang merupakan solusi dari PDP, x dan t disebut variabel bebas yang secara berurutan biasanya menyatakan varibel posisi dan variabel waktu. Namun secara umum arti dari variabel-variabel tersebut dapat disesuaikan dengan fenomena fisisnya. Pada (1.1) A, B, dan C diasumsikan tidak semuanya nol. Bentuk (1.1) dapat dituliskan dalam bentuk lain menjadi ( ) Auxx + Buxt + Cutt = R x, t, u, ux , ut , (1.2) A
dimana uxx = ∂ 2 u/∂x2 , uxt = ∂ 2 u/∂x∂t, dan utt = ∂ 2 u/∂t2 . Bentuk (1.2) berdasarkan tanda dari diskriman B 2 − 4AC dapat dikelompokan menjadi 2 Tipe Parabolik, B − 4AC = 0; (1.3) Tipe Eliptik, B 2 − 4AC < 0. Tipe Hiperbolik, B 2 − 4AC > 0; Dalam tabel berikut diberikan contoh-contoh PDP berserta tipenya;
1
AY
G
Bentuk Tipe Nama 2 ut = c uxx Parabolik Persamaan panas/Persamaan difusi uxx + uyy = 0 Eliptik Persamaan Laplace 2 utt = c uxx Hiperbolik Persamaan gelombang uxx + xuyy = 0 x < 0: Hiperbolik Persamaan Tricomi x > 0: Eliptik
Persamaan panas 1. Penurunan persamaan. Perhatikan suatu batang homogen dengan rapat massa tetap ρ dan penampang tetap A yang diletakkan pada sumbu x dari titik x = 0 sampai dengan x = L. Batang diasumsikan cukup tipis sedemikian sehingga aliran panas hanya dalam arah sumbu x saja dan tidak ada panas yang hilang sepanjang batang. Misalkan u(x, t) menyatakan suhu penampang batang pada posisi x setiap waktu t dan c menyatakan panas jenis batang (besarnya panas yang diperlukan untuk menaikkan suhu 1 derajat dari sebuah batang yang memiliki 1 satuan massa). Besarnya panas setiap saat yang terjadi diantara penampang batang pada posisi x dan penampang pada posisi x + △x adalah ∫ x+△x Q(t) = cρAu(s, t)ds. (2.1) x
Pada saat yang sama, laju aliran panas yang menembus penampang pada posisi x adalah sebanding dengan luas penampang dan gradien suhu pada penampang (Hukum Fourier untuk konduksi panas), yaitu −κAux (x, t),
(2.2)
dimana κ menyatakan konduktifitas panas batang. Tanda negatif pada (2.2) menyatakan bahwa panas mengalir ke arah suhu yang lebih rendah. hal yang sama juga untuk penampang pada posisi x + △x: −κAux (x + △x, t).
(2.3)
Selisih dari jumlah panas yang masuk pada penampang di x dan jumlah panas yang keluar dari penampang di x + △x harus sama dengan perubahan panas yang terjadi pada segmen x ≤ s ≤ x + △x. Dengan demikian kita peroleh [ ] ∫ x+△x ∂Q ∂u(x + △x, t) ∂u(x, t) ∂u(s, t) = cρA ds = κA − . (2.4) ∂t ∂t ∂x ∂x x
G
Asumsikan bahwa integran dari (2.1) adalah fungsi kontinu terhadap s, kita peroleh dengan menggunakan teorema nilai rata-rata untuk integral, ∫ x+△x ∂u(ξ, t) ∂u(s, t) ds = cρA △x, x < ξ < x + △x. (2.5) cρA ∂t ∂t x Selanjutnya, (2.4) menjadi
[ ] ∂u(ξ, t) ∂u(x + △x, t) ∂u(x, t) cρ △x = κ − . ∂t ∂x ∂x 2
AY
2
(2.6)
Dengan membagi kedua ruas (2.6) oleh cρ△x kemudian mengambil △x → 0, kita peroleh persamaan panas satu dimensi (disebut pula persamaan difusi) ut = α2 uxx ,
(2.7)
dengan tetapan difusifitas α2 = κ/(cρ). Catatan: • Jika terdapat panas tambahan dari luar diberikan pada batang dengan laju f (x, t) per satuan volum per satuan waktu, maka kita harus menambahkan ∫ x+△x pada turunan terhadap waktu Persamaan (2.4) suku x f (s, t)ds. Selanjutnya, jika proses yang sama dilakukan akan diperoleh ut = α2 uxx + F (x, t),
(2.8)
dengan F (x, t) = f (x, t)/(cρ). Persamaan (2.8) disebut persamaan panas tak homogen. Jika yang terjadi adalah kehilangan panas, maka Persamaan (2.8) masih tetap berlaku dengan menggantikan F (x, t) oleh −F (x, t). • Jika batang bergerak dalam arah sumbu x positif dengan kecepatan tetap v, maka kita harus menambahkan pada ruas kanan Persamaan (2.4) suku [vcρAu(x, t) − vcρAu(x + △x, t)]. Selanjutnya, jika proses yang sama dilakukan akan diperoleh ut = −vux + α2 uxx ,
(2.9)
Suku vux disebut suku konveksi sehingga Persamaan (2.9) biasa disebut Persamaan konveksi-difusi. 2. Kondisi batas. Terdapat tiga jenis kondisi batas: • Kondisi Dirichlet; pada setiap ujung selang nilai solusi diketahui, sebagai contoh u(0, t) = u0 , u(L, t) = uL .
G
• Kondisi Neumann; pada setiap ujung selang nilai dari turunan normal solusi diketahui, sebagai contoh ux (0, t) = u0 , ux (L, t) = uL . Jika kedua ujung batang diisolasi sehingga tidak ada panas yang keluar ataupun masuk, maka kondisi Neumann menjadi ux (0, t) = 0 = ux (L, t).
3
AY
• Kondisi Robin; pada setiap ujung selang diketahui nilai dari kombinasi linear dari dua tipe kondisi batas sebelumnya, yaitu ux (0, t) − hu(0, t) = A dan ux (L, t) + hu(L, t) = B dengan A, B, dan h konstanta positif. Kondisi ini biasa dipakai untuk menyatakan bahwa di ujung-ujung selang panas diradiasikan.
Untuk melengkapi persamaan panas, selain kondisi batas di atas perlu ditambahkan kondisi awal, yaitu u(x, 0) = g(x), (2.10) untuk suatu fungsi g(x) yang diketahui. 3. Pemisahan variabel. Teknik penyelesian yang umum dan sangat populer untuk masalah persamaan panas dalam daerah yang terbatas adalah pemisahan variabel. Gagasan dari teknik ini adalah menuliskan solusi u(x, t) menjadi perkalian dua buah fungsi yang masing-masing hanya bergantung pada satu peubah (x saja atau t saja). Berikut prosedur pemisahan variabel untuk menyelesaikan persamaan panas (2.7) dengan kondisi awal (2.10) dan kondisi batas Dirichlet; u(0, t) = 0 = u(L, t): I. Tuliskan u(x, t) = T (t)X(x). II. Substitusikan (I) ke (2.7), diperoleh T˙ (t)X(x) = α2 T (t)X ′′ (x), dimana T˙ (t) menyatakan turunan fungsi T (t) terhadap t dan X ′′ (x) menyatakan turunan kedua fungsi X(x) terhadap x. Dalam bentuk lain dapat dituliskan menjadi T˙ (t) X ′′ (x) = . α2 T (t) X(x)
(2.11)
III. Ruas kiri pada Kesamaan (2.11) hanya bergantung pada t sedangkan ruas kanannya hanya bergantung pada x. Kesamaan ini hanya dipenuhi oleh fungsi konstan. Misalkan, T˙ (t) X ′′ (x) = = −λ, (2.12) α2 T (t) X(x) dengan −λ suatu konstanta yang akan ditentukan kemudian. Tanda negatif dicantumkan untuk kemudahan saja. IV. Dari (2.12), kita peroleh dua persamaan T˙ (t) + α2 λT (t) = 0,
(2.13)
X ′′ (x) + λX(x) = 0, X(0) = 0 = X(L).
(2.14)
dan
1
AY
G
Persamaan (2.14) diselesaikan terlebih dahulu untuk menentukan konstanta λ. Persamaan (2.14) tidak lain merupakan masalah Sturm-Liouville Regular1 , yang dipenuhi oleh tak berhingga nilai eigen λn yang berpadanan dengan fungsi eigen Xn (x). Nilai dan fungsi eigen yang memenuhi (2.14) sangat bergantung pada kondisi batas yang diberikan. Silakan baca catatan saya: Masalah Sturm-Liouville Regular.
4
V. Misalkan sudah diperoleh nilai λn dan Xn (x) dari (2.14). Substitusikan nilai λn pada Persamaan (2.13) sehingga diperoleh T˙n (t) + α2 λn Tn (t) = 0, dengan solusi umumnya adalah Tn (t) = An exp(−α2 λn t). Nilai konstanta An akan ditentukan dari kondisi awal (2.10). VI. Solusi pada (I) dituliskan kembali menjadi un (x, t) = Tn (t)Xn (x). Persamaan panas (2.7) merupakan PDP linear. Oleh karena itu solusi umumnya adalah superposisi dari un (x, t), yaitu u(x, t) =
∞ ∑
un (x, t) =
n=0
∞ ∑
An Xn (x) exp(−α2 λn t).
(2.15)
n=0
Untuk t = 0, kita peroleh u(x, 0) = g(x) =
∞ ∑
An Xn (x).
(2.16)
n=0
Koefisien An pada (2.16) tidak lain adalah koefisien Fourier untuk fungsi g(x) dengan fungsi basisnya Xn (x). Setelah menentukan koefisien Fourier2 , solusi dari masalah persamaan panas diberikan oleh (2.15). 4. Teknik transformasi Laplace. Misalkan diberikan masalah berikut: ut = α2 uxx , a < x < b
(2.17)
dengan syarat awal u(x, 0) = u0 (x) dan syarat batas u(a, t) = u1 (t), u(b, t) = u2 (t). Misalkan pula transformasi Laplace untuk u(x, t) terhadap variabel t diberikan oleh ∫∞ L[u(x, t)] = U (x, s) =
u(x, t)e−st dt.
(2.18)
0
2
Silakan baca catatan saya: Deret Fourier.
5
AY
G
Dengan menerapkan (2.18) pada (2.17) akan diperoleh persamaan diferensial biasa orde dua 1 s (2.19) Uxx (x, s) − 2 U (x, s) = − 2 u0 (x) α α dengan kondisi batas U (a, s) = U1 (s), U (b, s) = U2 (s). Setelah menyelesaikan masalah syarat batas (2.19), u(x, t) dapat diperoleh dengan menghitung invers transformasi
Laplace U (x, s). Secara umum, menghitung invers dari transformsi Laplace memerlukan pengetahuan integral pada fungsi kompleks (di luar ruang lingkup kuliah ini). Pada bahasan kali ini akan ditampilkan contoh-contoh dimana invers transformasi Laplace nya sudah tersedia pada buku-buku yang memuat Tabel transformasi Laplace. Contoh. Persamaan panas pada daerah setengah bidang, x > 0 (diambil dari H.S. Carslaw & J.C. Jaeger, Conduction of heat in solids, Oxford Univ. Press, 1959 ). I. Misalkan temperatur awalnya nol, u(x, 0) = 0, dan syarat awalnya adalah u(0, t) = f (t). Persamaan (2.19) menjadi Uxx (x, s) − q 2 U (x, s) = 0, dimana q 2 = s/α2 dan kondisi batas U (0, s) = F (s). Solusi yang mensyaratkan agar solusi terbatas untuk x → ∞ dipenuhi oleh U (x, s) = F (s)e−qx . i. Jika u(0, t) = f (t) = C0 , konstan. Diperoleh U (x, s) =
C0 −qx e . s
Dari tabel transformasi diperoleh ( u(x, t) = C0 erfc
x √ 2α t
) ,
√ ∫u 2 dimana erfc(u) = 1 − erf(u) dengan erf(u) = (2/ π) 0 e−y dy. Fungsi erf (u) dan erfc (u) masing-masing secara berurutan disebut Error function and Complementary error function. ii. untuk sebarang syarat awal f (t), solusinya diberikan oleh x u(x, t) = √ 2α π
∫t 0
e−x /(4α (t−y)) f (t) dy. (t − y)3/2 2
2
G
II. Misalkan temperatur awalnya nol, u(x, 0) = 0 dan di x = 0 diberikan radiasi panas oleh suatu mediaum dengang temperatur tetap u0 . Kondisi batasnya diberikan oleh ux (0, t) = h(u(0, t) − u0 ),
6
AY
dengan h konstanta positif. Transformasi Laplace untuk syarat batas tersebut diberikan oleh hu0 . Ux (0, s) − hU (0, s) = − s
Solusinya diberikan oleh U (x, s) =
hu0 e−qx , s(q + h)
yang bersesuaian dengan ) ( ( ) √ x x hx+h2 α2 t √ − u0 e √ + hα t . u(x, t) = u0 erfc erfc 2α t 2α t 5. Metode beda hingga. Metode beda hingga adalah metode numerik yang memberikan nilai-nilai diskrit pada suatu kordinat (xj , tn ), yang disebut titik grid. Nilainilai numerik ini merupakan suatu nilai hampiran untuk solusi kontinu pada selang (xj − △x/2, xj + △x/2) dan (tn − △t/2, tn + △t/2), dimana △x dan △t masingmasing merupakan jarak diantara dua titik grid dalam arah sumbu x dan arah sumbu t. Untuk memudahkan penulisan, kita misalkan u(xj , tn ) = unj . Misalkan kita ingin mencari solusi numerik dari masalah berikut: ut = α2 uxx , (0 < x < L, t > 0);
(2.20)
dengan kondisi awal dan batas: u(x, 0) = g(x), u(0, t) = p(t), u(L, t) = q(t).
(2.21)
Kita bagi selang L menjadi N bagian yang masing-masing panjangnya sama, sehingga △x = L/N dan definisikan xj = j△x, j = 0, 1, 2, · · · , N . Hal yang sama untuk kordinat t, tn = n△t, n = 0, 1, 2, · · · . Tahap pertama dalam metode numerik adalah menghampiri solusi kontinu dengan hampiran beda hingganya. Alat yang akan digunakan adalah hampiran deret Taylor. Misalkan kita fokuskan untuk menghitung deret Taylor u(x, t) terhadap x di titik (xj , tn ), yaitu unj+1
=
unj
△x ∂unj (△x)2 ∂ 2 unj (△x)3 ∂ 3 unj + + + + ··· . 1! ∂x 2! ∂x2 3! ∂x3
(2.22)
Hal yang sama juga dapat diperoleh untuk unj−1 = unj −
(△x)2 ∂ 2 unj (△x)3 ∂ 3 unj △x ∂unj + − + ··· . 1! ∂x 2! ∂x2 3! ∂x3
(2.23)
Ada beberapa pilihan untuk menghampiri ∂unj /∂x: ∂unj unj+1 − unj ≈ ; disebut beda maju. ∂x △x ∂unj unj − unj−1 • dari (2.23): ≈ ; disebut beda mundur. ∂x △x ∂unj unj+1 − unj−1 • dari (2.22) dan (2.23): ≈ ; disebut beda pusat. ∂x 2△x 7
AY
G
• dari (2.22):
Untuk pendekatan dalam kordinat x, kita selanjutnya akan menggunakan beda pusat. Dengan cara yang sama hampiran beda pusat untuk turunan kedua terhadap x adalah ∂ 2 unj unj+1 − 2unj + unj−1 (2.24) ≈ ∂x2 (△x)2 Sedangkan untuk turunan pertama terhadap waktu kita akan gunakan beda maju. Dengan pendekatan beda hingga, (2.20) menjadi un+1 − unj un − 2unj + unj−1 j 2 j+1 =α △t (△x)2
(2.25)
= runj+1 + (1 − 2r)unj + runj−1 , un+1 j
(2.26)
atau dengan r = a2 △t/(△x)2 . Bentuk beda hingga (2.26) dikenal dengan bentuk FTCS (Forward Time-Centered Space). Bentuk (2.26) dilengkapi dengan kondisi awal u0j = g(xj ) = g(j△x) = gj (j = 1, 2, · · · , N − 1) dan kondisi batas un0 = p(tn ) = p(n△t) = pn , unN = q(tn ) = q(n△t) = qn (n = 1, 2, · · · ).
3
Persamaan Laplace
4
Persamaan gelombang
5
Latihan 1.
8
AY
G
Referensi: D.G. Duffy, Advanced Engineering Mathematics, CRC, 1998.