1. Home
  2. Archives
  3. Vol 6 (1987) Issue 1&2
  4. Articles

Komputasi Numerik Aerodinamika Kendaraan Dengan Metoda Identitas Green

Abstract

Suatu program komputasi untuk perhitungan aerodinamika benda sebarang dalam aliran potensial telah disusun dan dikembangkan dengan dasar metoda panel elemen kwadrilateral dan dengan menggunakan prinsip identitas Green, dengan salah satu tujuan penggunaan perhitungan aerodinamika kendaraan bermotor. Program yang telah disusun diuji pada benda-benda dengan geometri sederhana guna menilai keandalan program serta tingkat ketelitian yang dapat dicapai. Selanjutnya dilakukan perhitungun aerodinamik pada benda dengan tahanan minimum dan kendaraan darat lain (secara dwi-matra). Dari hasil perhitungan ini terlihat bahwa program komputasi yang telah disusun mampu meramalkan sifat aerodinamika (potensial_kecepatan dan distribusi tekanan) benda-benda dengan hasil yang cukup baik.

1. PENDAHULUAN

Pengaruh aerodinamika pada kendaraan, terutama pada kecepatan yang reatif tinggi (di atas 75 km/jam), memegang peranan praktis yang penting, ditinjau dari persentase daya yang diperlukan guna mengatasi tahanan aerodinamik, dan pengaruh aerodinamika pada stabilitas kendaraan [1][2].

dinamika Seiak komputasional mulai diperkenalkan pada akhir dekade 1950-an dan awal dekade 1960-an, sejalan dengan perkembangan efisien dan teknologi komputer yang berkecepatan tinggi dan metoda komputasi, berbagai penerapan praktis peramalan aerodinamika guna berkembang. Kebutuhan yang mendesak guna terselesaikannya berbagai masalah dinamika fluida -terutama aerodinamikajuga membantu pengembangan dinamika fluida komputasional secara pesat dan luas. Dinamika fluida komputasional pemecahan herbagai memungkinkan persoalan dinamika fluida yang belum dapat dipecahkan secara analitik, antara lain pola aliran yang sederhana yang vang menyangkut geometri aliran sebarang, yang merupakan suatu masalah yang cukup rumit, sebelum era dinamika fluida komputasional. Pada saat ini, dinamika fluida komputasional dapat dikatakan telah memasuki generasi ketiga [3], dan merupakan saat yang baik sekali dinamika fluida/aerodinamika teroretik, eksperimental maupun sebagai penghubung keduanya. Hubungan antara ketiganya dapat dilihat pada gambar 1.

penyelesaian cara Salah satu masalah aerodinamika dengan komputasi numerik adalah penggunaan metoda panel, yang telah dipergunakan secara luas oleh banyak peneliti, dan umumnya telah memberikan hasil yang baik [4][5][6][7][]. Metoda panel digunakan pada penelitian ini untuk perhitungan aerodinamika benda/kendaraan dengan menerapkan metoda identitas merupakan Pekerjaaan ini Green. kelanjutan serta pengembangan metoda yang telah dikembangkan dari oleh dan Grose[], & Djojodihardjo[9][10]. Hasil yang diperoleh diharapkan dapat memberikan gambaran awal tentang sifat aerodinamik suatu benda/ kendaraan serta dapat membantu perancangan pengujian terowongan angin:

Gambar 1 Hubungan antara Dinamika Fluida Teoritik, Eksperimental, dan Komputasional

* Jurusan Mesin ITB

2. PERSAMAAN LAPLACE

Medan aliran sekitar suatu benda/kendaraan yang bergerak relatif terhadap lingkungan fuida sekelilingnya, dalam pendekatan dasar dipandang tak viskos, inkompresibel dan irotasional. Dengan demikian, kita dapat mendefinisikan potensial kecepatan sebagai

\[V = \nabla \Phi \tag{1}\]

Karena di sini berlaku persamaan kontinuitas

\[\nabla \cdot \mathbf{V} = 0 \tag{2}\]

maka potensial kecepatan \(\Phi\) memenuhi persamaan Laplace

