PEMODELAN DAN SOLUSI NUMERIK ALIRAN GAS DALAM SALURAN PIPA MENGGUNAKAN METODE CRANK-NICOLSON (MODELING AND NUMERICAL SOLUTION OF A GAS FLOW IN PIPELINE USING CRANK-NICOLSON METHOD) Oleh: Zusnita Meyrawati (1205 100 071) Pembimbing: 1. Prof. DR. Basuki Widodo, M.Sc. 2. Drs. Kamiran, M.Si. Jurusan Matematika Fakultas Matematika dan Ilmu Pengetahuan Alam Institut Teknologi Sepuluh Nopember Surabaya 2010 Abstrak Jaringan pipa adalah sistem dengan rangkaian pipa yang panjang. Panjang dari jaringan pipa ini bisa mencapai ribuan kilometer. Sedangkan distribusi gas melalui saluran pipa bisa memiliki arti menghubungkan pipa-pipa tersebut dari pusat produksi sampai ke tempat-tempat penyimpanan gas. Oleh karena itu, pada sistem jaringan pipa distribusi gas dilengkapi dengan stasiun kompresi beserta peralatan lain, seperti katup dan regulator. Dengan demikian sistem ini bekerja pada tekanan tinggi dan menggunakan stasiun kompresi, untuk menyediakan energi yang cukup bagi gas tersebut, agar bisa mengalir sepanjang jarak yang akan ditempuh. Ketika gas mengalir, gas akan mengalami kehilangan energi dan tekanan akibat gesekan antara gas dan dinding dalam pipa, dan juga disebabkan transfer panas antara gas dan lingkungan di sekitar pipa. Sedangkan saluran pipa adalah komponen terpenting dari jaringan pipa distribusi gas. Hal itu dikarenakan saluran pipa tersebut merupakan penentu karakteristik dinamik utama. Karena itu perlu diketahui perilaku gas dalam sistem jaringan distribusi ini. Pada Tugas Akhir ini dikaji tentang model aliran gas dalam suatu saluran pipa untuk mengetahui perilaku gas dalam sistem jaringan distribusi. Untuk itu dibuat asumsi dan batasan masalah serta digunakan persamaan konservasi massa, persamaan konservasi momentum, dan persamaan keadaan gas nyata. Model yang diperoleh berupa suatu sistem persamaan diferensial parsial (PDP) jenis parabolik linier orde dua. Model tersebut selanjutnya diselesaikan secara numerik menggunakan metode Crank-Nicolson dan disimulasikan dengan bantuan program Matlab 7.6. Hasil dari simulasi menunjukkan bahwa jika diberikan tekanan 50 bar di inlet, ketinggian di ujung pipa 3000 m, dan laju alir sebesar 50 m3/s, serta dengan menggunakan nilai parameterparameter lain yang telah ditentukan, maka tekanan tersebut akan menurun secara periodik. Hal tersebut terlihat dari nilai tekanan di outlet sebesar 33,485 bar. Nilai tekanan di outlet tersebut akan berubah jika nilai parameter diubah. Dari hasil simulasi beberapa contoh kasus, terlihat bahwa nilai tekanan di outlet dipengaruhi oleh nilai dari parameter-parameter, diantaranya laju alir dan juga ketinggian di ujung pipa. Kata kunci: aliran gas, saluran pipa, metode Crank-Nicolson. Peningkatan ini didukung oleh beberapa fakta, diantaranya, gas lebih bersih daripada sumber energi lain, gas relatif lebih murah, terutama jika dibandingkan dengan minyak atau batu bara, dan yang utama karena Indonesia mempunyai cadangan gas yang melimpah. Distribusi dan transportasi gas telah sangat umum dilakukan di banyak negara di
I.
Pendahuluan Gas merupakan salah satu sumber energi alternatif yang layak diperhitungkan, mengingat kenyataan bahwa cadangan minyak dunia saat ini telah menipis. Di Indonesia penggunaan sumber energi alternatif ini meningkat sejalan dengan perkembangan industri yang terjadi di berbagai daerah. 1
Pada Tugas Akhir ini, model aliran gas dalam suatu jaringan pipa gas dibuat dengan menitikberatkan pada saluran pipa. Hal ini karena saluran pipa adalah komponen terpenting dari sistem tersebut. Model diturunkan dari persamaan konservasi massa, persamaan konservasi momentum, dan persamaan keadaan gas nyata, kemudian model tersebut diselesaikan menggunakan metode Crank-Nicolson. Dan selanjutnya hasil penyelesaian dari model ini divisualisasikan menggunakan bantuan program Matlab 7.6.
dunia, terutama negara yang mempunyai sumber-sumber gas. Di Indonesia sendiri, pemanfaatan gas, yaitu gas alam, dimulai pada tahun 1960-an dimana produksi gas alam dari ladang gas alam PT. Stanvac Indonesia di Pendopo, Sumatera Selatan dikirim melalui saluran pipa ke pabrik pupuk Pusri IA, PT. Pupuk Sriwidjaja di Palembang. Perkembangan pemanfaatan gas alam di Indonesia meningkat pesat sejak tahun 1974. Gas alam sendiri adalah bahan bakar fosil berbentuk gas, dengan komponen utamanya adalah metana (CH4), yang merupakan molekul hidrokarbon rantai terpendek dan teringan. Gas alam lebih ringan dari udara sehingga cenderung mudah tersebar di atmosfer. Namun gas alam dapat berbahaya karena sifatnya yang sangat mudah terbakar dan menimbulkan ledakan. Apabila berada di dalam ruang tertutup, seperti di dalam pipa, konsentrasi gas dapat mencapai titik campuran yang mudah meledak. Sehingga cukup sulit untuk menyimpan gas alam karena hal ini sangat mahal dan berbahaya. Salah satu cara untuk mengatasi masalah tersebut adalah dengan mendapatkan konsumen terlebih dahulu sebelum mendistribusikan gas alam. Pengiriman gas alam dari daerah produksi ke konsumen dapat dilakukan dengan beberapa cara tergantung situasi dan kondisi, antara lain dengan sistem transmisi pipa atau dikonversi dahulu ke bentuk lain seperti Compressed Natural Gas (CNG) atau Liquid Natural Gas (LNG). Apabila pengiriman gas alam dilakukan dengan sistem transmisi pipa, maka perlu diperhatikan perilaku aliran gas alam dalam saluran pipa selama proses transportasi dari tempat produksi ke konsumen agar pengiriman gas alam dapat dilakukan dengan benar dan aman. Pada penelitian-penelitian sebelumnya, telah banyak dilakukan pemodelan pada aliran gas dalam saluran pipa yang dilakukan dengan mengasumsikan bahwa aliran gas dalam keadaan steady, dimana kondisi tekanan, laju alir, dan temperatur tidak berubah terhadap waktu. Pada beberapa kondisi, asumsi ini memberikan hasil yang cukup baik. Namun kenyataan di lapangan menunjukkan bahwa terdapat banyak situasi, dimana aliran gas dalam keadaan unsteady. Oleh karena itu perlu dibangun suatu model matematika yang merepresentasikan aliran gas dalam saluran pipa yang unsteady (transient).
II. 2.1
Pada tahun 1995, Junyang Zhou dan Adewumi telah melakukan penelitian untuk mensimulasikan aliran gas dalam saluran pipa menggunakan skema numerik Total Variations Diminishing (TVD) yang memiliki tingkat kesulitan yang tinggi dan bila dilihat dari segi numerik mempunyai waktu komputasi yang lama. Dan pada tahun 1996 mereka meneliti hal yang sama dengan menggunakan skema numerik Godunov. Model yang dikembangkan mengasumsikan bahwa aliran gas dalam kondisi isothermal. Faktor deviasi gas (Z) bernilai konstan dan dihitung dengan menggunakan persamaan Dranchuk dan Abou-Kaseem yang diambil secara rata-rata. Pada tahun 2000, Zhou dan Adewumi melaporkan lagi hasil penelitiannya pada kasus yang sama dengan menggunakan hybrid TVD/Godunov Scheme, hybrid TVD/Roe Scheme dan hybrid TVD/Lax-Wendrof. Dinyatakan bahwa hasil simulasi dengan menggunakan hybrid TVD/Godunov Scheme dan hybrid TVD/Roe Scheme memberikan hasil yang memuaskan. Sedangkan hasil simulasi dengan menggunakan hybrid TVD/Lax-Wendrof memberikan hasil yang cukup baik. Ketiga penelitian yang dilakukan oleh Zhou dan Adewumi tersebut mengasumsikan bahwa harga Z yang digunakan adalah faktor Z rata-rata yang dihitung dengan menggunakan metode Dranchuk dan Abou-Kassem. Dalam laporannya tidak disebutkan bahwa apakah harga Z diambil rata-rata dari harga inlet dan outlet ataukah harga Z dihitung dari kondisi awal yang masih steady state. Asumsi lainnya adalah bahwa kecepatan suara (c) adalah 2
Tinjauan Pustaka Studi dari Penelitian Sebelumnya
aliran gas dinamik satu dimensi yang melalui suatu saluran pipa.
konstan sepanjang pipa pada setiap waktu. Parameter yang diamati oleh Zhou dan Adewumi adalah perubahan tekanan dan laju alir pada pipa tunggal. Mereka membuktikan bahwa model yang dikembangkan sesuai dengan kondisi di lapangan yaitu dengan melakukan validasi dengan data lapangan yang diambil dari paper yang telah dipublikasikan oleh Demsey, J.R., dkk (1972). 2.2
Gas Alam Gas alam (natural gas) sering juga disebut sebagai gas bumi atau gas rawa, adalah bahan bakar yang berasal dari fosil kemudian menjadi gas, dimana susunan kimiawinya yang utama terdiri dari metana (CH4). Gas alam dapat ditemukan di ladang minyak, ladang gas bumi dan juga tambang batu bara. Komponen utama dalam gas alam adalah metana (CH4), yang merupakan molekul hidrokarbon rantai terpendek dan teringan. Gas alam juga mengandung molekul-molekul hidrokarbon yang lebih berat seperti etana (C2H6), propana (C3H8) dan butana (C4H10), selain juga gas-gas yang mengandung sulfur (belerang), susunan tersebut dapat dilihat di Tabel 2.1. Gas alam juga merupakan sumber utama untuk sumber gas helium.
Gambar 2.1: Volume Kendali dalam Suatu Saluran Pipa
2.4 2.4.1
Gas Nyata Persamaan Keadaan Gas Nyata Pada kenyataanya semua gas yang ada di alam tidak ada yang bersifat ideal. Namun, perilaku dari kebanyakan gas nyata tidak berbeda jauh dari perilaku gas ideal. Oleh karena itu, digunakan sebagai faktor pengkoreksi atau faktor deviasi persamaan gas ideal, sehingga persamaan keadaan gas nyata dapat dinyatakan sebagai: (2.1) pV = ZnRuTemp dengan: p : tekanan V : volume Z : faktor deviasi n : jumlah mol gas Ru : konstanta gas universal Temp : suhu Dengan mensubstitusi persamaan m = ρV dan m = nM pada Persamaan (2.1), adalah massa dengan adalah massa, molekul relatif gas, dan adalah massa jenis gas, maka persamaan keadaan gas tersebut dapat ditulis: ρ ZRuTemp (2.2) p= = ρ ZRg Temp
Tabel 2.1: Komponen Gas Alam (McCain, 1990) Komponen Gas Alam
%
Metana (CH4)
70 - 98
Etana (C2H6)
1 - 10
Propana (C3H8) and Butana (C4H10)
<5
2.3
Jaringan Pipa Jaringan pipa adalah sistem dengan rangkaian pipa yang panjang. Panjang dari jaringan pipa ini bisa mencapai ribuan kilometer. Sedangkan distribusi gas bisa memiliki arti menghubungkan pusat produksi dengan tempat-tempat penyimpanan gas. Salah satu persoalan pokok yang paling penting dalam sistem distribusi gas adalah perhitungan tekanan yang harus diberikan pada gas dari stasiun kompresi untuk mencapai titik-titik pemakaian. Dengan menerapkan hukum konservasi massa, momentum, energi dan persamaan keadaan gas pada perhitungan volume kendali seperti yang ditunjukkan pada Gambar 2.1, akan didapatkan suatu model yang menjelaskan
M
dengan Rg = Ru / M 2.4.2
Faktor Deviasi Faktor deviasi adalah perbandingan volume aktual n mol gas pada tekanan dan temperatur tertentu dengan volume n mol gas pada tekanan dan temperatur tertentu jika berperilaku sebagai gas ideal, yaitu: V (2.3) Z = actual Videal
3
dengan kuantitas yang ditransportasikan adalah massa, momentum dan energi. Oleh karena itu, model persamaan aliran gas dalam saluran pipa diturunkan dari hukum konservasi massa, momentum, energi dan persamaan keadaan gas nyata. 2.5.1 Persamaan Konservasi Massa Persamaan konservasi massa diturunkan dengan menggunakan prinsip hukum kekekalan massa. Dari penurunan tersebut, maka didapatkan persamaan konservasi massa sebagai berikut: ∂ρ ∂ (2.6) + ( ρv) = 0 ∂t ∂x
Gambar 2.2: Bentuk Tipikal Faktor Sebagai Fungsi Tekanan (McCain, 1990)
2.4.3
Faktor Gesekan Faktor gesekan merupakan penyebab terjadinya kerja yang hilang selama proses aliran. Faktor gesekan yang akan dipakai dalam penelitian ini dinyatakan oleh faktor Fanning , yang didefinisikan sebagai berikut: τ (2.4) f =
2.5.2
Persamaan Konservasi Momentum Persamaan konservasi momentum diturunkan dengan menggunakan prinsip hukum kekekalan momentum (Hukum Kedua Newton). Dari penurunan tersebut, maka persamaan konservasi momentum untuk volume kendali dengan luas penampang dan sudut inklinasi dapat ditulis sebagai berikut:
1 2 ρv 2
dengan: τ : tegangan geser dinding v : kecepatan aliran gas
∂ ∂ ( ρ vS ) + pS + ρ v 2 S + τ π D + ρ Sg sin θ = 0 ∂t ∂x
(
(2.7) Node Network atau Mesh Node adalah pembagian media menjadi daerah-daerah kecil, kumpulan dari node tersebut disebut node network atau mesh. Node network untuk metode beda hingga Crank-Nicolson dapat ditunjukkan pada Gambar 2.3.
2.4.4 Persamaan Kecepatan Gelombang Suara Kecepatan suara dalam gas disini didefinisikan sebagai kecepatan suatu gangguan kecil di dalam tabung aliran fluida.
2.6
Volume Kendali
S
)
S
Gambar 2.3: Skematik Penurunan Persamaan Kecepatan Suara (Streeter dan Wylie, 1990)
Gambar 2.4: Node Network untuk Metode CrankNicolson
Dengan menurunkan gerak gelombang tekanan akibat gangguan kecil dalam suatu volume kendali, didapatkan persamaan kecepatan suara sebagai berikut: c=
ZRuTemp = ZRgTemp M
2.7
Metode Crank-Nicolson Metode Crank-Nicolson diperoleh dari rata-rata metode Eksplisit dan metode Implisit. u j , n +1 − u j ,n ⎛ ∂u ⎞ (2.8) ⎜ ⎟ ≈ Δt ⎝ ∂t ⎠ j ,n
(2.5)
1 ⎛ u j +1,n − u j −1,n u j +1,n +1 − u j −1,n +1 ⎞ ⎛ ∂u ⎞ + ⎟ ⎜ ⎟ ≈ ⎜ 2Δx 2Δx ⎝ ∂x ⎠ j ,n 2 ⎝ ⎠ (2.9)
2.5 Persamaan Dasar Aliran Gas Unsteady Aliran gas dalam saluran pipa dapat dipandang sebagai masalah transportasi, 4
⎛ ∂ 2u ⎞ 1 ⎛ u j −1,n − 2u j ,n + u j +1,n u j −1,n +1 − 2u j ,n +1 + u j +1,n +1 ⎞ + ⎜ 2⎟ ≈ ⎜ ⎟ Δx 2 Δx 2 ⎝ ∂x ⎠ j ,n 2 ⎝ ⎠
∂ ∂ ( ρ vS ) + pS + ρ v 2 S + τ π D + ρ Sg sin θ = 0 ∂t ∂x
(2.10) dimana menyatakan titik grid spatial dan menyatakan step waktu.
ZR Temp p = ρ. u = ρ ZRg Temp M
(
IV. Pemodelan dan Penyelesaian Numerik 4.1 Pemodelan Matematika Aliran Gas Unsteady dalam Saluran Pipa Pemodelan aliran gas pada Tugas Akhir ini didasarkan pada model aliran gas unsteady yang dikembangkan oleh Zhou dan Adewumi (1995) yang mempunyai beberapa asumsi sebagai berikut: 1. Pipa transmisi lurus dan horizontal 2. Luas penampang pipa konstan 3. Aliran gas satu dimensi 4. Aliran berada pada kondisi isothermal (suhu di sepanjang saluran pipa konstan) 5. Fluida aliran adalah gas satu fasa (gas kering) 6. Ekspansi dinding pipa diabaikan 7. Berlaku gesekan untuk kondisi unsteady 8. Kecepatan suara ( ) dianggap konstan 9. Tidak ada kerja yang dilakukan gas selama aliran terjadi 10. Massa jenis gas dianggap konstan Namun pada Tugas Akhir ini, pipanya tidak horizontal, tetapi mempunyai sudut inklinasi. Kemudian, karena aliran diasumsikan isothermal yaitu tidak ada gangguan dari luar sehingga suhu di dalam saluran pipa konstan, maka persamaan konservasi energi tidak diperlukan dalam membangun model matematika disini. Dengan demikian, berdasarkan asumsi-asumsi yang telah dibuat, model matematika yang dikembangkan untuk menjelaskan perilaku aliran gas unsteady dalam saluran pipa terdiri dari persamaan konservasi massa, persamaan konservasi momentum, dan persamaan keadaan gas nyata, yang dapat dituliskan kembali sebagai berikut: ∂ρ ∂ (4.1) + ( ρv) = 0
disini, kecepatan ( v ) akan diganti dengan laju alir massa yang didefinisikan sebagai berikut: (4.5) q = ρ vS = ρ Q = ρ nQn Untuk penyederhanaan pertama, dengan menggunakan Persamaan (4.4) dan (4.5), persamaan konservasi massa dan momentum dapat dinyatakan sebagai fungsi dari laju alir dalam kondisi normal Qn ( x, t ) dan tekanan p ( x, t ) . Maka Persamaan (4.1) dapat ditulis kembali menjadi: c 2 ρ n ∂Qn ∂p (4.6) =− ∂t S ∂x Sedangkan untuk Persamaan (4.2), dengan menggunakan faktor Fanning ( ) pada Persamaan (2.4) dan persamaan aliran massa pada Persamaan (4.5), diperoleh: 2 2 f ρv fq (4.7) = τ = 2 2ρ S 2 Dari Persamaan (4.2) dan Persamaan (4.7), maka didapatkan: ∂q ∂ ⎛ q2 ⎞ fq 2 + ⎜ Sp + π D + ρ Sg sin θ = 0 ⎟+ ∂t ∂x ⎝ S ρ ⎠ 2ρ S 2 (4.8) Dengan menggunakan definisi luas 2 penampang pipa, S = 1 π D 2 dan karena q
4
Sρ
konstan terhadap jarak maka Persamaan (4.8) menjadi: ∂Qn S ∂p 2 fc 2 ρ n Qn2 Sg sin θ p (4.9) =− − − ρ n ∂x DS p ρnc 2 ∂t Sehingga dari Persamaan (4.6) dan (4.9) didapatkan penyederhanaan pertama untuk model adalah:
∂x
5
(4.2) (4.3)
dengan Rg = Ru / M Dengan menggunakan kecepatan suara ZRuTemp . Persamaan keadaan gas nyata c= M pada Persamaan (4.3) menjadi: (4.4) p = c2 ρ Permintaan gas oleh konsumen biasanya dinyatakan dengan laju alir massa ( q ) . Maka
III. Prosedur Kerja 1. Studi literatur. 2. Pemodelan aliran gas dalam saluran pipa. 3. Penyelesaian numerik. 4. Simulasi model.
∂t
)
⎧ ⎪ ⎪ ⎨ ⎪ ∂Qn = − S ⎪⎩ ∂t ρn
c 2 ρ n ∂Qn ∂p =− S ∂x ∂t (4.10) ∂p 2 fc 2 ρ n Qn2 Sg sin θ p − − ρn c 2 DS p ∂x
Akhirnya, dengan menggunakan Q ρ = Qn ρ n , maka kembali hubungan penyederhanaan kedua dari model dapat dinyatakan sebagai berikut: ∂ 2 p 2 2 g sin θ ∂p 2 16Qf ∂p 2 (4.14) + = ∂x 2 c2 ∂x c 2π D3 ∂t Dengan mengasumsikan bahwa dan β adalah suatu nilai konstan, maka Persamaan (4.14) menjadi suatu model linier dengan respek terhadap p 2 . Dan asumsi yang
Penyederhanaan selanjutnya dari model dapat diturunkan dari hipotesis bahwa kondisi batas tidak berubah secara cepat dan kapasitas dari saluran gas relatif besar sehingga bisa dianggap konstan terhadap jarak dan juga waktu. Pada kasus ini, kondisi pertama dari Persamaan (4.10) dapat juga dieliminasi dari model. Kemudian, dengan mengekspresikan lagi model tersebut sebagai fungsi dari laju alir dalam kondisi normal Qn ( x, t ) dan tekanan p ( x, t ) . Didapatkan
kedua yaitu Q ( x, t ) merupakan suatu rata-rata laju alir sepanjang saluran gas dalam tiap interval waktu ∆ . Dari asumsi tersebut, maka telah didapatkan model aliran gas dalam saluran pipa yang berupa PDP parabolik linier orde dua dengan respek terhadap dalam tiap ∆ adalah sebagai berikut: ⎧ ∂2 p2 ∂p 2 ∂p 2 +β =α ⎪⎪ 2 ∂x ∂x ∂t (4.15) ⎨ ⎪dengan α = 16Qf dan β = 2 g sin θ ⎪⎩ c 2π D3 c2
persamaan sebagai berikut: 1 ∂q 2 fc 2 ρn2 Qn2 g sin θ ∂p p =− − − S ∂t DS 2 p c2 ∂x Dari asumsi yang diberikan bahwa q konstan terhadap jarak dan waktu, maka: 2 fc 2 ρn2 Qn2 g sin θ ∂p ∂x
=−
Dengan
DS 2
p
−
c2
p
menggunakan
persamaan
4.2 Skema dan Penyelesaian Numerik 4.2.1 Syarat Awal dan Syarat Batas 1. Syarat Dirichlet (4.16) p ( 0, 0 ) = p ( 0 )
∂p ∂p = 2 p , didapatkan: ∂x ∂x 2
2 2 2 2 ∂p ⎞ ⎛ ∂p ⎞ ⎛ ∂p ⎞ ⎛ 2 fc ρ n Qn g sin θ ⎞ ⎛ = − − p ⎟⎜ 2 p ⎟ ⎟ ⎜ ⎜ ⎟⎜ ∂x ⎠ DS 2 p c2 ⎝ ∂x ⎠ ⎝ ∂x ⎠ ⎝ ⎝ ⎠
p ( 0, t ) = p ( 0 ) ; t > 0
4 fc 2 ρ n2Qn2 2 g sin θ 2 ∂p 2 (4.11) =− − p DS 2 c2 ∂x Dan ketika Persamaan (4.11) diturunkan parsial terhadap x , didapatkan: ∂ 2 p 2 ∂ ⎛ 4 fc 2 ρn2Qn2 2 g sin θ 2 ⎞ = ⎜− − p ⎟ ∂x 2 ∂x ⎝ DS 2 c2 ⎠
2.
(4.13) hubungan
∂p |x =0 = 0 ∂t ∂Q |0≤ x ≤ L = 0 ∂t
p ( x) =
∂p ∂p = 2 p , Persamaan (4.13) menjadi: ∂t ∂t ∂ 2 p 2 16Qn f ρ n ∂p 2 2 g sin θ ∂p 2 = 2 − c ρπ D3 ∂t c2 ∂x 2 ∂x 2
(4.18)
(4.19) (4.20)
p (0) − μ x 2
(4.22)
dengan μ = f ⎛ 2 ρ cQ ⎞ dan σ = 2 g sin θ ⎜ ⎟ 2 2
D⎝
S
⎠
c
Sedangkan nilai Z dihitung dari suatu syarat awal dengan menggunakan algoritma yang ditunjukkan dalam Gambar 4.1. Pada penelitian ini, berdasarkan persamaan empiris 6
Q ( x, t ) = konstan; 0 ≤ x ≤ L Syarat Neumann
Untuk perhitungan tekanan pada tiap titik pada saat t = 0 sepanjang pipa dapat dinyatakan oleh Persamaan (4.21) untuk kasus dimana saluran pipa mempunyai sudut inklinasi dan oleh Persamaan (4.22) untuk saluran pipa horizontal. μ 2 ⎛ ⎞ (4.21) p ( x ) = ⎜ p ( 0 ) − ( e μ x − 1) ⎟ e − μ x σ ⎝ ⎠
8Qn fc 2 ρ n2 ∂Qn 2 g sin θ ∂p 2 ∂2 p2 (4.12) = − − ∂x 2 ∂x ∂x DS 2 c2 Dari Persamaan (4.6) dimasukkan ke Persamaan (4.12), maka diperoleh: ∂ 2 p 2 8Qn f ρ n ∂p 2 g sin θ ∂p 2 = − DS ∂t c2 ∂x 2 ∂x Dengan menggunakan
(4.17)
yang telah diusulkan oleh Segeler, Ringler, dan Kafka, nilai Z dihitung dengan menggunakan persamaan: p (4.23) Z = 1 − rata − rata 390 dengan prata − rata = 1 ⎡ p ( 0 ) + p ( L ) ⎤ 2⎣
Menghitung
hasil pendiskritan pada Persamaan (4.26) menjadi: ⎛ ( −2 a − 2 ) ⎜ ⎜ (a − b) ⎜ 0 ⎜ 0 ⎜ ⎜ # ⎜ ⎜ 0 ⎝
⎦
,
atau (4.22)
Menghitung Z dengan Z0 = Z
menggunakan Persamaan (4.23)
Benar Selesai
% 0
% 0
⎛ ( 2a − 2 ) ( − a − b ) 0 ⎜ ⎜ ( − a + b ) ( 2a − 2 ) ( − a − b ) ⎜ 0 ( − a + b ) ( 2a − 2 ) =⎜ 0 0 ( −a + b ) ⎜ ⎜ # % % ⎜ ⎜ 0 0 0 ⎝ ⎛ ( −a + b ) ( u0,n + u0,n +1 ) ⎞ ⎜ ⎟ 0 ⎜ ⎟ ⎜ ⎟ 0 ⎟ +⎜ ⎜ ⎟ 0 ⎜ ⎟ # ⎜ ⎟ ⎜ ( −a − b ) ( u J +1,n + u J +1,n +1 ) ⎟ ⎝ ⎠
dengan
menggunakan Persamaan (4.21)
Z0 = Z
0 0 0 (a + b) 0 0 ( −2a − 2 ) ( a + b ) 0 ( a − b ) ( −2 a − 2 ) ( a + b ) 0 ( a − b ) ( −2 a − 2 ) ( a + b ) % 0 0 0
% 0 0 0
0 ( −a − b ) ( 2a − 2 ) ( − a − b )
" %
0 0
% %
0 0
% "
(a − b)
%
" %
0 0
% %
0 0
%
%
%
%
0
0
"
( −a + b )
⎞ ⎛ u1,n +1 ⎞ ⎟⎜ ⎟ ⎟⎜ u2,n +1 ⎟ ⎟⎜ u3,n +1 ⎟ 0 ⎟⎜ ⎟ 0 ⎟⎜ u4,n +1 ⎟ % ⎟⎜ # ⎟ ⎟⎜ ⎟ ( −2a − 2 ) ⎟⎠ ⎜⎝ u J ,n+1 ⎟⎠ 0 0
⎞⎛ u1,n ⎞ ⎟⎜ ⎟ ⎟⎜ u2,n ⎟ 0 ⎟⎜ u3,n ⎟ ⎟⎜ ⎟ 0 ⎟⎜ u4,n ⎟ % ⎟⎜ # ⎟ ⎟⎜ ⎟ ( 2a − 2 ) ⎟⎠ ⎜⎝ uJ ,n ⎟⎠ 0 0
(4.27)
Salah |Z – Z0| < Tol
Gambar 4.1: Algoritma untuk Menghitung Z dan
V. 5.1
4.2.2
Untuk menyelesaikan Persamaan (4.27), disusun algoritma penyelesaian sebagai berikut: 1. Mendefinisikan parameter-parameter yang dibutuhkan. 2. Menghitung syarat batas dan faktor deviasi Z menggunakan algoritma seperti yang ditunjukkan pada Gambar 4.1. 3. Dari Z yang telah diperoleh pada langkah 2, dihitung nilai syarat awal menggunakan Persamaan (4.21) untuk kasus dimana saluran pipa mempunyai sudut inklinasi atau Persamaan (4.22) untuk saluran pipa horizontal. 4. Memasukkan syarat awal dan syarat batas ke dalam skema numerik. 5. Skema numerik yang berupa matrik tridiagonal diselesaiakan menggunakan algoritma eliminasi Gauss. 6. Selanjutnya dihitung nilai error berdasarkan hasil perhitungan numerik dan perhitungan eksak.
p ( L, t )
Skema Numerik Dengan memisalkan p 2 = u , maka model aliran gas pada Persamaan (4.15) dapat ditulis kembali sebagai berikut: ∂ 2u ∂u ∂u (4.24) +β =α 2 ∂x ∂x ∂t Pada persamaan tersebut, tiap kondisi dapat didekati dengan skema beda hingga. Dengan menerapkan pendekatan metode Crank-Nicolson untuk model aliran gas seperti yang ditunjukkan dalam Persamaan (2.8), (2.9) , dan (2.10), maka akan didapatkan skema numerik sebagai berikut: ⎛ −4Δt + 4αΔx 2 ⎞ ⎛ 2Δt − β ΔxΔt ⎞ ⎛ 2Δt + β ΔxΔt ⎞ ⎟ u j ,n + ⎜ ⎜ ⎟ u j −1,n + ⎜ ⎟ u j +1,n + 2 2 2Δx 2 Δt ⎝ 2 Δx Δ t ⎠ ⎝ 2 Δx Δt ⎠ ⎝ ⎠ ⎛ −4Δt − 4αΔx 2 ⎞ ⎛ 2Δt − β ΔxΔt ⎞ ⎛ 2Δt + β ΔxΔt ⎞ ⎟ u j ,n +1 + ⎜ ⎜ ⎟ u j −1,n +1 + ⎜ ⎟ u j +1, n +1 = 0 2 2 2 Δx 2 Δ t ⎝ 2 Δx Δ t ⎠ ⎝ 2 Δx Δ t ⎠ ⎝ ⎠
(4.25) ∆
Persamaan (4.25) dikalikan dengan , maka selanjutnya didapatkan skema numerik untuk j = 1, 2, 3, … , J dan n = 0,1, 2,…, N adalah: ( a − b ) u j −1,n+1 + ( −2a − 2 ) u j ,n+1 + ( a + b ) u j +1,n+1 =
5.2
Program Algoritma pada Subbab 5.1 tersebut diterapkan ke dalam program dengan menggunakan bantuan program 7.6. Listing program selanjutnya disajikan dalam lampiran.
( −a + b ) u j −1,n + ( 2a − 2 ) u j ,n + ( −a − b ) u j +1,n (4.26) Δt 2 Δt 2 dengan = a= 2αΔx 2 Δt αΔx 2 β Δ x Δt 2 β Δt b= = 2 2αΔx Δt 2αΔx
dan
5.3
Simulasi Pada penelitian ini disimulasikan beberapa skenario yang berbeda untuk memperoleh gambaran yang lebih lengkap
Dengan memasukkan syarat awal dan syarat batas pada skema numerik, selanjutnya 7
Simulasi dan Pembahasan Algoritma Program
tentang perilaku aliran gas transient dalam saluran pipa. Contoh kasus yang disajikan didasarkan pada beberapa perkiraan kondisi operasi yang ada di lapangan. Berikut adalah tampilan dari interface simulasi:
Gambar 5.3: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 50 m3/s dan H = 7800 m
Dari simulasi kedua didapatkan 19,1293 bar 0,00185209.
Gambar 5.1: Interface Simulasi
T = 3600 s
L = 10 m v = 28 m / s Rg = 392 m2 / s 2 K
ρ = 0.73 kg / m3
5
nilai
error
sebesar
c. Untuk Qo = 50 m3/s dan H = 500 m, didapatkan hasil sebagai berikut:
Untuk kebutuhan simulasi, digunakan parameter-parameter sebagai berikut: D = 0.6 m Temp = 278 o K p ( 0, t ) = 50 bar
dan
p ( L, t ) =
f = 0.003
Pada Kasus 1 dan Kasus 2, digunakan M=N=50. Sedangkan pada Kasus 3, nilai M dan N tersebut diubah. 1.
Kasus 1 Pada kasus ini disimulasikan beberapa skenario dengan memasukkan Qo = 50 m3/s dan nilai H yang berbeda-beda. a. Untuk Qo = 50 m3/s dan H = 3000 m, didapatkan hasil simulai sebagai berikut:
Gambar 5.4: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 50 m3/s dan H = 500 m
Dari simulasi ketiga didapatkan
44,1291 bar dan nilai error sebesar 0,00013069. Dari ketiga hasil yang diperoleh di atas, terlihat bahwa semakin kecil nilai H, maka nilai syarat batas p(L) semakin besar dan grafik distribusi tekanan akan semakin mendekati linear. Dan terlihat juga bahwa nilai error semakin kecil. 2.
Kasus 2 Sedangkan pada kasus kedua ini, nilai Qo dinaikkan dan diberikan dua nilai H yang berbeda. a. Untuk Qo = 100 m3/s dan H = 4750 m, didapatkan hasil sebagai berikut:
Gambar 5.2: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 50 m3/s dan H = 3000 m
Dari simulasi pertama didapatkan p ( L, t ) = 33,485 bar dan nilai error sebesar 0,00075822. b. Untuk Qo = 50 m3/s dan H = 7800 m, didapatkan hasil simulai sebagai berikut: 8
p ( L, t ) =
Gambar 5.5: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 100 m3/s dan H = 4750 m
Gambar 5.7: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 50 m3/s, H = 3000 m, J=N=30
Dari simulasi keempat didapatkan p ( L, t ) = 10,734 bar dan nilai error sebesar 0,00060289.
Dari simulasi keenam didapatkan nilai error sebesar 0,00076542.
b. Untuk Qo = 100 m3/s dan H = 2000 m, didapatkan hasil sebagai berikut:
b. Untuk Qo = 50 m3/s ,H = 3000 m, M=N=70, didapatkan hasil simulai sebagai berikut:
Gambar 5.6: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 100 m3/s dan H = 2000 m
Gambar 5.8: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 50 m3/s, H = 3000 m, J=N=70
Dan dari simulasi kelima didapatkan p ( L, t )
Dari simulasi ketujuh didapatkan nilai error sebesar 0,00075046. Dari kedua contoh pada Kasus 3 dan dengan melihat kembali contoh pada Kasus 1 untuk Qo = 100 m3/s, H = 3000 m tetapi dengan M=N=50, terlihat bahwa semakin banyak jumlah pendiskritan, maka nilai error semakin kecil. Jadi hasil perhitungan numerik semakin mendekati hasil perhitungan eksak.
= 24,3797 bar dan nilai error sebesar 0,00026281. Dari kedua hasil pada kasus kedua di atas menunjukkan bahwa ketika nilai Qo dinaikkan, maka hal tersebut akan mengakibatkan nilai syarat batas p ( L, t ) semakin kecil. Jadi, perubahan nilai H dan Qo sangat mempengaruhi besar kecilnya nilai syarat batas p ( L, t ) .
4.
Kasus 4 Untuk mengetahui bagaimana perilaku tekanan gas ketika saluran gasnya menurun atau sinus sudut inklinasinya bernilai negatif, pada kasus ini disimulasikan dua contoh skenario. a. Untuk Qo = 50 m3/s , H = -1000 m, didapatkan hasil simulai sebagai berikut:
3.
Kasus 3 Untuk melihat hubungan banyaknya jumlah pendiskritan dengan besarnya nilai error, akan disimulasikan dua contoh sebagai berikut. a. Untuk Qo = 50 m3/s ,H = 3000 m, J=N=30, didapatkan hasil simulai sebagai berikut:
9
VI. 6.1
Simpulan Kesimpulan
Berdasarkan pembahasan dari bab sebelumnya maka dapat disimpulkan bahwa: 1. Model matematika yang menggambarkan perilaku aliran gas dalam suatu saluran pipa dapat dinyatakan sebagai berikut:
⎧ ∂2 p2 ∂p 2 ∂p 2 β α + = ⎪⎪ ∂x 2 ∂x ∂t ⎨ ⎪dengan α = 16Qf dan β = 2 g sin θ ⎪⎩ c 2π D3 c2
Gambar 5.9: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 50 m3/s, H = -1000 m
Dari simulasi kedelapan didapatkan p ( L, t ) = 52,2059 bar 0,00026687.
dan
nilai
error
2.
sebesar
b. Untuk Qo = 50 m3/s , H = -3000 m, didapatkan hasil simulai sebagai berikut:
Dari simulasi dan penyelesaian numerik yang diperoleh dengan bantuan program Matlab 7.6, terlihat bahwa nilai tekanan di outlet dipengaruhi oleh nilai parameterparameter, diantaranya laju alir dan ketinggian di ujung pipa. Sedangkan besarnya error, selain dipengaruhi oleh nilai parameter, juga dipengaruhi oleh banyaknya pendiskritan. Semakin banyak pendiskritan, maka nilai error semakin kecil.
6.2
Saran Untuk pengembangan penelitian selanjutnya, disarankan: 1. Pada Tugas Akhir ini menggunakan asumsi bahwa laju alir gas konstan sepanjang saluran pipa, selanjutnya dapat dikembangkan penelitian untuk laju alir gas yang tidak konstan agar lebih mendekati kondisi di lapangan. 2. Tugas Akhir ini masih bersifat analitis pada tahap pemodelan dan numerik untuk penyelesaiannya, belum ada data laboratorium yang dipakai sebagai pembanding. Diharapkan kedepannya bisa dilakukan uji laboratorium sehingga model ini dapat diterapkan di lapangan.
Gambar 5.10: Grafik Distribusi Tekanan Gas dalam Saluran Pipa untuk Qo = 50 m3/s,H -3000 m
Dari simulasi kesembilan didapatkan p ( L, t ) = 65,4484 bar dan nilai error sebesar 0,00082373. Dari kedua contoh kasus pada Kasus 4, terlihat bahwa ketika saluran pipa menurun, semakin lama tekanan gas akan semakin meningkat. Hal ini terbukti dari hasil simulasi, jika diberikan tekanan awal sebesar 50 bar, maka tekanan akhirnya menjadi 52,2059 bar untuk H = -1000 m. Terlihat juga bahwa nilai error-nya semakin besar jika sinus sudutnya semakin negatif. Melihat nilai error pada masingmasing simulasi, terlihat bahwa nilai error cukup kecil. Maka dapat disimpulkan bahwa penyelesaian model aliran gas menggunakan metode Crank-Nicolson memberikan hasil yang baik. Besarnya nilai error selain dipengaruhi oleh nilai parameter-parameter, juga dipengaruhi oleh banyaknya jumlah pendiskritan. Untuk nilai dari matriks tekanan dari simulasi dapat dilihat dalam lampiran.
Daftar Pustaka Gonzalez, Alberto H., dkk. 2008. Modeling and Simulation of a Gas Distribution Pipeline Network. Journal of Applied Mathematical Modelling, 1584–1600. Lurie, Michael V. 2008. Modeling of Oil Product and Gas Pipeline Transportation. Weinheim: WileyVCH Verlag GmbH & Co. KGaA. McCain, William D. Jr. 1990. The Properties of Petroleum Fluids. Oklahoma: PennWell Publishing Company. 10
Segeler, C.G., Ringler, M.D., dan Kafka, E.M. 1969. Gas Engineers’ Handbook. New York: AGA. Smith, G.D. 2005. Numerical Solution of Partial Differential Equations. Oxford: Clarendon Press. Streeter, V.L. dan Wylie, E.B. 1990. Diterjemahkan oleh Prijono, A. Mekanika Fluida Jilid I Edisi 8. Jakarta: Erlangga. Sulistyarso, Harry B. 2007. Aplikasi Suatu Model Aliran Gas Transient pada Kasus Line Packing untuk Lapangan Gas. Disertasi, Departemen Teknik Perminyakan, Institut Teknologi Bandung. Sulistyarso, Harry B., dkk. 2004. Solusi Model Aliran Gas Dalam Pipa pada Kondisi Line Packing Menggunakan Skema Richtmyer. Proceedings ITB Sains & Teknik, Vol. 36A, No. 2, 159177. Yang, Won Y., dkk. 2005. Applied Numerical Methods Using MATLAB. John Wiley & Sons, Inc. Zhou, Junyang dan Adewumi, M. A. 2000. Simulation of Transient in Natural Gas Pipelines Using Hybrid TVD Schemes. International Journal of Numerical Method for Fluids, 32: 407 – 437.
11
