MODEL MATEMATIK DAN NUMERIK DARI ALIRAN FLUIDA DALAM JARINGAN PIPA BERBASIS MER DENGAN FORMULASI AKAR TERKECIL Ade Jamal,Agus Sainjati,Aris Suwatjono,LebongAndalaluna'
ABSTRAK MODEL MATEMATIK DAN NUMERIK DARI ALmAN FLUffiA DALAM JARINGAN PIPA BERBASIS MER DENGAN FORMULASI AKAR TERKECIL. PersamaanmatematisdaTi analisis aliran fluida dalamjaringan pipa dibahas untuk diterapkandalam pemodelannumerik yang sesuai.Dalam hal ini dipilih model numerik berbasisMetode ElemenHingga (MEH). Dalam rangka menjagamodel numerik MEH ini dapatdiaplikasiseluas-luasnya untuk segalakasusaliran fluida, maka formulasiakar terkecil dipilih sebagaipenggantiformulasi Galerkin atau formulasi Residu. Persarnan aliran fluida yang dikaji dibatasipada alirantunak dan tak-mampumampat.Selainitu, artikel ini dibatasi hanya sampaipengembanganmodel numeriknya,tidak memasukkankajian lanjutan seperti prosedur penyelesaain modelyangdibuat.
ABSTRACT MATHEMATICAL AND NUMERICAL MODEL OF FLUm FLOW IN A PIPE NETWORK BASED ON FEM WITH LEAST ROOT FORMULATION. Mathematicalformulations for analysis of fluid flow in a pipe network has been studied to be implementedin an appropriate numericalmodel. In this case,a numerical model based on Finite ElementMethod (FEM) has been chosen.In orderto keepthe range of applicability as wide as possible,a finite elementmodelhas been derivedusingthe LeastSquare formulationin place of the frequentlyused WeightedResidual(Galerkin) formulation.The fluid flow studiedhere is assumedto be incompressibleand steady.Furthermore,this article was limited itself to the developmentof the numericalmodel. Hence, further studies,such as solutionprocedureof the numericalmodel,areexcluded.
PENDAHULUAN Jaringanpipa untuk aliranfluida, baik cairanmaupungas,dirancangsedemikian rupa hingga tidak terjadi hambatanaliran yang menyebabkantidak ekonomisnya energi tekananyang dibutuhkan. Sistemrancangandistribusi didasarkandua faktor utama yaitu keperluanjumlah aliran (debit) dari fluida yang ditransportasikandan besar tekanan (energi) yang dibutuhkan. Kegiatan rekayasa ini dikenal sebagai .Pusat PengkajiandanPenerapan Teknologi lnformasi danElektronika,BPPTeknologi
217
RisalahLokakaryaKornputasidalamgains daDTeknologiNuklir XIV, Juli 2003 (217-228)
analisis dan perancanganjaringan pipa. Analisis jaringan pipa juga dibutuhkan untuk operasional clan pengontrolan, akuisisi suplai, optimisasi kinerja jaringan terhadap biaya, clan lain-lain. Sistem tersebut dapat dimodelkan sebagai sambungan seri mapun paralel daTi elemen-elemen pipa yang berhubungan satu sarna lain. Analisis clan rancanganjaringan pipa menimbulkan masalah yang relatif kompleks, terutama sekali jika jaringan terdiri daTi banyak pipa. Analisis tersebut umumnya menggunakan program komputer. Tahapan yang kritis adalah penentuan individual penurunan tekanan clan debit aliran dalam elemen pipa. Model numerik yang dibuat di sini membatasi masalah aliran fluida yang tunak (steadyflow), tak mampu mapat (incompressible flow), tapi tetap mengikut sertakan peranan gesekan(viscos flow) baik untuk daerah laminar maupun turbulen.
Model Matematis Aliran Dalam Jaringan Pipa Persamaan fundamentaluntuk aliran fluida secaraumumdiatur olehtiga hukum kekekalan,yaitu persamaankekekalanmasa, kekekalanmomentumdan kekekalan energi[I], yaitu:. persamaan kekekalanmasa:
-8p +div(p'V)=O,
(1)
at
(2) danpersarnaan kekekalanenergi
p
Dh
= Dp +div(k.VT)+'rij Dt Dt
-avi ax.J
(3)
Tiga persamaanini mengaturtiga variabeldasar:kecepatanE, tekananp, clan temperatur T. Variabel sekunder yang diikutsertakan adalah entalpi h, densitas fluida p. Untuk aliran dalam saluranpipa, tiga persamaandi atas dapatdisederhanakan menjadi problem satu dimensi. Selain itu permasalahandibatasi lebih jauh untuk kondisi tunak (steady state) clan tak-mampu mampat (incompressible),sehingga persamaan umumnyaadalah[2]:
218
Model Matematik dan Numerik daTi Aliran Fluida dalamjaringan Pipa Berbasis MEH
(Ade Jamal, et.al.»
persamaan kekekalanmasaataupersamaan kontinuitas
~=o, ox
(5) danpersamaan kekekalanenergi 1 OQh+ oW = dE + -dp + V .dV + g .dZ p
(6)
Variabel dasarnya adalah kecepatan aliran searahpipa u, tekanan dalam pipa p dan temperatur T yang implisit dinyatakan dalam fungsi energi dalam E, dan energi panas kedalam fluida Qh.Kerja yang dilakukan oleh fluida seperti gesekandinyatakan dalam W merupakan variabel sekunder. Viskositas .u merupakan parameter fisis daTi fluida. Kecepatan aliran rata-rata V didapat daTi persamaanberikut:
v =11U .dA,
(7)
Persarnaankesetimbangan gaya (5) yang mengikutsertakan viskositas (gesekan fluida), juga disebut persamaanNavier-Stokes, secara umum sangatjarang diketahui solusinya, kecuali untuk kasus-kasus yang sangatkhusus clan sederhana.Hal ini terjadi karena persamaanini menghasilkan persarnaandiferensial parsial clantak-linear. Salah satu pendekatan yang sering dipakai adalah pendekatan empiris di mana penurunan tekanan (pressure drop) atau kerugian tekanan (pressure loss) didekatkan dengan persamaanhasil percobaan seperti kerugian tekanan untuk aliran laminar daTi Hagen-Poiseuille:
(8)
Q=V.A
(9)
Risalah Lokakarya Komputasi dalam Sains dan Teknologi Nuklir XIV, Juli 2003
Untuk aliran turbulensi penuh digunakan persamaanempiris daTi DarcyWeisbach: (10)
Metode Solusi Aliran Fluida Dalam Jaringan Pipa Terdapatbeberapametodeyang tersediadalamliteraturyang seringdigunakan untuk menyelesaikan permasalahaliranfluida dalamjaringanpipa [3]. Padaumumnya metode solusi ini berdasarkananalogi hukum Kirchhoff yang biasa dipakai untuk jaringanelektrik,yaitu: Penjumlahansecaraaljabar debit aliran yang masuk clan keluar daTi semua cabang atau simpul (node) hams nolo Aturan ini diturunkan daTi persamaan kontinuitas(4) L Q masuk= LQkeluar
(11)
Penjumlahansecaraaljabarkerugiantinggi tekananpadasirkuit tertutup(closed loop) barnsnolo Aturan ini diturunkan dan persamaankekekalanenergi (6) setelah menerapkan hukum-hukumtermodinamikamenjadi: 2
-
-
p'g
2.g
~+Zl+~
/
-:-:-+Z2 + V; P2
2.g
=1
!i) loss
(12)
atau (13)
~ -~ =hL di mana v.2 I
2.g
hL =
~ ,p.g
loss
Metode yang menggunakan dua persamaan dasar ini antara lain metode "koreksi debit" atau lebih dikenal dengan nama metode Hardy-Cross [3,4] di mana debit untuk pipa awalnya diasumsikan untuk setiap pipa sehingga persamaan kontinuitas terpenuhi di setiap titik simpul, lalu secara berturut-turut diterapkan kesetimbangan energi pada setiap sirkuit untuk mendapatkan koreksi debit. Variasi dari metode ini yang lebih mudah untuk diprogramkan dengankomputer, yaitu dengan
220
Model Matematik dan Numerik daTi Aliran Fluida dalamjaringan Pipa Berbasis MEH
(Ade Jamal, et.al.»
cara persamaan untuk kesetimbangan diselesaikan secara serentak. Metode lain, di mana debit aliran Q dalam persamaan kontinuitas dieliminasi dengan bantuan persamaanenergi, sehingga persamaan yang harus diselesaikan hanya dalam variabel tekanan tapi tak-linear sepenuhnya.Metode matriks [3] diperkenalkan oleh Steffer, di mana persamaan kontinuitas dan energi diselesaikan secara bersamaan dalam bentuk matriks berukuran besar. Metode matriks ini sederhana dan persamaan yang dicari langsung didapat. Hanya saja, matriks ini sering terlalu besar dan tidak sempurna (ill-conditioned). Dari semua metode yang disebut di atas, yang paling populer adalah Metode Hardy-Cross, karena dapat dilakukan secara hitungan tangan untuk jaringan pipa yang kecil. Perlu dicatat di sini bahwa tidak ada satu metode pun yang menggunakan persamaanmomentum dalam proses penyelesaiannya. Metode Elemen Hingga (MEH) telah sangat populer sejak tiga dasawarsa terakhir. Berangkat daTi mekanika solid dan struktur di dunia teknik pesawat terbang, saat ini teknik MEH (FEM-Finite Element Method) telah diaplikasikan untuk bermacama-macamanalisa media kontinum, termasuk fluida. Yang menarik adalah kebanyakan literatur dan referensi yang ada hanya membahas MEH untuk aliran fluida pada domain dua dimensi dan tiga dimensi [5-7]. Untuk masalah fluida satu dimensi seperti aliran fluida dalam pipa tertutup, jarang yang mebahasnya [8], dan jika ada dalam literatur hanya terbatas pada aliran laminar [9].
Model ElemenHingga Formulasi Akar Terkecil (LeastSquare Formulation) Persamaandasar yang digunakan dalam model elemen hingga adalah persamaankontinuitas (4) dan persamaanmomentum yang diturunkan dengan pendekatanberdasarkankesetimbangangaya pada fluida dalam pipa seperti pada gambarberikut:
\-'-
Gambar
PersamaanMomentum BerdasarkanKeseimbanganGaya pada Fluida dalamPipa
221
RisalahLokakaryaKomputasidalamgainsdan TeknologiNuklir XN, Juli 2003
Suku terakhir daTi persamaan (14) adalah perubahan momentum aliran masuk clankeluar, yang mana untuk pipa lurus menjadi hilang karena persamaankontinuitas. Suku kedua merupakan gaya karena gravitasi, biasanya digabungkan dengan gaya tekanan dengan memperkenalkan variabel tinggi tekanPhyaitu: Ph =p+p.g.Z,
sehinggapersamaanmomentumuntukpipa lurus adalah:
-dph
4 D
(16) dx Berawaldan persamaan kontinuitas(4), yangkita tulis ulang denganpersamaan :7)menjadi -=-'l",
av ax
= 0,
bersama-samadengan persamaanmomentum(16), kita terapkan teknik elemen hingga dengan formulasi akar terkecil, yaitu dengan cara mengasumsikan fungsi Ph dan V untuk setiap elemen pipa sebagai fungsi linear dari variabel bebas x sebagaiberikut: Ph
= Ph/N
v
= V IN
./N
+ PhOUT
. OUT
':II IN + V OUT
\f OUT
Substitusipersamaan(18) kedalampersamaan (17), lalu denganformulasi akar terkecildidapat: Tf ~N
=
y OUT
=
V
(20)
e'
yang tidak lain adalahkontinuitas dalam elemenpipa. Dengan menerapkan formulasiakarterkecilpadapersamaan(17) clan(16), didapat:
Le
222
1
]{ Pl
-l
1
pz
4
'l'"e(Ve) -'l'"e(Ve)
}.
Model Matematik dan Numerik daTi Aliran Fluida dalamjaringan Pipa Berbasis MEH
(Ade Jamal,et.al.»
Persamaanini diturunkan denganmenerapkansecaraimplisit persamaan(20) clan asumsi bahwa faktor gesekan Te penyebab turunnya tekanan dalam pipa tergantungpada kecepatantapi tidak padajarak, yang mana hal ini berlaku untuk elemenpipa dengandiameterkonstan.Dua buahpersamaan dalammatriks (21) saling ketergantungan sehinggahanyasatupersamaan yangdipakaiyaitu:
}
[1 -ll!Pl ~P2
= 4.Le ~l"e(Ve
e
Persamaan(20) clan(22) merupakanpersamaanuntuk satuelemenpipa. Untuk kasus aliran laminar, persamaanini sarna dengan model elemen hingga yang tradisionalditurunkandenganformulasi Galerkin[9], yaitu: 1l".D4 128.p"Le
1 -1
~1]{;~}={~~}
di manaPi adalahtekananhidrostatisdi tiap ujung elemen,clanQi adalahdebit yang masuk(positif) atau keluar (negatif) di ujung elemen.Pada saatpenggabungan (assembly)semuamatriks elemenke dalam matriks sistemjejaring pipa, persamaan kontinuitaspadatiap titik simpul N di terapkanpadasukukananpersamaan (23). Persamaan kontinuitasitu adalah:
di mana qN adalah debit aliran yang keluar daTi sistem pipa di titik simpul N dan jumlah pipa yang berhubungan dengan titik simpul adalah m. Dalam metode elemen hingga tradisional ini tekanan Pi didetinisikan sebagai variabel primer, yang biasanya harus dicari, dan debit aliran keluar/masuk qi adalah variabel sekunder, yang biasanya diketahui. Jumlah variabel yang dicari sarna denganjurnlah persamaannya yaitu sarna denganjumlah titik simpul. Untuk kasus yang tidak dibatasi pada aliran laminar, gesekan 'ie mungkin merupakan fungsi tak-linear daTi debit Qe (atau kecepatan fluida Ve=Q/Ae dalam pipa), sehingga persamaan sistem jejaring pipa juga tidak linear dengan variabel yang harus dicari selain tekanan pada setiap titik simpul Ph juga kecepatan aliran V atau debit aliran Q dalam pipa. Sistem persamaanumum untukjejaring pipa yang terdiri M elemen pipa dan N titik simpul adalah: sejumlah M persamaan (25) dan sejumlah N persamaan(24). Sistem persamaan ini mengatur sejumlah N variabel tekanan Ph dan
RisalahLokakaryaKomputasidalamSainsclanTeknologi Nuklir XN, Juli 2003
sejumlah M variabel debit aliran Q:
~K] [T(Q)nJWh} }= { {O}} l[O] [q Jl {Q} {q}
.
di mana matriks perbedaantekanan [Ke] berukuran(M*N), matriks gesekan [Te(Q)] adalah matriks diagonal berukuran (M*M) clan matriks kontinuitas [C] berukuran(N*M). Perlu dicatatbahwa sistemmatriks persamaan(25) tidak simetris sepertibiasanyapadasistemmatriksMEH yangtradisional. Setidaknyasatu persamaandalam matriks kontinuitas [C] bergantungsecara linear denganyang lainnya, karenaitu minimal satudebit keluarq dijadikan variabel menggantikan satu variabel Ph yang diketahui nilainya sebagai kondisi batas (boundarycondition)agarpersamaan(25) diatasdapatdiselesaikan:
di mana variabel tekanan Ph diuraikan menjadi PhAyang diketahui dan PhByang dicari serta matriks perbedaantekanan [K] diuraikan menjadi [KA] dan [KB] masing-masing sesuai denganvariabel tekanan Ph yang diketahui dan dicari. Jika semua tekanan Ph di semua titik simpul diketahui, maka sistem persamaanmatriks (26) menjadi:
~T(Q)] [0]If{!2J}= {{[-KJ {Ph}
L[C]
[-J]Jl{q}
{O}
}
Untuk kasusumum di mana gesekanmernpakanfungsi tidak linear daTivaribel debit Q, maka persamaan(25), (26) clan (27) menjadi tidak linear clan harns diselesaikandenganmetodeiterasi sepertipadametodeNewtonRaphson.
224
Model MatematikdanNumerik daTiAliran FluidadalamjaringanPipaBerbasisMEH
(Ade Jamal,et.al.»
KESIMPULAN DAN REKOMENDASI Suatupendekatanyang lain untuk analisis aliran fluida tunak clan tak-mampu mampattelah dikembangkandenganteknik ElemenHingga dengan formulasi akar terkecil. Dengan metode MEH ini debit dalam pipa dapat langsung dihitung bersamaandengan variabel lain di titik simpul yaitu tekanan atau debit yang keluar/masukdi titik tersebut.Metode elemenhingga yang menggunakanformulasi Galerkinsebagaianalogidari model elemenhinggadari strukturbatang,penghitungan debit dalampipa dilakukansebagaianalisislanjut (post-analysis).Debit aliran dalam pipa denganMetode Cross-Hardyjuga didapatsetelahiterasikoreksidebit. Dibanding dengan Metode Matriks Steffer, keunggulanmetode ini adalah terhindamya kemungkinan sistem matriks yang tidak sempurna (ill-conditioned matrix)dari matrikskontinuitas. Model elemenhingga untuk elemenpipa ini dapatdikembangkanlebih lanjut untuk komponen-komponen jaringan pipa lainnya sepertiuntuk elbow, klep, strainer clan komponenlainya dengan tetap konsistenpada sistem pengembanganelemen matriksnya.Selain itu algoritma clan prosedursolusi dari sistempersamaanmatriks dapatdikaji lebih jauh untuk mendapatkanefisiensidari segi waktu perhitunganclan jumlah memory komputer yang diperlukan denganmemperhatikanbentuk sistem matriksyangkhususdari persamaan(26) clan(27).
DAFTARPUSTAKA WHITE, F.M., ViscousFluid Flow, 2ndEd. MCGraw-Hill (1991)
2. BENEDICT, R. P., Fundamentalof PipeFlow, JohnWiley andSons(1980) 3
KODOATIE, R. J., Hidrolika Terapan-Aliran pada saluranTerbuka dan Pipa, PenerbitAndi, Yogyakarta,(2002)
4. STREETER,V.L., WYLIE, E.B., Fluid Mechanics, 8ili Edition, McGrawHill (1985) 5 DHA1T, G., Finite Element Modeling of Fluids in Computational Fluid Dynamics, Lecture Series 1992-04, yon Karman Institute for Fluid Dynamics (1992)
6. CONNOR, J.J., BREBIA, C.A., Finite Element Technicquesfor Fluid Flow, Newnes-Butterworths, London(1976)
225
RisalahLokakaryaKomputasidalamSainsdaDTeknologiNuklir XIV, Juli 2003
7. FERZIGER, J. H., PERIC, M., ComputationalMethods for Fluid Dynamics, Springer2ndPrinted (1997) 8. ARASU, K., Analysis of Pipe Networks, EngineeringDevelopmentDivision,
hldira GandhiCentrefor Atomic Research(1993)
9. REDDY, J.N., An Introdcution to The Finite Element Method, McGrawHil1 (1985)
DISKUSI
HUDIHASTOWO Apa ada rencana untuk mengembangkansoftware yang sudah dibuat untuk menyelesaikanmasalahfluida dalam2 phasa?Permasalahan ini banyak dijumpai di bidang teknik nuklir, tetapi memang disadari tidak mudah. Bagaimana kita menyelesaikan3 persamaankontinuitas momentumclan energi dalam 2 phaseyang berbeda?
ADE JAMAL Rencana2 phasebarn bisa dimulai jika model untuk rasegas (incompressible)telah selesai.Tahapawal ini barn sampaialirancompressible.
ENDANG ROSADI Apakah pengaruh temperatur fluida dilibatkan dalam formulasi/pemodelan ini?
ADE JAMAL Pemodelanini hanya dilaksanakanuntuk air clanminyak saja.Denganasumsiyang dipakai untuk aliran incompressible(tidak mampumampat),makafungsi temperatur tidak dilibatkan. Untuk tahap pengembanganmodel selanjutnya di mana aliran dihitung, makatemperaturharnsdiperhitungkan.
226
Model Matematika dan Numerik daTi Aliran Fluida dalam Jaringan Pipa Berbasis MEH
(Ade Jamal,et.al.
Dalam sebuahjaringan pipa, biasanyaada sebuahkomponenyang mengandung tekanan clan volume besar. Apakah komponen seperti itu sudah termodelkan denganmodelyangsekarang? 2. Model solusialiran fluida, untuk setiaptitik ada variabeltekanandan debit volume yang salingbergantungan,tapi tidak dipengaruhinilai besarannya.Model numerik selalu membutuhkansalah satu nilai daTi variabel tekanan atau debit volume. Variabelyangbelumdiketahuinilainya dihitungdengankomputasi. M. SYAMSA ARDISASMITA Apa alasanSaudaramencari penyelesaianpemodelanmatematikadenganmetode Kirchoff yangsederhanadibandingkandenganMetodeElemenHingga?
ADEJAMAL Metode Kirchoff tidak kita gunakan dalam metode solusi yang kita kembangkan. Metode Kirchoff ini dinaikkan sebagai gambar "state of the art" teknologi yang digunakan oleh design manual dan textbook yang ada. Model ini dikembangkan murni Metode Elemen Hingga dengan formulasi Kuadrat Terkecil.
227
Home