\[\nabla^2 \Phi = 0 \tag{3}\]

3. KONDISI BATAS

Pada pendekatan aliran tak viskos. kondisi batas yang harus dipenuhi agar persamaan (3) dapat dipecahkan adalah kondisi batas aliran menyinggung permukaan, yaitu:

\[\mathbf{y} \cdot \mathbf{n} = 0 \tag{4}\] dengan n - vektor normal satuan (tegak lurus) pada permukaan benda. Dengan menggunakan variabel potensial kecepatan, kondisi batas ini dapat dituliskan sebagai:

\[\nabla \Phi \cdot \mathbf{n} = 0 \tag{5}\] yang dikenal pula sebagai kondisi batas. Neumann.

Kondisi batas yang lain menyatakan kondisi radiasi, yaitu bahwa fungsi potensial & bernilai konstan (yang dapat ditetapkan sama dengan nol) pada titik medan aliran yang terletak sangat jauh dari permukaan benda.

Persamaan Laplace (3) merupakan persamaan diferensial parsial eliptik, yang memerlukan penyelesaian secara simultan pada semua titik dalam medan aliran. Penyelesaian persamaan (3) yang dikenai syarat batas Neumann permukaan benda dan syarat batas radiasi di tak berhingga dapat sebagai persamaan integral dirumuskan Fredholm. Berbagai literatur [11][12] menunjukkan bahwa agar dapat dipecahkan tanpa kesulitan kestabilan numerik, persamaan integral Fredholm tersebut berbentuk jenis kedua. Metoda hendaknya identitas Green yang dirumuskan oleh Bristow dan Grose[8] memenuhi persyaratan tersebut.

4. PERUMUSAN PERSAMAAN KERJA DALAM BENTUK PERSAMAAN INTEGRAL IDENTITAS GREEN

Karena persamaan Lapace bersifat linear, maka setiap jawab dapat disusun dari berbagai jawab elementer yang juga memenuhi persamaan Laplace: prinsip superposisi berlaku.

Oleh karena itu, guna kemudahan, potensial aliran \(\Phi\) di sekitar benda dapat dianggap sebagai superposisi dari potensial aliran bebas (aliran dasar) dan potensial gangguan \(\phi\):

\[\Phi = V \times + \phi \tag{6}\]

Karena sifat persamaan Laplace yang linear, berdasarkan prinsip superposisi, secara mudah dapat ditunjukkan bahwa \(\phi\) juga memenuhi persamaan Laplace:

\[\nabla^2 \phi = 0 \tag{7}\]

Medan kecepatan aliran juga dapat dibedakan atas vektor kecepatan mutlak aliran dan vektor kecepatan gangguan:

\[V = \nabla \Phi = V + v \tag{8}\]

\[\mathbf{v} = \nabla \phi \tag{9}\]

Untuk memudahkan, kerangka acuan dipilih berimpit dengan benda sedemikian sehingga aliran tak terganggu mempunyai kecepatan V dalam arah sumbu-x positif, seperti diperlihatkan pada Gambar 2.

Gambar 2 Sistem Koordinat Setelah V dan v diketahui, tekanan dapat dihitung dengan menggunakan persamaan Bernoulli:

\[C_{p} = \frac{p - p}{\frac{1}{2} \rho V^{2}} = 1 - \left\{ \frac{V}{V} \right\}^{2}\] (10)

Gaya aerodinamik dapat segera diperoleh dengan integrasi tekanan pada seluruh permukaan benda, yang kemudian dapat diuraikan dalam komponennya: gaya tahanan dalam (melawan) arah gerakan dan gaya angkat tegak lurus pada arah gerakan.

\[D = \frac{1}{2} \rho V^2 \int \int C_{\rho} n \cdot i dS \qquad (11)\]

\[L = \frac{1}{2} \rho V^2 \int \int C_p n \cdot k dS \qquad (12)\]

Demikian pula momen aerodinamik dapat langsung diperoleh pula dari:

\[M = \frac{1}{2} \rho V^2 \int \int C_p(n \cdot iz + n \cdot kx) dS\] (13)

