Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Solusi PERSAMAAN DIFERENSIAL PARSIAL (PDP) dengan HARGA AWAL dan KONDISI BATAS dalam PEMODELAN dan MODEL MATEMATIS
5
Bentuk umum :
Persamaan Diferensial Biasa (PDP) linier order 2 yang paling umum dalam 2 variabel bebas dapat dituliskan sebagai,
∂2 ∂2 ∂2 a 2 f(x, y) + b f(x, y) + c 2 f(x, y) ∂x ∂y ∂x ∂y
=
g( x , y )
dapat diklasifikasikan sebagai : 1. PDB ‘Hiperbolik’,
yaitu jika
b2 − 4 a c > 0
2. PDB ‘Parabolik’,
yaitu jika
b2 − 4 a c = 0
3. PDB ‘Eliptik’,
yaitu jika
b2 − 4 a c < 0
Beberapa Notasi Penulisan : 1. f x =
∂ f(x,y) ∂x
2. f y =
∂ f(x,y) ∂y
3. f xy =
∂2 f(x,y) ∂x ∂y
4. f xx =
∂2 f(x,y) ∂x 2
5. f yy
∂2 f(x,y) = ∂y 2
Beberapa contoh PDP yang terbentuk : 1. Persamaan Laplace (eliptik)
: f xx + f yy = 0
2. Persamaan gelombang (hiperbolik)
: f xx − c 2 f yy = 0
3. Persamaan difusi (parabolik)
: f x − D f yy = 0
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (1)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Notasi atau kemungkinan konfigurasi PDP lain :
• Persamaan Tricomi (mixed type) –
berbentuk hiperbolik, bila :
–
berbentuk parabolik, bila :
–
berbentuk eliptik, bila :
: y f xx + f yy = 0
y < 0 y = 0 y > 0
Perbedaan PDB dengan PDP : PDB
PDP
• Primitif persamaan merupakan persamaan aljabar dengan satu variabel bebas dan beberapa variabel tak bebas.
• Primitif persamaan berupa persamaan aljabar dengan beberapa variabel yang dapat saling berkaitan.
• Persamaan berbentuk suatu turunan hanya terhadap satu variabel bebas.
• Persamaan berupa turunan parsial terhadap variabel-variabelnya.
• Domain solusi atau penyelesaiannya berdimensi satu yang dapat sangat presisi : hanya bergantung pada pemilihan langkah (step) variabel bebas.
• Domain solusi atau penyelesaiannya berdimensi dua : bergantung pada pemilihan langkah (step) variabel-variabelnya.
• Solusi numerik berupa IVP.
• Berbentuk solusi IVP dan atau BVP.
Gambar 1. (a). BVP-IVP dan (b). BVP
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (2)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Strategi Solusi PDP : A.
Persamaan-persamaan yang berbentuk eliptik dapat diselesaikan dengan pemilihan beberapa metode yang sesuai : diskretisasi ‘Rectangular Mesh’ dengan diferensial hingga (finite differences), finite volume dan finite element.
B.
Solusi persamaan-persamaan yang berbentuk hiperbolik dapat berupa metode : diskretisasi ‘Rectangular Mesh’ dalam diferensial hingga (finite differences), persamaan-persamaan diferensial-hingga pendekatan (Skema Padé dan Skema Crank-Nicolson), bentuk reduksi menjadi sistem PDB, dan metode karakteristik.
C.
Persamaan-persamaan yang berbentuk parabolik dapat diselesaikan dengan pemilihan beberapa metode : diskretisasi dalam ‘grid segi-4’ atau pendekatan diferensial hingga eksplisit, MOL, metode Crank-Nicolson implisit, dan eliminasi Gauss (persamaan implisit).
BEBERAPA CONTOH JENIS SOLUSI PDP
Kasus #1 (Solusi PDP bentuk Parabolik) :
∂f ( x , t ) ∂2 f ( x , t ) = D , ∂t ∂x 2 D = konstanta
0 < x < 1, t > 0
Secara ringkas, beberapa metode praktis yang dapat digunakan untuk Kasus #1 ini adalah sbb : (a). MOL (Method Of Lines) h = x i +1 − x i , i = 1 ,K , N - ‘Grid’ pada arah-x : x 1 = 0 , x N +1 = 1
-
D du i = 2 [u i +1 − 2 u i + u i −1 ] Diskretisasi : dt h u ( t ) ≅ f ( x , t ) i
-
• Dirichlet : f ( 0 , t ) = g1 ( t ) Jenis kondisi batas : • Neumann : f x ( 1 , t ) = g 2 ( t ) • Robin : α f ( 0 , t ) + β f ( 0 , t ) = g ( t ) x 3
(b). Diferensial Hingga (Low-Order Time Approximation) h = x i +1 − x i , i = 1 ,K , N - ‘Grid’ pada arah-x : u 1 = 0 , u N +1 = 0
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (3)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
-
D du i dt = h 2 [u i +1 − 2 u i + u i −1 ], i = 2 ,K , N u i , j ≅ f ( x j , t j ), t j = j ∆t , j = 0 ,1 ,K Diskretisasi : sehingga u i , j +1 − u i , j D = 2 [u i +1 , j − 2 u i , j + u i −1 , j ] ∆t h
-
Implementasi : Eliminasi Gauss atau TDM.
(c). Metode Crank-Nicolson
-
∂ 2 f ( x ,t ) ∂f ( x , t ) = Pendekatan diskretisasi: 2 ∂t i , j + 12 ∂x i , j + 12 dengan u i , j +1 − u i , j k
=
1 2
u i +1 , j − 2 u i , j + u i −1 , j u i +1 , j +1 − 2 u i , j +1 + u i −1 , j +1 + 2 h h2
memberikan
− r u i −1 , j +1 + ( 2 + 2 r ) u i , j +1 − r u i +1 , j +1 = r u i −1 , j + ( 2 − 2 r ) u i , j + r u i +1 , j dengan
r = k/ h2 -
Secara umum, sisi kiri dari persamaan di atas mengandung 3 besaran tak diketahui, dan di sisi kanan 3 besaran diketahui, dari besaran u , seperti digambarkan di bawah ini :
Catatan : langkah waktu (time step) δt = k
haruslah sangat kecil, karena
proses tersebut hanya berlaku untuk 0 < k / h
h = δx harus dipertahankan kecil. Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
2
<
1 2
dan
Halaman (4)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Kasus # 2 (Solusi PDP bentuk Hiperbolik) : ∂f ( x , t ) ∂f ( x , t ) + a = 0, x > 0, t > 0 ∂t ∂x a = konstanta (positif dan nyata)
• kondisi awal (initial condition) :
f ( x ,0 ) = g ( x ),
x ≥ 0
• kondisi batas (boundary condition) : f ( 0 , t ) = b( t ), t > 0 Secara ringkas, salah satu metode praktis yang dapat digunakan untuk Kasus #2 ini adalah : (#). Backward Difference dan Sistem PDB :
-
Pendekatan ‘backward difference’ untuk turunan x (x-derivative) pada titik (x,t) dengan formula : ∂f f ( x ,t ) − f ( x − h,t ) = + O( h ) ∂x h
-
Dengan menganggap x tetap pada suatu posisi tertentu, PDP di atas dapat dituliskan sebagai suatu PDB : ∂f ( t ) a = − [f ( x , t ) − f ( x − h , t ) ] + O ( h ) ∂t h
-
Dengan menuliskan PDB di atas pada grid dengan N buah titik ( x i = i h , i = 1 ,K , N ), sepanjang bidang-waktu t, akan diperoleh suatu harga V i ( t ) yang merupakan harga pendekatan dari f i ( t ) yang berjumlah N buah persamaan : dV 1 a = − [V 1 − b( t ) ] dt h dV 2 a = − [V 2 − V 1 ] dt h M dV N a = − [V N − V N −1 ] dt h
Kasus #3 (Solusi PDP bentuk Eliptik) : ∂2 f ( x , y ) ∂2 f ( x , y ) + = 0, ∂x 2 ∂y 2
0 ≤ x ≤ 0, 0 ≤ y ≤ 0
Secara ringkas, salah satu metode praktis yang dapat digunakan untuk Kasus #1 ini adalah : Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (5)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
(#). Diskretisasi dalam grid segi-4 :
-
Pendekatan ‘grid bujur-sangkar’ sehingga diperoleh :
∆x = ∆y = h -
Jika Nh = 1 , makam jumlah ‘titik grid internal’ adalah ( N − 1 ) 2
-
Maka diskretisasi diferensial hingga untuk persamaan order-2 di atas adalah : 1 ( ∆x )2
[u i +1 , j − 2 u i , j + u i −1 , j ] +
1 ( ∆y )2
[u i , j +1 − 2 u i , j + u i , j −1 ] = 0
dengan
u i,j ≅ f ( x i , y j ) xi = i h yj = jh -
Karena ∆x = ∆y , maka persamaan diskretisasi di atas dapat dituliskan sebagai :
u i , j −1 + u i +1 , j − 4 u i , j + u i −1 , j + u i , j +1 = 0 dengan ketelitian atau akurasi sebesar O ( h 2 ) .
-
Membentuk SISTEM PERSAMAAN ALJABAR LINIER (SPAL).
-
Penyelesaian akhirnya akan sangat bergantung pada jenis kondisi batasnya : Dirichlet Problem, Neumann Problem atau Robin Problem. BEBERAPA CONTOH PENYELESAIAN PDP
Contoh #1 (Penyelesaian PDP bentuk Parabolik) : Selesaikan PDP di bawah ini :
∂U ∂ 2U = ∂t ∂x 2 •
yang memenuhi kondisi awal berikut : U =1
•
untuk
0 ≤ x ≤1
pada saat
t =0
dan kondisi batas berikut :
∂U = U pada x = 0 , untuk setiap t ∂x ∂U = −U pada x = 1 , untuk setiap t ∂x menggunakan ‘metode diferensial-hingga eksplisit’ mengimplementasikan ‘central-differences’ untuk kondisi batasnya. Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
dan Halaman (6)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Jawab : Representasi ‘diferensial-hingga’ dari PDP di atas adalah : u i , j +1 − u i , j
=
δt
u i −1 , j − 2 u i , j + u i +1 , j ( δx ) 2
,
x a i dan t a j ’
atau bentuk eksplisitnya :
u i , j +1 = u i , j + ρ (u i −1 , j − 2 u i , j + u i +1 , j )
dengan ρ =
(1)
δt ( δx )2
Pada saat x = 0 ,
u 0 , j +1 = u 0 , j + ρ (u −1 , j − 2 u 0 , j + u 1 , j )
(2)
Penulisan formulasi kondisi batas pada x = 0 dalam term ‘central-differences’ adalah sebagai berikut : u 1 , j − u −1 , j 2 δx
= u 0,j
(3)
Penyisihan ‘nilai imejiner’
u −1 , j
dalam persamaan (2) dengan (3) di atas :
u 0 , j +1 = u 0 , j + 2 ρ (u 1 , j − [1 + δx ]u 0 , j )
(4)
Jika diambil δx = 0 ,1 . Sehingga, pada saat x = 1 , persamaan (1) di atas menjadi :
u 10 , j +1 = u 10 , j + ρ (u 9 , j − 2 u 10 , j + u 11 , j )
(5)
dan formula kondisi batasnya adalah : u 11 , j − u 9 , j 2 δx
= − u 10 , j
Penyisihan ‘nilai fiktif’
(6)
u 11 , j
pada persamaan (5) dengan bantuan persamaan
(6) di atas, diperoleh :
u 10 , j +1 = u 10 , j + 2 ρ (u 9 , j − [1 + δx ] u10 , j )
(7)
Dari persamaan (4) dan (7) di atas, dapat disimpulkan bahwa terdapat suatu ‘simetri’ pada saat x = 12 . Salah satu syarat yang perlu diketahui, dalam skema solusi ini adalah bahwa :
r≤
1 2 + h δx
Jika dipilih r = ¼, maka persamaan-persamaan (4) dan (1) di atas dapat dituliskan sebagai : Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (7)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
u 0 ,1 =
1 2
u i , j +1 =
(0 ,9 u 1 (u 4
0,j
i −1 , j
+ u 1 , j ), + 2 u i , j + u i +1 , j )
(i = 1 ,2 ,3 ,4 ),
dan pada titik simetri, dapat dituliskan (tanpa perlu menuliskan formula penuh seprti persamaan (7)) :
u 5 , j +1 =
1 4
2 u 4 , j + 2 u 5 , j
Karena harga awal adalah 2 pertama bila t = ρ (δx ) =
u 0 ,1 = u 1 ,1 =
(0 ,9 + 1 ) = 0 ,95 1 (1 + 2 + 1 ) = 1 = 4
u =1, 1 400
maka nilai-nilai , adalah :
u
pada akhir ‘time-step’
1 2
u 2 ,1 = u 3 ,1 = u 4 ,1 = u 5 ,1 ,
dan nilai-nilai pada akhir ‘time-step’ kedua adalah :
u 0 ,2 = u 1 ,2 = u 2 ,2 =
(0 ,9 × 0 ,95 + 1 ) = 0 ,9275 1 (0 ,95 + 2 + 1 ) = 0 ,9875 4 1 (1 + 2 + 1 ) = 1 = u 3 ,2 = u 4 ,2 4 1 2
= u 5 ,2 .
Dan seterusnya, dan hasil-hasilnya disajikan sebagai berikut :
t = 0,0000 0,0025 t = 0,0050 0,0075 0,0100 0,0125 0,0150 0,0175 0,0200 … … 0,1000 0,2500 0,5000 1,0000
i=0 x=0
1 0,1
2 0,2
3 0,3
4 0,4
5 0,5
1,0000 0,9500 0,9275 0,9111 0,8978 0,8864 0,8764 0,8673 0,8590
1,0000 1,0000 0,9875 0,9756 0,9648 0,9549 0,9459 0,9375 0,9296
1,0000 1,0000 1,0000 0,9969 0,9923 0,9872 0,9818 0,9762 0,9708
1,0000 1,0000 1,0000 1,0000 0,9992 0,9977 0,9956 0,9931 0,9902
1,0000 1,0000 1,0000 1,0000 1,0000 0,9998 0,9993 0,9985 0,9974
1,0000 1,0000 1,0000 1,0000 1,0000 1,0000 0,9999 0,9996 0,9991
0,7175 0,5542 0,3612 0,1534
0,7829 0,6048 0,3942 0,1674
0,8345 0,6452 0,4205 0,1786
0,8718 0,6745 0,4396 0,1867
0,8942 0,6923 0,4512 0,1917
0,9017 0,6983 0,4551 0,1933
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (8)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Representasi ‘Grid’ atau ‘Mesh’ Segi-4 dari PDP Parabolik di atas :
j (0,0)
(0,1)
(0,2)
(0,3)
(0,4)
(0,5)
(0,6)
(0,7)
(0,8)
(0,9)
(0,10 )
(1,0)
(1,1)
(1,2)
(1,3)
(1,4)
(1,5)
(1,6)
(1,7)
(1,8)
(1,9)
(1,10 )
∆t
(2,0)
(2,1)
(2,2)
(2,3)
(2,4)
(2,5)
(2,6)
(2,7)
(2,8)
(2,9)
(2,10 )
∆t
(N,0) (N,1) (N,2) (N,3) (N,4) (N,5) (N,6) (N,7) (N,8) (N,9)
(N,10 )
i
∆x
∆x
∆x
Bidang Simetri
Solusi Analitis dari PDP Parabolik dengan kondisi-kondisi batas dan kondisi awal seperti di atas adalah : U = 4
dengan
αn
sec α n −4 α n2 t e cos 2 α n (x − 12 ) 2 n =1 3 + 4 α n ∞
∑ (
)
(0 < x < 1 ),
adalah akar-akar positif dari persamaan berikut :
α tan α =
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
1 2
Halaman (9)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Harga-harga U dari hasil solusi analitis seperti di atas dapat ditabelkan sebagai berikut : x=0
x = 0,1
x = 0,2
x = 0,3
x = 0,4
x = 0,5
0,9460 0,9250 0,9093 0,8965 0,8854 0,8755 0,8666 0,8585
0,9951 0,9841 0,9730 0,9627 0,9532 0,9444 0,9362 0,9286
0,9999 0,9984 0,9950 0,9905 0,9855 0,9802 0,9748 0,9695
1,0000 0,9999 0,9994 0,9984 0,9967 0,9945 0,9919 0,9891
1,0000 1,0000 1,0000 0,9998 0,9994 0,9988 0,9979 0,9967
1,0000 1,0000 1,0000 1,0000 0,9999 0,9996 0,9992 0,9985
0,7176 0,5546 0,3619 0,1542
0,7828 0,6052 0,3949 0,1682
0,8342 0,6454 0,4212 0,1794
0,8713 0,6747 0,4403 0,1875
0,8936 0,6924 0,4519 0,1925
0,9010 0,6984 0,4558 0,1941
t= 0,0025 0,0050 0,0075 0,0100 0,0125 0,0150 0,0175 0,0200 … … 0,1000 0,2500 0,5000 1,0000
Harga-harga U pada saat x = 0,2 antara hasil solusi dengan diferensialhingga (DH) dan solusi analitis dapat dibandingkan seperti ditabelkan di bawah ini : Solusi DH (pada x = 0,2)
Solusi Analitis (pada x = 0,2)
Persentase Galat
1,0000 0,9126 0,8345 0,6452 0,4205 0,1786
0,9984 0,9120 0,8342 0,6454 0,4212 0,1794
0,16 0,07 0,04 -0,03 -0,16 -0,45
t = 0,005 0,050 0,100 0,250 0,500 1,000
Contoh #2 (Penyelesaian PDP bentuk Parabolik) : Selesaikan PDP di bawah ini (sama seperti soal sebelumnya) :
∂U ∂ 2U = ∂t ∂x 2 •
yang memenuhi kondisi awal berikut : U =1
•
untuk
0 ≤ x ≤1
pada saat
t =0
dan kondisi batas berikut :
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (10)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
∂U =U ∂x ∂U = −U ∂x
pada pada
x = 0 , untuk setiap t x = 1 , untuk setiap t
menggunakan ‘metode Crank-Nicholson’ dan aplikasi ‘central-differences’ pada kondisi-kondisi batasnya. Jawab : Representasi ‘diferensial-hingga’ dari PDP di atas adalah : u i , j +1 − u i , j δt
=
1 u i +1 , j +1 − 2 u i , j +1 + u i −1 , j +1 u i +1 , j − 2 u i , j + u i −1 , j + 2 ( δx ) 2 ( δx ) 2 ’
x a i dan t a j
dapat dituliskan kembali sebagai :
− ρ u i −1 , j +1 + (2 + 2 ρ )u i , j +1 − ρ u i +1 , j +1 = ρ u i −1 , j + (2 − 2 ρ )u i , j + ρ u i +1 , j dengan : ρ =
(1)
δt
(δx )2
Representasi ‘central difference’ dari kondisi batas pada x = 0 : u 1 , j − u −1 , j 2
δx
= u 0,j
sehingga kedua persamaan berikut memenuhi kondisi di atas,
u −1 , j = u 1 , j − 2 δx ⋅ u 0 , j dan
u −1 , j +1 = u 1 , j +1 − 2 δx ⋅ u 0 , j +1 .
Kedua persamaan terakhir di atas dapat digunakan untuk mengganti termterm u −1 , j dan u −1 , j +1 pada persamaan (1), pada saat i = 0 . Dengan cara yang sama, kondisi batas pada x = 1 dapat diturun-kan seperti di atas, namun akan lebih mudah untuk menggunakan ‘prinsip simetri’ pada x = 12 , dalam hal ini : u 6,j = u 4,j Skema penyelesaian metode ini secara umum berlaku untuk sembarang nilai nyata ρ , namun demikian akan lebih baik jika harga ρ ini dipertahankan masih cukup kecil agar supaya solusinya lebih dekat pada solusi analitis PDP. Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (11)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Pilih ρ = 1 dan δx = 0,1 . Maka akan diperoleh suatu SPAL yang representatif untuk ‘nilai-nilai pivotal’ di atas, yaitu yang diwakili oleh besaran-besaran u 0 , j +1 , u 1 , j +1 , K , u 5 , j +1 , sebagai berikut : 2 ,1
u 0 , j +1 − u 1 , j +1 = − 0 ,1 u 0 , j + u 1 , j ,
− u i −1 , j +1 + 4 u i , j +1 − u i +1 , j +1 = u i −1 , j + u i +1 , j
(i = 1 , 2 , 3 , 4 )
− u 4 , j +1 + 2 u 5 , j +1 = u 4 , j Sehingga pada ‘time-step’ yang pertama akan diperoleh SPAL : 2 ,1 u 0 − u 1 = 0 ,9 − u 0 + 4 u 1 − u 2 = 2 ,0 − u 1 + 4 u 2 − u 3 = 2 ,0 − u 2 + 4 u 3 − u 4 = 2 ,0 − u 3 + 4 u 4 − u 5 = 2 ,0 − u 4 + 2 u 5 = 1 ,0
Bentuk SPAL yang didapat adalah (Bentuk TDM) :
0 0 0 u 0 0 ,9 2 ,1 − 1 0 2 ,0 − 1 4 − 1 0 0 0 u 1 2 ,0 0 −1 4 −1 0 0 u 2 = ⋅ 0 − 1 4 − 1 0 u 3 2 ,0 0 2 ,0 0 0 0 − 1 4 − 1 u 4 u 0 0 0 0 1 2 − 1 ,0 5 Hasil dari penyelesaian SPAL untuk beberapa ‘time-step’ disajikan pada tabel berikut : i=0
i=1
i=2
i=3
i=4
i=5
1,0000 0,8908 0,8624 0,7179 0,5547 0,3618 0,1540
1,0000 0,9707 0,9293 0,7834 0,6054 0,3949 0,1680
1,0000 0,9922 0,9720 0,8349 0,6458 0,4212 0,1793
1,0000 0,9979 0,9900 0,8720 0,6751 0,4404 0,1874
1,0000 0,9994 0,9964 0,8944 0,6929 0,4520 0,1923
1,0000 0,9997 0,9979 0,9018 0,6989 0,4559 0,1940
t= 0,00 0,01 0,02 0,10 0,25 0,50 1,00
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (12)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Harga-harga U pada saat x = 0,2 antara hasil solusi dengan metode CrankNicolson (C-N) dan solusi analitis bila dibandingkan akan diperoleh data seperti ditabelkan di bawah ini : Solusi C-N (pada x = 0,2)
Solusi Analitis (pada x = 0,2)
Persentase Galat
0,9922 0,9131 0,8349 0,6458 0,4212 0,1793
0,9905 0,9120 0,8342 0,6454 0,4212 0,1794
0,17 0,12 0,08 0,06 0,00 -0,06
t = 0,01 0,05 0,10 0,25 0,50 1,00
Contoh #3 (Penyelesaian PDP bentuk Parabolik) : Selesaikan PDP berikut (masih sama seperti soal sebelumnya) :
∂U ∂ 2U = D ∂t ∂x 2 •
dengan kondisi awal berikut : U =1
•
untuk
0 ≤ x ≤1
pada saat
t =0
dan kondisi batas berikut :
∂U = U pada x = 0 , untuk setiap t ∂x ∂U = −U pada x = 1 , untuk setiap t ∂x menggunakan ‘Method Of Lines (MOL)’ dan aplikasi ‘central-differences’, bila D = 1,0 , pada kondisi-kondisi batasnya. Jawab : Representasi bentuk MOL (solusi dalam ‘sistem PDB’) dengan formula ‘diferensial-hingga’, pada ‘time step j’ dari PDP di atas :
du i , j dt
=
D ( ∆x )2
[u
i −1 , j
− 2 u i , j + u i +1 , j ],
x a i dan t a j
atau dalam formula lebih sederhana, untuk sembarang du i = dt
•
1 ( ∆x )2
[u i −1 − 2 u i
+ u i +1 ]
j
: (1)
Penerapan kondisi awal : Pada saat t = 0 , atau j = 0 , semua harga u = 1
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (13)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
•
Penerapan kondisi batas : Representasi ‘central difference’ pada x = 0 , atau i = 0 :
u 1 − u −1
= u0 ,
2 ∆x
sehingga diperoleh PDB pada kondisi batas ini : du 0 = dt
2 ( ∆x )2
[− (1 +∆x ) u 0 + u 1 ]
Jika diambil ∆x = 0 ,1 , maka sama seperti pada problem solusi sebelumnya, yaitu terdapatnya ‘bidang simetri’ pada x = 0 ,5 , atau i = 5 , dalam hal ini :
∂u ∂x
u 6,j − u 4,j
=
2 ∆x
x = 0 ,5
= 0
sehingga
u 6,j = u 4,j sehingga diperoleh PDB pada ‘bidang simetri ini’ adalah: du 5 = dt
2 ( ∆x )2
[u 4 −
u5 ]
Secara keseluruhan keenam PDB yang terbentuk dapat dituliskan sebagai berikut :
du 0 dt du 1 dt du 2 dt du 3 dt du 4 dt du 5 dt
= = = = = =
2 ( ∆x ) 2 1 ( ∆x ) 2 1 ( ∆x ) 2 1 ( ∆x ) 2 1 ( ∆x ) 2 2 ( ∆x ) 2
[− (1 + ∆x ) u 0 + u1 ] [u 0 − 2 u1 + u 2 ] [u1 − 2 u 2 + u 3 ] [u 2 − 2 u 3 + u 4 ] [u 3 − 2 u 4 + u 5 ] [u 4 −
u5 ]
Bentuk sistem PDB di atas dapat diselesaikan dengan metode Runge-Kutta yang telah dipelajari (termasuk RKG, RKM, RKF45).
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (14)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
Hasil integrasi dengan metode RKG ( h = 0,00001 ) disajikan sebagai berikut :
t = 0,0000 0,0025 0,0050 0,0075 0,0100 0,0125 0,0150 0,0175 0,0200 … … 0,1000 0,2500 0,5000 1,0000
i=0 x=0
1 0,1
2 0,2
3 0,3
4 0,4
5 0,5
1,0000 0.9614 0.9358 0.9173 0.9038 0.8926 0.8828 0.8733 0.8645
1,0000 0.9956 0.9867 0.9765 0.9674 0.9587 0.9505 0.9421 0.9341
1,0000 0.9997 0.9980 0.9950 0.9914 0.9873 0.9828 0.9778 0.9725
1,0000 1.0000 0.9998 0.9991 0.9981 0.9967 0.9949 0.9926 0.9898
1,0000 1.0000 1.0000 0.9999 0.9996 0.9990 0.9982 0.9970 0.9952
1,0000 1.0000 1.0000 0.9997 0.9993 0.9984 0.9972 0.9953 0.9928
0.7069 0.5077 0.3090 0.1145
0.7697 0.5528 0.3365 0.1247
0.8156 0.5858 0.3566 0.1321
0.8434 0.6058 0.3688 0.1367
0.8527 0.6125 0.3729 0.1382
0.8433 0.6057 0.3687 0.1366
Harga-harga U pada saat x = 0,2 antara hasil solusi dengan ‘Method Of Lines’ dan solusi analitis bila dibandingkan akan diperoleh data seperti ditabelkan di bawah ini :
t = 0,01 0,05 0,10 0,25 0,50 1,00
Solusi MOL (pada x = 0,2)
Solusi Analitis (pada x = 0,2)
Persentase Galat
0,9914 0,9105 0,8156 0,5858 0,3566 0,1321
0,9905 0,9120 0,8342 0,6454 0,4212 0,1794
0,09 -0,16 -2,28 -10,17 -18,11 -35,80
Contoh #4 (Penyelesaian PDP bentuk Eliptik) : Suatu pelat logam bujur-sangkar R = {(x , y ) : 0 ≤ x ≤ 1 , 0 ≤ y ≤ 1} dipanaskan kedua sisinya, dengan suhu berbeda) sampai tercapai kondisi tunak. Pada keadaan ini persamaan konduksi panas tersebut dapat dinyatakan sebagai : ∂ 2T ∂ 2T = 0 + ∂y 2 ∂x 2 Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
Halaman (15)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
•
dengan kondisi-kondisi batas berikut : T (0 , y ) = T1
pada sisi pemanasan (suhu tetap)
T (1 , y ) = T 2
pada sisi lawan/seberangnya (suhu tetap)
∂T (x ,0 ) = 0 (permukaan yang diisolasi secara sempurna) ∂y ∂T (x ,1 ) = k [T (x ,1 ) − T 2 ] (panas terkonveksi sepanjang y = 1 ) ∂y dengan T1 , T 2 , dan k adalah tetap dan T 1 ≥ T (x , y ) ≥ T 2 . Selesaikan PDP eliptik di atas, untuk menghitung perkembangan suhu atau pemanasan permukaan sepanjang sumbu x dan y , dengan memperhatikan kondisi-kondisi batasnya. Jawab : Representasi ‘grid’ dari PDP perpindahan panas pada pelat logam di atas dapat digambarkan dalam suatu domain R , sedemikian rupa sehingga, x i = i h y j = j h (dalam hal ini ∆x = ∆y ) dan N h = 1 . Dalam hal ini, untuk setiap tempat kedudukan ‘titik interior’ dalam grid berlaku persamaan :
u i , j −1 + u i +1 , j − 4 u i , j + u i −1 , j + u i , j +1 = 0 dengan
•
(1)
u i , j ≈ T (x i , y j )
Pada daerah-daerah batas x = 0 dan x = 1 berlaku kondisi-kondisi batas menurut Dirichlet, yaitu :
•
u 0 , j = T1 ,
untuk
j = 0 ,1 ,K , N − 1
(2)
u N ,j = T 2 ,
untuk
j = 0 ,1 ,K , N − 1
(3)
Pada daerah y = 0 , permukaan yang diisolasi sempurna mengarahkan pada kondisi Neumann, yang berarti dapat didiskretisasi dengan ‘central difference’ sebagai berikut :
u i , −1 − u i ,1 = 0 , •
untuk
i = 0 ,1 ,K , N − 1
(4)
Dan pada daerah y = 1 , merupakan kondisi batas menurut Robin, yaitu : u i ,N −1 − u i ,N +1 2
h
= k [u i ,N − T 2 ], untuk i = 0 ,1 ,K , N − 1
Property of Setijo Bismo: Seri Perkuliahan Pemodelan dan Mater - Solusi PDP
(5)
Halaman (16)
Seri Mata Kuliah : PEMODELAN dan MATEMATIKA TERAPAN
y
∂T = k (T − T 2 ) ∂y
0,4
0,3 T = T1
0,2
0,1
1,4
2,4
3,4
4,4
1,3
2,3
3,3
4,3
1,2
2,2
3,2
4,2
1,1
2,1
3,1
4,1
0,0
1,0
2,0 ∂T = 0 ∂y
3,0
T = T2
4,0
x
Persamaan-persamaan (1), (2), (3), (4), dan (5) di atas, jika disusun sebagai solusi untuk titik-titik interior di dalam representasi grid di atas maka akan diperoleh SPAL sebagai berikut : − 4 1 0 1
1 0 2 −4 1 0 1 −4 0 0 1
0 −4 0 1 1 0 1
2 0 1 −4 1 0 1
2 0 1 −4 0 0 1
1 0 0 −4 1 0 1
1 0 1 −4 1 0 1
1 0 1 −4 0 0 1
1 0 0 −4 1 0 2
1 0 1 −4 1 0 2
1 0 1 −4 0 0 2
− T1 u 1 ,0 u 0 2 ,0 u 3 ,0 −T 2 − T1 u 1 ,1 u 2 ,1 0 − T2 u 3 ,1 u − T 1 1 ,2 ⋅ u 2 ,2 = 0 −T 2 u 3 ,2 u 1 ,3 − T1 1 0 0 1 0 u 2 ,3 u − T 0 0 1 2 3 ,3 1 0 ( 4 +2 h k ) − (T1 + 2 h k T 2 ) u 1 , 4 u 2h kT −( 4 + 2 h k ) 1 1 2 2 ,4 − (T1 + 2 h k T 2 ) − ( 4 + 2 h k ) u 3 ,4 0 1
⇒ Periksalah jumlah elemen matriks yang berharga 0 !
Property of Setijo Bismo