MATEMATIKA
TEKNIK KIMIA
Dr. Ir. Setijo Bismo, D.E.A.
Dr.rer.nat. Ir. Yuswan Muharam, M.T.
SILABUS
Pendahuluan
Formulasi problem fisikokimia
Teknik penyelesaian model persamaan
diferensial biasa (PDB)
Teknik penyelesaian model persamaan
diferensial parsial (PDP)
REFERENSI
Applied Mathematics and Modeling for Chemical
Engineers, Rice, 1995
Numerical Methods for Chemical Engineers with
MATLAB Applications, Constantinides, 1999.
Numerical Methods for Engineers and Scientists,
2nd Edition, Hoffman, 2001
Applied Numerical Methods Using Matlab, Yang,
2005
Numerical Analysis Using MATLAB and
Spreadsheets.2ed Ed, Karris, 2004
EVALUASI
UTS = 20 %
UAS = 30 %
Tugas = 30 %
Proyek = 20 %
PENDAHULUAN
Dr.rer.nat. Ir. Yuswan Muharam, M.T.
Departemen Teknik Kimia
Fakultas Teknk
Universitas Indonesia
DEFINISI MODEL (TERMINOLOGI)
“Sebuah objek M (benda, sistem fisika atau
kimia, atau proses) adalah model apabila
terdapat analogi antara objek M dan objek
lain O sehingga kesimpulan mengenai O
dapat dibuat”.
DEFINISI MODEL (TERMINOLOGI)
Model M
Representasi objek O;
Taksiran objek O yang diisolasi dari seluruh
realitas,
Menggambarkan kenyataan atau bagian dari
kenyataan.
Dapat disederhanakan menjadi bagian dari
kenyataan jika perlu kesimpulan tertentu saja.
DEFINISI MODEL (TERMINOLOGI)
Keterbatasan analogi model M dan objek O
Keterbatasan kesesuaian fungsi,
Keterbatasan lesesuaian struktur dan perilaku,
Keterbatasan akurasi.
Model M dan objek O boleh berbeda skala.
Hasil model bagus apabila variabel dan
fenomena pentingnya direpresentasikan
secara benar dalam konteks atau investagi
tertentu.
DEFINISI MODEL (TERMINOLOGI)
Analogi antara model M dan objek O dapat
dibuat dalam bentuk persamaan matematis.
Model matematis menggambarkan
seperangkat persamaan aljabar dan/atau
diferensial dan/atau integral yang digunakan
untuk menjelaskan perilaku objek O.
TUGAS CHEMICAL ENGINEER
Mengoperasikan dan mengoptimalkan proses
yang ada;
Merancang pabrik baru dan memodifikasi
pabrik yang ada.
APLIKASI MODEL MATEMATIS
DI INDUSTRI KIMIA
Percobaan
Simulasi
Analisis sensitivitas
Kendali dan operasi
Optimisasi
Eksplorasi
KETERBATASAN MODEL
MATEMATIS
1. Jenis, jumlah serta keakuratan data;
2. Perkakas matematis;
3. Interpretasi hasil model.
INTERPRETASI HASIL MODEL
PENYUSUNAN DAN
KLASIFIKASI MODEL
Dr.rer.nat. Ir. Yuswan Muharam, M.T.
Departemen Teknik Kimia
Fakultas Teknk
Universitas Indonesia
PENYUSUNAN MODEL
MATEMATIKA
Penyusunan model matematika adalah
pengesetan seperangkat persamaan
matematika.
Persamaan matematika adalah hubungan
antara variabel proses.
TAHAP-TAHAP PEMODELAN
1. Formulasi persoalan, pengumpulan objektif
dan kriteria keputusan;
2. Pengamatan terhadap proses dan
klasifikasinya untuk membagi proses
menjadi beberapa subsistem (elemen
proses);
3. Penentuan hubungan antara subsistem;
4. Analisis variabel dan hubungan antar
variabel pada setiap elemen proses;
TAHAP-TAHAP PEMODELAN
5. Pembentukan persamaan matematika
dengan menggunakan variabel dan
parameter; Pengumpulan data;
6. Pengamatan representasi proses oleh
model; perbandingan hasil simulasi dengan
data proses nyata;
7. Instalasi model; interpretasi dan
pemeriksaan hasil.
TAHAP-TAHAP PEMODELAN
8. Analisis sensitivitas model untuk
mengidentifikasi parameter yang
berpengaruh kuat dan lemah terhadap
respons model;
9. Penyederhanaan model.
10. Tahap 4 – 9 diulang, sampai interpretasi
hasil model sesuai dengan kriteria objektif
dan solusi yang diharapkan.
KEGUNAAN MODEL
Untuk memformulasikan fenomena fisika dan
fisikokimia, yaitu perpindahan panas,
perpindahan massa dan perpindahan
momentum, serta reaksi kimia di dalam
sistem homogen dan heterogen.
Untuk mendesain operasi perpindahan
massa, menghitung penukar panas,
merekayasa reaksi kimia, dan
mengendalikan proses.
KLASIFIKASI
MODEL MATEMATIKA
MODEL BERDASARKAN
PRINSIP FISIKOKIMIA
Digunakan untuk memformulasi fenomena
perpindahan.
Proses dibagi menjadi sejumlah elemen
proses yang dijelaskan dengan hukum
kekekalan massa, momentum, dan energi.
MODEL BERDASARKAN
PRINSIP FISIKOKIMIA
Model deterministik atau elemen model:
Nilai atau seperangkat nilai setiap variabel atau parameter
model pada kondisi tertentu telah ditentukan.
Model statistik atau elemen model statistik
Variabel dan parameter model merupakan besaran
statistik, berupa probabilitas atau momen dari fungsi
densitas probabilitas.
Misalnya
Jika fungsi densitas probabilitas P(Y ) berlaku untuk
variabel statistik Y, maka P(Y) dY adalah probabilitas
variabel tersebut yang berada dalam rentang dY di sekitar
Y.
MODEL BERDASARKAN
PRINSIP FISIKOKIMIA
Klasifikasi berdasarkan jenis persamaan
Tingkat kesulitan metode penyelesaian
berkurang dari kanan ke kiri.
MODEL PDF
Model berbasis persamaan transport dalam bentuk fungsional
P(1, . . . , n).
Probabilitas menemukan variabel terikat (1, . . . , n) dalam
rentang d1, . . . , dn di sekitar fungsi 1(x, t), . . ., n(x, t) adalah
P(1, . . . , n)d1, . . . , dn.
Memberi informasi statistik proses statistik.
Memberi fungsi distribusi variabel proses.
Contoh:
mekanika statistik, teori kinetik gas, campuran makro dalam
distribusi waktu tinggal, distribusi ukuran kristal, distribusi
aktivitas pada pelet katalis, dan distribusi umur dan ukuran
biakan mikrobiologi.
MODEL EMPIRIS
Korelasi respons proses terhadap perubahan satu
atau beberapa variabel proses.
Contoh:
Fitting polinomial pada data eksperimen, respons proses
pada pengendalian proses dalam bentuk fungsi transfer
pada domain waktu atau frekuensi.
Merupakan model statistik karena data diperoleh
secara eksperimen dan berisi kesalahan statistik.
Memiliki makna terbatas dalam menjelaskan proses
atau elemen proses;
Misal: prediksi berada di luar rentang percobaan.
MODEL BERDASARKAN
PRINSIP FISIKOKIMIA
Dr.rer.nat. Ir. Yuswan Muharam, M.T.
Departemen Teknik Kimia
Fakultas Teknik
Universitas Indonesia
ILUSTRASI PROSES PEMODELAN
Proses pendinginan fluida yang mengalir di
dalam pipa berpenampang lingkaran.
Dimulai dengan model yang paling sederhana.
Menambah tingkat kesulitan untuk meningkatkan
keakuratan.
MODEL 1 –
ALIRAN SUMBAT
Buat sketsa sistem.
Plug flow:
Profil kecepatan fluida
berbentuk plug (merata
pada posisi radial).
Elemen fluida
bercampur sempurna
ke arah radial sehingga
temperatur fluida
merata pada bidang
normal terhadap bidang
aliran (arah radial).
MODEL 1 –
ALIRAN SUMBAT
Jika tube tidak panjang atau perbedaan temperatur tidak besar,
maka sifat fisik fluida tidak banyak berubah.
Asumsi:
1. Keadaan tunak;
2. Sifat fisik fluida (, Cp, k dll) konstan;
3. Temperatur dinding konstan dan merata (tidak berubah ke arah z
atau r) dengan nilai Tw;
4. Temperatur inlet konstan dan merata (tidak berubah ke arah r)
dengan nilai T0, dimana T0 > Tw;
5. Profil kecepatan berbentuk plug atau datar sehingga merata ke
arah z atau r;
6. Fluida bercampur sempurna (turbulen Re > 2100) sehingga
temperatur merata ke arah radial;
7. Konduksi termal sepanjang sumbu relatif kecil dibandingkan
konveksi.
MODEL 1 –
ALIRAN SUMBAT
Buat sketsa elemen volume diferensial
sistem (fluida alir) atau “volume kontrol."
MODEL 1 –
ALIRAN SUMBAT
Kembangkan hukum kekekalan energi umum
Keadaan tunak akumulasi nol.
Tidak ada sumber kimia, nuklir atau listrik tidak ada
pembangkit panas.
Panas hanya berpindah melalui perimeter elemen akibat
perbedaan temperatur antara fluida dan dinding.
Laju pengambilan panas menggunakan hukum
pendinginan Newton (+)
MODEL 1 –
ALIRAN SUMBAT
Kembangkan hukum kekekalan energi umum
Luas kontak = keliling x panjang.
Koefisien perpindahan panas, h konstan.
Bar di atas T menyatakan nilai rata-rata antara
T(z) dan T (z + z)
MODEL 1 –
ALIRAN SUMBAT
Kembangkan hukum kekekalan energi umum
Sepanjang sumbu, panas masuk dan keluar
hanya melalui konveksi (aliran) sehingga
Dua suku pertama: laju alir massa x entalpi lokal
(temp. rujukan untuk entalpi = 0).
MODEL 1 –
ALIRAN SUMBAT
Disusun kembali dan dibagi z, diperoleh
Dengan
menjadi
MODEL 1 –
ALIRAN SUMBAT
Pengelompokan parameter menjadi satu suku
(parameter lumping)
menjadi
.
dimana
.
MODEL 1 –
ALIRAN SUMBAT
Persamaan diferensial biasa orde pertama.
MODEL 2 –
KECEPATAN PARABOLIK
Jika aliran lebih lambat (Re < 2100), kecepatan berbentuk
parabola.
v0 = kecepatan rata-rata
vz = kecepatan lokal (variatif).
Asumsi 5, 6, dan 7 dimodifikasi:
5. Profil kecepatan arah z berbentuk parabola dan tergantung pada
posisi r.
6. Fluida tidak tercampur sempurna ke arah radial sehingga
konduksi panas radial diperhitungkan.
7. Karena konveksi lebih kecil, konduksi panas aksial
dipertimbangkan.
MODEL 2 –
KECEPATAN PARABOLIK
Volume kontrol berbentuk cincin dengan tebal r dan panjang
z;
Panas melewati dua permukaan, area anular yang normal
terhadap aliran fluida, dan area sepanjang keliling cincin;
Fluks panas (laju per satuan luas normal) menggunakan
konduksi molekular.
MODEL 2 –
KECEPATAN PARABOLIK
Laju bersih pembentukan (pelepasan) panas
oleh konduksi = fluks x luas area normal
terhadap arah fluks.
Hukum kekekalan panas elemen volume
MODEL 2 –
KECEPATAN PARABOLIK
Dua koordinat posisi proses diferensiasi parsial,
misalnya
disusun kembali dan dibagi dengan 2zr
.
.
MODEL 2 –
KECEPATAN PARABOLIK
Dengan limit, diperoleh
Turunan terhadap z menunjukkan nilai r konstan,
sehingga r dapat ditempatkan di luar suku; dengan
membagi dengan r dan menata kembali, diperoleh
MODEL 2 –
KECEPATAN PARABOLIK
Substitusi hukum Fourier dan uz
ke
diperoleh
Persamaan diferensial parsial orde dua
GABUNGAN LAJU DAN
KESETIMBANGAN
Contoh: adsorpsi menggunakan unggun padat granular.
Adsorpsi lebih cepat dibandingkan difusi internal, sehingga pada
dan dekat partikel terjadi kesetimbangan lokal.
q = komposisi rata-rata fasa padat (mol solut teradsorpsi per
satuan volume partikel),
C* = komposisi solut (mol solut per satuan volume fluida), yang
setimbang.
Asumsi:
Pengontrol laju: laju perpindahan antara fasa mengalir dan fasa
diam (padat).
GABUNGAN LAJU DAN
KESETIMBANGAN
Konsep aliran sumbat profil
kecepatan fluida datar.
Adsorbat di dalam fluida encer
efek panas diabaikan (isotermal).
Partikel sangat kecil efek difusi
aksial diabaikan transportasi
fasa fluida disebabkan aliran
konveksi.
GABUNGAN LAJU DAN
KESETIMBANGAN
Transportasi antarfasa mengikuti hukum laju yang
berangkat dari keadaan kesetimbangan
termodinamika.
Luas antarfasa total tidak diketahui koefisien
perpindahan volumetrik (kca); a = luas antarfasa
total per satuan volume kolom paking.
Persamaan laju inkremental
.
.
.
GABUNGAN LAJU DAN
KESETIMBANGAN
Neraca solut di kedua fasa
Vo: kecepatan superfisial fluida (terjadi jika tube kosong);
: fraksi volume kosong di antara partikel (volume kosong interstitial)
(1 - ): fraksi volume fasa padat;
Laju akumulasi: fasa fluida (C) dan fasa padat (q).
Pembagian dengan Az dan limit menghasilkan
GABUNGAN LAJU DAN
KESETIMBANGAN
Neraca solut di fasa diam saja
Tidak ada reaksi kimia;
Laju akumulasi sama dengan laju perpindahan ke padatan
Dibagi dengan A z
Jika kesetimbangann dicapai C C*
.
GABUNGAN LAJU DAN
KESETIMBANGAN
Substitusi
menghasilkan
Kondisi batas
.
.
.
PROSEDUR PEMODELAN
1. Gambar sketsa sistem dan definisikan besaran kimia, fisika dan
geometri.
2. Pilih variabel terikat (respons).
3. Pilih variabel bebas (misal z, t).
4. Buat daftar parameter (konstanta fisik, ukuran dan bentuk); buat
pula daftar parameter tak konstan (misal viskositas yang berubah
terhadap temperatur).
5. Gambar sketsa perilaku variabel terikat, seperti profil temperatur
yang diharapkan.
6. Buat “volume kontrol" untuk elemen diferensial atau berhingga
sistem (misal CSTR); buat sketsa elemen dan indikasikan semua
lintasan masuk dan keluarnya.
PERSAMAAN
DIFERENSIAL BIASA -
PROBLEM NILAI AWAL
Dr.rer.nat. Ir. Yuswan Muharam, M.T.
PERSAMAAN
DIFERENSIAL BIASA (PDB)
Persamaan diferensial untuk fungsi yang
hanya tergantung pada satu variabel
Ruang (x, y, z, r)
Waktu (t).
Solusi PDB:
Kondisi awal (problem nilai awal);
Kondisi batas (problem nilai batas).
PERSAMAAN
DIFERENSIAL BIASA (PDB)
Problem nilai awal:
jika semua kondisi berada pada satu titik dan
dapat diintegrasi mulai dari titik tersebut.
Problem nilai batas dua titik:
jika pada satu titik terdapat satu atau lebih kondisi
dan pada titik lain terdapat satu atau lebih kondisi
yang lain.
Contoh problem PDB:
kontrol parameter, kinetika di dalam reaktor batch,
reaktor alir sumbat.
KLASIFIKASI PDB
Dasar klasifikasi:
Orde,
Kelinearan,
Kondisi batas.
KLASIFIKASI BERDASARKAN
ORDE
Orde persamaan diferensial = orde tertinggi
dari derivat (turunan).
Orde pertama:
dy
y kx
dx
Orde kedua: d2y dy
2
y kx
dx dx
Orde ketiga: 2
3
d y d y dy
2
a 2 b kx
dx
3
dx dx
KLASIFIKASI BERDASARKAN
KELINEARAN
Linear: tidak mengandung perkalian variabel
terikat, derivatnya atau keduanya.
Tak linear: mengandung perkalian variabel
terikat atau derivatnya atau keduanya.
Linear: dy
y kx
dx
Tak Linear: d2y dy
2
y kx
dx dx
2
3
d y d y dy
2
a 2 b kx
dx
3
dx dx
KLASIFIKASI BERDASARKAN
KONDISI BATAS
Problem nilai awal:
Semua nilai variabel terikat dan/atau turunanya
diketahui pada nilai awal variable bebas.
Problem nilai batas:
Variabel terikat dan/atau turunannya diketahui
pada lebih dari satu variabel bebas.
KLASIFIKASI BERDASARKAN
KONDISI BATAS
PDB orde ke-n:
dny d n 1 y
b0 x n b1 x n1 ... bn 1 x bn x y Rx
dy
dx dx dx
R(x) = 0 homogen.
R(x) 0 tak homogen.
Koefisien {bi | i = 1, 2, …, n}
koefisien variabel jika fungsi dari x;
koefisien konstan jika skalar.
KLASIFIKASI BERDASARKAN
KONDISI BATAS
Untuk mendapatkan solusi sebuah PDB orde
ke-n atau sejumlah n PDB orde pertama,
diperlukan spesifikasi n nilai variabel terikat
(turunannya) pada nilai-nilai tertentu variabel
bebasnya.
SOLUSI PDB -
PROBLEM NILAI
AWAL
Dr.rer.nat. Ir. Yuswan Muharam, M.T.
KUADRATUR
Hanya satu PDB (linear atau tidak linear)
f y
dy
dt
y 0 y0
Pemisahan variabel:
dy
dt
f y
y t
dy
y0
dt
f y 0
● Jika dapat diselesaikan secara analitik solusi eksak.
KUADRATUR
Misal: problem kinetika untuk reaksi orde dua
dc
kc 2
dt
c0 c0
Pemisahan variabel dan integrasi:
dc
2
kdt
c
1
kt D
c
Kondisi batas menghasilkan:
1 1
kt
c c0
METODE EKSPLISIT
Jika nilai y pada tn diketahui, maka
perhitungan vektor y pada waktu berikutnya
tn +1 hanya memerlukan nilai vektor y yang
diketahui tersebut serta turunannya dy/dt =
f(y) pada waktu tn (dan waktu sebelumnya).
METODE EKSPLISIT
Integrasi numeris PDB dapat dilakukan jika sistem
terdiri dari n PDB orde pertama simultan dalam
bentuk:
f1 y1 , y2 ,..., yn , x
dy1
dx
f 2 y1 , y2 ,..., yn , x
dy2
dx
. Bentuk kanonis
.
.
f n y1 , y2 ,..., yn , x
dyn
dx
METODE EKSPLISIT
Apabila kondisi awal diberikan pada titik x0:
y1 x0 y1,0
y2 x0 y2,0
f1 y1 , y2 ,..., yn , x
. dy1
.
dx y1 F1 x
f 2 y1 , y2 ,..., yn , x y2 F2 x
dy2
.
dx
yn x0 yn ,0 .
.
.
. Solusinya:
. .
yn Fn x
f n y1 , y2 ,..., yn , x
dyn
dx
METODE EKSPLISIT
Dalam bentuk matriks
f1 y1 , y2 ,..., yn , x
dy1
dx
dy2
f 2 y1 , y2 ,..., yn , x y1 x0 y1,0 y1 F1 x
y2 F2 x
dx
y2 x0 y2,0
.
. .
.
. . .
. .
f n y1 , y2 ,..., yn , x
dyn
dx yn x0 yn ,0 yn Fn x
f x, y
dy
dx y x0 y 0 y F x
METODE EKSPLISIT
Persamaan diferensial orde tinggi
d nz dz d 2 z d n 1 z
n
G z , , 2 ,..., n 1 , x
dx dx
dx dx
dapat diubah menjadi seperangkat
persamaan orde satu.
Caranya?
METODE EKSPLISIT
d nz dz d 2 z d n 1 z
G z , , 2 ,..., n 1 , x
dx dx
dx n
dx z y1
dz dy1
y2
dx dx
d 2 z dy2
2
y3
dx dx
.
.
Transformasi
.
d n 1 z dyn 1
n 1
yn
dx dx
d n z dyn
n
dx dx
METODE EKSPLISIT
d nz dz d 2 z d n 1 z
G z , , 2 ,..., n 1 , x
dx dx
z y1 dx n
dx
dz dy1
y2 dy1
dx dx y2
d 2 z dy2 dx
y3 dy2
y3
2
dx dx
. dx
.
. substitusi
.
.
n 1
.
d z dyn 1
yn
G y1 , y2 , y3 ,..., yn , x
n 1 dyn
dx dx
dx
d n z dyn
n
n persamaan orde pertama
dx dx bentuk kanonis
METODE EKSPLISIT
dy1
y2
dx
dy2
y3
dx
.
.
.
G y1 , y2 , y3 ,..., yn , x
dyn
dx
Jika sisi kanan PDB bukan fungsi variabel bebas, maka disebut
persamaan otonom.
f y
dy
dx
Jika f(y) linear terhadap y, maka dapat ditulis: y’ = Ay
METODE EKSPLISIT
Ubah persamaan berikut ke bentuk
kanonisnya!
d 4z d 3z d 2z dz
4
5 3 2 2 6 3z 0
dt dt dt dt
d 4z d 3z d 2z dz
5 3 2 2 6 3 z e t
dt 4 dt dt dt
3
d 3z d 2 z dz
z 2 2 2z 0
dx 3 dx dx
d 3z d 2z 2 dz
3
z 3
2
z 5z 0
dx dx dx
METODE EKSPLISIT
Metode Euler
Metode Adam-Bashford
Runge-Kutta
METODE EULER
f t , y
dy
Bentuk kanonis:
dt
Diferensial:
Nilai rata-rata f pada h adalah f(y(tn)).
METODE EULER - CONTOH
Contoh:
Dekomposisi nitrogen dioksida di dalam reaktor alir sumbat
dengan laju reaksi
Konstanta laju reaksi pada 383°C = 5030 ml/mol/detik.
Asumsi:
Difusi aksial sangat kecil sehingga diabaikan,
Profil kecepatan berbentuk plug.
Hitung profil konsentrasi keadaan tunak pada temperatur
konstan!
METODE EULER - JAWAB
Neraca massa
u = kecepatan, S = luas penampang lintang reaktor.
METODE EULER - JAWAB
Bagi dengan z dan susutkan elemen menjadi nol
(limit)
Kondisi awal:
Solusi analitik:
METODE EULER - JAWAB
Kalikan sisi kiri dengan S/S
Jadikan persamaan tak-berdimensi
METODE EULER - JAWAB
Metode Euler:
Jika h = 0,2
METODE EULER - JAWAB
METODE EULER - LATIHAN
Selesaikan PDB di bawah dengan
menggunakan metode Euler!
dy
t y
dt
y 0 1
METODE EULER - JAWABAN
tn yn f(yn) t f(yn)
METODE ADAM-BASHFORD
Orde kedua:
Orde keempat:
METODE ADAM-BASHFORF -
LATIHAN
Selesaikan PDB di bawah dengan
menggunakan metode Adam-Bashford orde-
keempat!
dy
t y
dt
y 0 1
METODE EKSPLISIT
Metode eksplisit orde tinggi perlu solusi (sisi kanan) yang
dievaluasi pada waktu-waktu sebelumnya.
Evaluasi mudah dilakukan kecuali pada permulaan evaluasi
gunakan metode Euler dengan ukuran tahap yang sangat kecil
selama beberapa tahap untuk mendapatkan nilai-nilai
permulaan.
Keuntungan metode Adam – Bashford orde keempat:
Hanya menggunakan satu evaluasi fungsi per tahap,
Akurasi orde tinggi.
Kelemahan metode Adam – Bashford orde keempat:
Perlu metode lain untuk memulai.
METODE RUNGE-KUTTA
Skema titik tengah:
titik tengah digunakan untuk menghitung titik tak
diketahui pada tn + 1;
Argument yn + (h/2)fn = slope pada tn + (h/2), titik
tengah antara tn dan tn + 1.
METODE RUNGE-KUTTA
Skema korektor predictor-trapezoid Euler.
METODE RUNGE-KUTTA-GILL
Orde ke-empat;
Paling banyak digunakan karena
memerlukan sedikit memori komputer;
Ditulis dalam bentuk vektor untuk sistem
PDB;
METODE RUNGE-KUTTA-GILL
METODE RUNGE-KUTTA-
FELDBERG
Orde ke-enam
Nilai yn+1 – zn+1 merupakan taksiran error untuk yn+1
LATIHAN METODE
RUNGE-KUTTA
Selesaikan dengan menggunakan metode
Runge-Kutta-Gill!
y 0 1
dy
8. y x2 ,
dx
Gunakan h = 0,01!