Jawab dari persamaan Laplace (3) dapat dituliskan sebagai integral dari fungsi Green sebagai berikut [8][13]

\[\phi = \iint_{S_{B+W}} \left\{ -\frac{\nabla \phi \cdot \mathbf{n}}{4\pi \mathbf{r}} + \frac{\phi}{4\pi} - \frac{\partial}{\partial \mathbf{n}} \cdot (\frac{1}{\mathbf{r}}) \right\} dS\] (14)

yang merupakan persamaan identitas Green yang ketiga, dan memiliki bentuk persamaan integral Fredholm jenis kedua. S merupakan permukaan batas antara fluida dengan benda pejal yang bergerak dan titik-titik di tak berhingga. Dengan menidentifikasikan \(\phi\) sebagai distribusi doblet dan \(\nabla \phi\).n sebagai distribusi sumber pada seluruh permukaan, dapat ditunjukkan [5][8][14], bahwa

\[\phi = \int \int \left\{ -\frac{\sigma}{4\pi r} + \frac{\mu}{4\pi} \frac{\partial}{\partial n} \left( \frac{1}{r} \right) \right\} dS\] (15)

dengan S permukaan benda / kendaraan berikut ulakannya. Dapat ditun jukkan bahwa untuk ulakan yang tanpa ketebalan, integral dari suku pertama ruas kanan persamaan (15) sama dengan nol [8][14]. Dalam persamaan ini:

σ - kekuatan sumber per satuan luas pada permukaan pembatas pada titik sebarang Q

\(\mu\) - kekuatan doblet per satuan luas pada permukaan pembatas pada titik sebarang Q

dS- elemen permukaan diferensial pada Q

\(r - |x_p - \xi_Q|\), jarak antara titik medan P terhadan titik sumber Q

terhadap titik sumber Q. n - vektor satuan normal (tegak-lurus) pada permukaan S di titik Q

Bentuk (14) atau (15) dapat dipergunakan sebagai dasar penyusunan persamaan kerja /

5. PERSAMAAN KERJA

Persamaan Laplace (3) untuk \(\phi\) dengan syarat batas (4) pada permukaan benda S dipenuhi oleh \(\phi\) yang dinyatakan dengan persamaan integral Fredholm jenis

kedua (14) yang merupakan identitas Green. Di sini, \(\nabla \phi\) n pada integran pertama pada ruas kanan , yang memiliki pengertian fisik distribusi kekuatan sumber (sama dengan \(\sigma\) pada persamaan (15)) merupakan syarat batas Neumann yang diketahui, karena

\[V = \nabla \Phi = V + \nabla \phi\]

\[\mathbf{V}_{\mathbf{n}} = \nabla \Phi \cdot \mathbf{n} = \mathbf{V} \cdot \mathbf{n} + \nabla \phi \cdot \mathbf{n}\]

sehingga

\[\nabla \phi \cdot \mathbf{n} = \nabla - \nabla \cdot \mathbf{n} \tag{16}\]

Syarat batas (4) menyatakan bahwa

\[\mathbf{V}_{\mathbf{n}} = \nabla \Phi + \mathbf{n} = 0 \tag{17}\]

Jadi, dari (16), diperoleh:

\[\sigma = \nabla \phi \cdot \mathbf{n} = -\mathbf{V} \cdot \mathbf{n} \tag{18}\]

Bila (18) dimasukkan ke dalam persamaan (14), diperoleh persamaan kerja

\[\phi = \iint_{S_{B+V}} \left\{ -\frac{V \cdot n}{4\pi r} + \frac{\phi}{4\pi} \frac{\partial}{\partial n} \left( \frac{1}{r} \right) \right\} dS\] (19)

yang akan dipergunakan sebagai dasar metoda panel yang akan dikembangkan di sini. Dengan menggunakan aliran fiktif di dalam benda yang menghasilkan persamaan integral untuk potensial yang seupa dengan (19), tetapi ruas kiri sama dengan nol, dan integran pertama pada ruas kanan juga nol karena di dalam airfoil tidak ada aliran , maka diperoleh persamaan kerja

\[\phi = \iint_{S_{\mathbf{n}}} \left\{ -\frac{\mathbf{v} \cdot \mathbf{n}}{4n\mathbf{r}} + \frac{\phi}{4n\mathbf{\partial n}} \frac{\partial}{\partial \mathbf{r}} \right\} dS\] \[+ \iint_{S_{\mathbf{v}}} \left\{ \frac{\Delta \phi}{4n\mathbf{\partial n}} \cdot \frac{\partial}{\partial \mathbf{n}} \cdot \left(\frac{1}{\mathbf{r}}\right) \right\} dS \qquad (20)\]

Di sini \(\phi\) menyatakan besarnya distribusi doblet pada permukaan \(S_{_{\rm B}}\) sedangkan \(\Delta\phi\) beda kekuatan doblet di atas dan di bawah permukaan atas ulakan \(S_{_{\rm W}}\).

6. PERUMUSAN METODA PANEL

Permukaan benda diwakili oleh sekumpulan panel; pada tiap panel kekuatan sumber dan doblet dianggap konstan. Hal ini merupakan pendekatan orde ke nol. Kekuatan sumber pada tiap panel sama dengan nilai V<sub>o</sub>n pada titik

beratnya, sedangkan kekuatan doblet pada tiap panel, \(\phi\), merupakan besaran yang tidak diketahui dan sama dengan besarnya potensial gangguan di titik tersebut. Dengan membagi permukaan benda \(S_{\rm B}\) menjadi n buah panel, akan diperoleh n buah kekuatan doblet yang harus dicari dengan n buah kondisi batas, sehingga persamaan (20) diubah menjadi himpunan n buah persamaan simultan. Jenis panel yang dipergunakan adalah permukaan datar, dengan bentuk kwadrilateral. Perincian pendiskritan persamaan (20) ini diperinci pada pustaka [10]. Persamaan Fredholm jenis kedua (20) diubah menjadi sitem persamaan simultan berikut ini:

\[\phi_{i} = \sum_{j=1}^{n} - V_{co} \cdot n_{j} C_{ij} + \sum_{j=1}^{n} \phi_{j} A_{ij}\] (21)

dengan

\[C_{ij} = \frac{1}{4\pi} \iint \frac{1}{r(x_i, \xi_j)} dS(\xi_j) (22)\]

\[A_{ij} = \frac{1}{4\pi} \int \int \frac{\phi}{\phi_j} \frac{\partial}{\partial n_j} \left\{ \frac{1}{r(x_i, \xi_j)} \right\} dS(\xi_j)\] (23)

Bila didefinisikan

\[B_{i} = \sum_{j=1}^{n} - V_{\infty} \cdot n_{j} C_{ij}\] (24)

maka persamaan (21) dapat ditulis kembali menjadi:

\[\phi_i = B_i + \sum_{j=1}^n \phi_j A_{ij}\] (25)

atau

\[\tilde{\mathbf{B}}_{i} = \sum_{j=1}^{n} \phi_{j} \tilde{\mathbf{A}}_{ij} \tag{26}\] dengan

\[\bar{A}_{ij} = A_{ij} - \delta_{ij} \tag{27}\]

\[\bar{B}_{ij} = -B_{ij}' \tag{28}\]

Di sini \(\delta_{ij}\) adalah delta Kronecker. Evaluasi dari koefisien pengaruh \(A_{ij}\) dan \(C_{ij}\) diperinci pada [4],[14] dan [15].

7. GARIS BESAR PROGRAM KOMPUTASI

Prosedur program komputasi yang telah disusun mencakup empat langkah utama, yaitu:

1) Geometri benda didefinisikan dengan penentuan letak titik sudut panel permukaan, yang dinyatakan dalam sistem acuan global.

2) Penghitungan vektor penyusun panel dan berbagai atributnya, seperti koordinat titik kontrol dan vektor normal satuan. Komponen vektor satuan normal dan tangensial juga merupakan matriks transformasi dari koordinat lokal ke koordinat global.

3) Penghitungan koefisien pengaruh. Penghitungan ini didahului oleh transformasi vektor penyusun panel dan atributnya dari koordinat global ke koordinat lokal (panel).

4) Penyelesaian persamaan simultan linear, dengan metoda eliminasi Gauss, disertai kondensasi tumpuan, agar kesalahan pembulatan sekecil mungkin dan diperoleh ketelitian yang cukup baik.

5) Penghitungan berbagai besaran aliran pada permukaan seperti kecepatan, tekanan dan sebagainya. Penghitungan koefisien tekanan di sini dilakukan dengan cara diferensiasi numerik, dengan lebih dahulu melakukan penyesuaian lengkungan (curve fitting) dengan metoda kesalahan pangkat dua terkecil.

8. PEMAPARAN HASIL DAN PEMBAHASAN

Metoda panel dengan penerapan identitas Green yang telah disusun lebih dahulu diterapkan pada berbagai benda dengan geometri sederhana yang sudah dikenal jawabnya guna menilai keandalan dan ketepatannya.

Pada Gambar 3 dan 4 diperlihatkan perbandingan hasil komputasi dengan jawab eksak untuk potensial kecepatan serta koefisien tekanan aliran fluida di sekitar silinder. Pembandingan ini menunjukkan hasil yang amat baik.

25

Gambar 3 Distribusi Potensial Kecepatan Sekitar Silinder

0

Gambar 4 Distribusi Tekanan Sekitar Silinder

Selanjutnya diperlihatkan hail yang diperoleh untuk aliran sekitar elips, airfoil NACA 0012, airfoil Joukowki, seperti terlihat pada Gambar 5, 6 dan 7.

3

Gambar 5 Distribusi Potensial Kecepatan Sekitar Elips

5

Gambar 6 Distribusi Kecepatan Sekitar Airfoil NACA 0012

Untuk aliran di sekitar bola (aliran tiga dimensi, Gambar 3) terlihat bahwa penambahan jumlah panel memperbaiki hasil komputasi secara konsisten. Selanjutnya diperlihatkan

8

Gambar 7 Hasil Komputasi untuk Profil Joukowski

hasil komputasi potensial kecepatan aliran sekitar spheroid, spheroid dalam arah melintang arah airan dan ellipsoid ( Gambar 9, 10 dan 11), yang

11

Gambar 8 Distribusi Potensial Kecepatan Sekitar Bola

13

Gambar 9 Distribusi Potensial Kecepatan Sekitar Spheroid

dibandingkan dengan jawab ekaknya. Hasil yang diperoleh cukup baik, walaupun penambahan jumlah panel ternyata belum memperbaiki ketelitian yang diperoleh. Dengan demikian perlu dilakukan evaluasi

0

Gambar 10 Distribusi Potensial Kecepatan Sekitar Spheroid dalam Arah Melintang

2

Gambar 11 Distribusi Potensial Kecepatan Sekirar Elipsoid

kembali pada program komputasi. hakekat pendekatan teoretik yang dipakai dan penerapannya pada program komputasi.

Selanjutnya metoda ini diterapkan pada dua buah model kendaraan darat, serta benda dengan tahanan aerodinamik minimum (yang dapat dipergunakan sebagai dasar / rujukan / base-line pada perancangan kendaraan darat). Secara kwalitatif diperlihatkan kecenderungan yang baik, seperti diperlihatkan pada Gambar 12, 13, 14 dan 15.

6

Gambar 12 Distribusi Tekanan Sekitar Kendaraan Versi 1

8

Gambar 13 Distribusi Tekanan Sekitar Kendaraan Versi 2

10

Gambar 14 Distribusi Potensial Kecepatan Sekitar Body of Minimum Drag

12

Gambar 15 Distribusi Tekanan Sekitar Body of Minimum Drag

Pada Gambar 16 diperlihatkan potensial perbandingan antara nilai kecepatan untuk aliran sekitar silinder yang diperoleh menurut metoda yang disusunt dengan metoda panel yang menggunakan sumber saja. Hasil yang diperoleh menun jukkan bahwa metoda identitas Green memiliki ketelitian yang lebih baik dari metoda sumber.

0

Addistrict Control (A) the property of the con-

1

Gambar 16 Distribusi Potensial Kecepatan Sekitar Silinder

9. KESIMPULAN DAN RENCANA SELANJUTNYA

Suatu program komputasi telah dikembangkan guna penentuan besaran - besaran aerodinamik sekitar benda / kendaraan, dengan memakai metoda panel dan penerapan metoda identintas Green ketiga. Pengujian metoda ini dilakukan pada benda benda dengan geometri sederhana dengan hasil yang cukup baik, dengan meninjau ketelitian yang diperoleh maupun konvergensinya. Penerapan pada model kendaraan darat juga telah memperlihatkan kecenderungan yang baik. Pada program yang disusun di sini masih diperlukan pengujian program yang menyertakan jejak / ulakan, untuk menghitung gaya angkat yang dialami benda. Pola program demikian dapat dilihat pada karya terdahulu untuk serupa [9]. Peningkatan metoda ketelitian dan efektifitas prosedur komputasi juga akan disempurnakan dengan menggunakan pendekatan orde yang lebih tinggi, misalnya penggunaan permukaan lengkung atau distribusi singularitas per panel dengan orde yang lebih tinggi.

References

  1. H.Djojodihardjo, Pertimbangan aerodinamika untuk penghematan penggunaan bahan bakar pada kendaraan bermotor, Konferensi kedua, IATO, Oktober 1979.
  2. H.Djojodihnrdjo, Aerodynamic Apects of Pasenger Vehicles - Review and Analysis, The Third International Pacific Conference on Automotive Engineering, Jakarta, November 1985.
  3. John D,Anderson. Introduction to Computational Fluid Dynamics. Lecture Series l986-01. Von Karman Institute for Fluid Dynamics. 1986.
  4. J.L.Hess dan A.M.O.Smith, Calculation of Potential Flow About Arbitrary Bodies. Progre in Aeronautical Sciences. Vol.8. Pergamon Press, New York. 1966.
  5. N.Djojodihardjo dan S.E.Widnal1, A Numerical Method for the Calculation of Nonlinear. Unsteady Lifting Potential Flow Problems, AIAA Journal, Vol.7, Oct.l969.
  6. B.Hunt. The Panel Method for Subsonic Aerodynamic Flows: A Survey of Mathematics Formulations and Numerical Models and an Outline of the New British Aeropace Scheme. in Computational Aerodynamics, Von Karman Institute for Fluid Dynamics, 1978.
  7. C.C.Kuo and L.Morino, Steady Subonic Flow Around Finite Thickness Wings. NASA CR"2616. National Aeronautics and Space Administration. Washington, D.C., 1978.
  8. D.R.Bristow and G.G.Grose, Modification of Douglas Neumann Program to Improve the Efficiency of Predicting Component Interference and High Lift Characteristics. NASA Contractor Report 3020. National Aeronautics and Space Administration, Washington, D.C., 1978.
  9. H Djojodihardjo, Development of Aerodynamic Computations Routine for Road Vehicles Using Panel Method, Proc.The Second International Pacific Conference on Automotive Engineering, Tokyo. November 1983.
  10. H.Djojodihardjo. H.Suhartono, dan M.Karnadi. Penggunaan Identitas Green dengan Metoda Panel Guna Komputasi Aerodinamik ubonik Tiga-Dimensi, Aeronotika dan Astronotika - Indonesia Vol.1. no.1. 1986.
  11. S.Mikhlin. Integral Equations, The MacMillan Company. New York, 1964.
  12. E.Martensen. Potential Theorie. Teubner Verlag. Stuttgart, 1963.
  13. O.D.Kellog. Foundations of Potential Theory. Dover Publications, Inc. 1953.
  14. H Djojodihardjo, A Numerical Method for the Calculation of Nonlinear Unsteady Lifting Potential Flow Problems. Doctor of Science Thesis, Massachusetts Institute of Technology, 1965.
  15. S.Bunyamin. H.Djojodihardjo, Sudarmipto dan H.Suhartono. Komputasi Aerodinamik dengan metoda Green, Seminar Komputasi dan Identifikasi Sistem Aero-Elasto-Mekanik, PAU Ilmu Rekayasa ITB, Oktober 1987.