## METODE OPTIMASI PADA SISTEM PENGENDALIAN PROSES TANGKI PEMANAS BERPENGADUK
Rudy Agustriyanto
Jurusan Teknik Kimia, Fakultas Teknik, Universitas Surabaya, Jl. RayaKalirungkut, Surabaya, 60293, Indonesia Email: rudy.agustriyanto@staff.ubaya.ac.id

## Abstrak
Pada penelitian ini, penyetelan pengendali PID (Proportional Integral  Derivative)  dilakukan  dengan  metode optimasi dan hasil akhirnya dibandingkan dengan metode penyetelan lain (yaitu metode sintesa langsung dan Ziegler-Nichols). Tujuan penelitian adalah mendapatkan response output terbaik ketika terjadi perubahan set point dan terjadi gangguan. Penyetelan pengendali adalah proses yang berguna dalam menetapkan parameter pengendali  yang  dapat  mempengaruhi  response  output  dari  perubahan  setpoint  (masalah  servo)  maupun ketika  terjadi  gangguan  (masalah  regulatori).  Kinerja  pengendali  dapat  dievaluasi  dari  harga  SSE  (Sum Squares of Error) ketika sistem dikendalikan oleh pengendali PID yang disetel dengan metode tertentu. Hasil penelitian menunjukkan bahwa metode optimasi tampak lebih baik dari metode  lainnya  karena  memberikan kinerja terbaik.
Kata kunci: Optimasi, Penyetelan, Pengendali PID

## Abstract
In this research, PID (Proportional Integral Derivative) controller tuning was done using optimization method and the result  was  compared  to  the  other  well-known  tuning  methods  (Direct  synthesis,  Ziegler–  Nichols).  The  purpose  of this research is to obtain the best output response when the set point and disturbance are changed. Tuning is a useful process to set the controller parameter that can affect the output response from the changing set point (servo problem) as well as from the disturbance (regulatory problem). The controller performance can be evaluated from the value of Sum  Squares  of  Error  when  the  system  is  controlled  by  PID  controllers  that  is  tuned  with  chosen  method.  The research  result  showed  that  optimization  method  was  better  than  the  other  methods  because  it  could  give  the  best result.
Keywords :  Optimization, Tuning, PID control

## Pendahuluan
Optimasi bisa digunakan dalam berbagai bidang yaitu  baik  dalam  bidang  ilmu  pengetahuan,  teknik, dan  bisnis.  Dalam  bidang  statistik  berbagai  prinsip optimasi  akan  menghasilkan  keuntungan  maksimum, dan  kerugian  minimum,  sementara  di  bidang  bisnis metode  optimasi  digunakan  untuk  memaksimalkan keuntungan  dan  meminimalkan  biaya  dengan menggunakan  sumber  daya  semaksimum  mungkin menggunakan usaha seminimal mungkin dalam upaya untuk  meningkatkan  keuntungan  bisnis.  Dalam bidang  teknik  metode  ini  dapat  digunakan  untuk menafsirkan beberapa persamaan.
Penggunaan  metode  optimasi  ini  mulai  banyak di  gunakan  pada  pengendali  PID.  Pada  optimasi, pengendali  PID  termasuk  pengendali  konvensional atau klasik yang masih populer di industri hingga saat ini  karena  efektifitasnya  luar  biasa,  implementasinya mudah  dan  aplikasinya  sangat  luas.  Seiring  dengan perkembangan teknologi, mulai dikembangkan pengendali modern, seperti nonlinear control, optimal control,  variable  structure  control,  adaptive  control. Masalah utama pengendali PID adalah penyetelan,yaitu  menentukan  nilai  parameter  kc,    τi
dan  τd agar diperoleh performa sistem yang optimal [6].
Penelitian  ini  digunakan    untuk  mendapatkan nilai  penyetelan  pengendali  proses  orde  satu  dari sistem  tangki  pemanas  berpengaduk  yang  terbaik menggunakan  metode  optimasi  jika  dibandingkan dengan metode konvensional lain.
Penelitian yang dilakukan Tavakoli menggunakan  analisis  dimensi  dan  teknik  optimasi numerik  membuktikan  bahwa  metode  yang  optimal untuk  penyetelan  pengendali  PID  untuk  proses  orde satu  dengan  time  delay  ditampilkan  dengan mempertimbangkan nilai Integral Square Error (ISE), Integral  Absolute  Error  (IAE)  dan  Integral  Time Absolute  Error  (ITAE).  Hasil  simulasi  menunjukkan bahwa  metode  yang  diusulkan  memiliki  keunggulan besar  jika  dibandingkan  dengan  teknik  konvensional [6].
Metode penyetelan dapat menawarkan perkiraan awal  nilai  parameter  pengendali  PID.  Namun  tujuan perancangan  yang  digunakan  untuk  pengembangan dari metode penyetelan dapat menjadi sangat berbeda dari  tujuan  kinerja  spesifik  untuk  proses  yang diselidiki.  Akibatnya,  insinyur  pengendali  sering menyesuaikan  kebutuhan  untuk  merancang  suatu sistem  lintas  tertutup  yang  memuaskan.  Diusulkan
sebuah  reformulasi  dari  rumusan  matematika,  yang mengubah  rumusan  awalnya  non  –  linear  untuk menjadi linear yang dapat menyelesaikan permasalahan optimasi secara global [5].
Algoritma  auto  –  tuning  untuk  pengendali  PI dan  PID  berdasarkan  pengulangan  percobaan  yang pernah  dilakukan  diusulkan  untuk  meminimalkan kesalahan  gangguan integral  error (IE)  dengan memaksimalkan integral gain, keuntungan subyek  ke fasa yang diinginkan dan meminimalkan batasan gain yang  diperlukan.  Keuntungan  utama  dari  algoritma auto  –  tuning  untuk  penyetelan  pengendali  PID dengan sistem lintas tertutup mirip dengan pengendali  PID  yang  dirancang  menggunakan optimasi  numerik  off–line.  Selain  itu  algoritma  ini berlaku  untuk  setiap  struktur  model  linier,  termasuk dead time dan sistem fase non–minimum [3].

## Teori
Optimasi merupakan salah satu metode kuantitatif  utama  sebagai  sarana  pengambilan keputusan  [2].  Optimasi  berfokus  pada  pemilihan sesuatu  yang  paling  baik  dari  keseluruhan  rangkaian dengan  menggunakan  suatu  metode  kuantitatif. Berbagai  macam  masalah  pada  desain,  konstruksi, operasi,  dan  analisa  pada  pabrik  kimia  dapat diselesaikan  dengan  optimasi.  Masalah  yang  terjadi dalam  bidang  teknik  dapat  direpresentasikan  oleh suatu persamaan atau mungkin sebuah data dari hasil eksperimen.  Optimasi  digunakan  dengan  tujuan menentukan  nilai  variabel  dalam  proses  yang  dapat menghasilkan kriteria performa kerja terbaik.
Dalam pabrik kimia, keuntungan dapat diperbesar  dengan  menaikkan  performa  kerja  pabrik tersebut,  misalnya  dengan  meningkatkan  nilai  yield dari  produk  yang  diinginkan  dan  berdaya  jual, mengurangi  konsumsi  energi  yang  digunakan, mempercepat  laju  pemrosesan,  dan  pengupayaan waktu  yang  lebih  lama  antara  tiap  proses  shutdown. Optimasi berfungsi untuk mengurangi biaya perawatan  alat,  mengurangi  jumlah  peralatan  yang digunakan,  dan  dapat  mereorganisasi  penggunaan sumber daya manusia agar lebih efisien.
Cakupan  penerapan  optimasi  dalam  suatu masalah  dapat  meliputi  keseluruhan  perusahaan, pabrik,  proses  ataupun  unit  operasi  tunggal,  satu bagian  dari  suatu  peralatan,  dan  sebagainya.  Dalam proses  optimasi  ini  dikenal  pula  istilah  program simulator  desain  proses  yang  merupakan  suatu  alat bantu  untuk  pemecahan  masalah  dalam  industri. Program  ini  dapat  melaksanakan  perhitungan  neraca massa  dan  energi  untuk  tiap-tiap  bagian  dari  suatu peralatan  dan  kemudian  mengkombinasikannya menjadi suatu unit produksi total. Penggunaan fungsi iterasi  untuk  simulasi  semacam  ini    biasanya  sangat diperlukan  terutama  untuk  mencapai  diagram  alir proses yang diinginkan.
Ada  2  data  penting  yang  diperlukan  dalam mengindetifikasi  suatu  proses  optimasi  yaitu:  data
tentang  keuntungan  dan  harga  dari  suatu  pabrik  atau masing-masing  unit  yang  berada  di  dalamnya,  dan data tentang catatan periodik dari sistem operasi yang bekerja dalam suatu pabrik tersebut. Catatan tersebut berisikan  informasi  mengenai  neraca  massa  dan energi sistem, efisensi, kualitas dan kuantitas produk yang  dihasilkan,  dan  jumlah  umpan  masuk  yang tersedia.
Tiap  masalah  optimasi  memiliki  3  komponen penting yaitu:
1. Setidaknya  terdapat  satu  fungsi  yang  akan dioptimalkan  (misalnya  fungsi  keuntungan,  fungsi biaya, dan sebagainya)
2. Batasan persamaan (persamaan)
3. Batasan pertidaksamaan (pertidaksamaan)
Langkah-langkah  untuk  memecahkan  masalah optimasi:
1.  Menganalisa  proses  itu  sendiri  sehingga  variabel proses  dan  karakter-karakter  spesifik  lainnya dapat terdefinisi.
2.  Menentukan  kriteria  optimasi  dan  membuat spesifik  fungsi  yang  diinginkan  untuk  dinyatakan dalam  variabel-variabel  dan  konstanta-konstanta yang  telah  didefinisikan  di  atas.  Langkah  ini menyediakan suatu performance model.
3.  Mengembangkan  model  matematis  dari  proses atau  peralatan  tersebut  yang  menunjukkan hubungan antara variabel input-output dari proses dan  juga  koefisien-koefisien  yang  berhubungan dengan variabel-variabel tersebut. Selain itu, perlu disertakan  juga  batasan  dalam  persamaan  dan pertidaksamaan  yang  digunakan.  Dapat  dipakai juga  konsep-konsep  neraca  massa  dan  energi, persamaan  empiris,  konsep  implisit  dan  batasan- batasan eksternal. Kemudian melakukan identifikasi  variabel  dependent  dan  variabel independent  untuk  memperoleh  jumlah  degree  of freedom.
4.  Jika  jangkauan    masalah  terlalu  luas  maka  dapat dilakukan  pemecahan  problem  tersebut  menjadi bagian-bagian  yang  lebih  kecil  atau  membuat fungsi  obyektif  menjadi  model  yang  lebih sederhana.
5.  Menerapkan  teknik  optimasi  yang  sesuai  dengan model matematika yang ada.
6.  Memeriksa  kembali  jawaban  yang  diperoleh,  dan menyelidiki  pengaruh  hasil  terhadap  koefisien- koefisien  dalam  masalah  serta  bentuk-bentuk asumsi yang mungkin digunakan.
Masalah  optimasi  dapat  menemui  kendala  atau mungkin kegagalan yang disebabkan oleh:
1. Fungsi  obyektif  dan  atau  fungsi  pembatas kemungkinan  mempunyai  daerah  diskontinu dalam parameter kontinu yang dipakai.
2. Fungsi  obyektif  dan  atau  fungsi  pembatas kemungkinan bukan suatu fungsi linear terhadap variabel-variabelnya.
3. Fungsi  obyektif  dan  atau  fungsi  pembatas kemungkinan didefinisikan dalam bentuk interaksi antar variabelnya yang terlalu rumit.
4. Nilai  dari  Fungsi  obyektif  dan  atau  fungsi pembatas  termasuk  tidak  sensitif,  atau  justru terlalu sensitif terhadap perubahan dari variabel- variabel yang ada.
5. Fungsi  obyektif  menunjukkan  suatu  kondisi optimum  lokal  saat  kondisi  optimum  global dicapai.

## Metodologi Penelitian
Matlab  dan  Simulink  digunakan  untuk melakukan  simulasi  proses  dan  implementasi  sistem pengendaliannya. Adapun sistem yang ditinjau adalah sistem tangki pemanas berpengaduk seperti ditunjukkan pada Gambar 1.
Gambar 1. Sistem tangki pemanas berpengaduk
Pada sistem tangki  pemanas  di  atas  air  masuk pada  suhu  To  (K)  dengan  laju  alir  volumetric  Fo (m3/menit).  Air  yang  masuk  dalam  tangki  tersebut diinginkan  memiliki  suhu  T  (K)  saat  keluar  tangki. Hal  ini  dapat  dilakukan  dengan  suplai  panas  yang berasal  dari  steam  yang  dialirkan  dalam  koil.  Pada sistem  dipasang thermocouple sebagai  sensor  untuk mengetahui  suhu  air  di  dalam  tangki,  untuk selanjutnya thermocouple akan  membandingkan apakah  suhu  dalam  tangki  sudah  sama  dengan  suhu set  point.  Selisih  dari  pembacaan thermocouple dan suhu set-point akan  mempengaruhi  kerja  pengendali yang akan mengatur bukaan valve supaya steam yang disuplai sesuai untuk mencapai suhu keluar (T) yang diinginkan.
Fungsi hantar untuk sistem ini telah diturunkan dengan  menggunakan  neraca  massa  dan  panas  [1], dan diperoleh hasil sebagai berikut [7]:
 Antara output T dengan gangguan To adalah:
         1 1 ) ( ) ( s Fo V s To s T
(1)
 Antara output T dengan Q adalah:
           1 1 ) ( ) ( s Fo V Fo cp s Q s T  (2)
Dimana :
 Fungsi transfer sensor (thermocouple) = 1
 Fungsi transfer pengendali valve = 1
 ρ = Massa jenis (kg/m 3 )
 A = Luas permukaan tangki (m 2 )
 h = Ketinggian tangki (m)
 t = Waktu (menit)
 Fo = Laju alir masuk (m 3 /menit)
 F = Laju alir keluar (m 3 /menit)
 T = Suhu keluar (K)
 To = Suhu masuk (K)
 Q = Kalor yang disuplai (kJ/menit)
 cp = Heat capacity (kJ/kg. K
Fungsi hantar tersebut kemudian disimulasikan dengan  Simulink  dan  dibuat  lintas  tertutup  dengan pengendali  PID.  Diagram  blok  untuk  sistem  lintas tertutup  dapat  dilihat  pada  Gambar  2.  Sedangkan skema  penyetelan  pengendali  PID  dengan  metode optimasi ditunjukkan pada Gambar 3.
Dari Gambar 2, dapat ditentukan rumus untuk gain  proses  (Kp)  dan  time  constant  proses  (  ) adalah:
Fo cp K p     1 (3)
Fo V   (4)
Sedangkan waktu tunda (θ) yang terjadi pada proses di atas , diasumsikan sebesar 1 menit.
Dalam penelitian ini dilakukan penyetelan pengendali  untuk  beberapa  variabel  percobaan seperti ditunjukkan pada Tabel 1.
Kp/( τs+1)
+ - PID 1 1 S Tsp E V Q To T 1/( τs+1) + +
Gambar 2. Diagram blok sistem lintas tertutup
Menentukan Objective Function
Menentukan Harga Trial Awal Parameter Pengendali
Mensimulasi pada Software MATLAB Menggunakan Metode Least Square Non Linear

## Menetapkan Harga Parameter Pengendali
Gambar 3. Skema penyetelan dengan metode optimasi
Tabel 1. Data Variabel Proses Simulasi No Variabel proses Harga Parameter Model 1 V = 25 m 3 Kp = 4,804.10 -5
       kJ menit K τ  = 5 menit θ = 1 menit Fo = 5 m 3 / menit To = 303,15 K
ρ = 995,68 kg/m 3 cp = K kg kJ  181 , 4 2
V = 35 m 3
Kp = 4,036.10 -5
       kJ menit K τ  = 5,8333 menit θ = 1 menit Fo = 6 m 3 / menit To = 303,15 K
ρ = 995,68 kg/m 3 cp = K kg kJ  181 , 4 3
V = 40 m 3
Kp = 4,804.10 -5
       kJ menit K τ  = 8 (menit) θ = 1 menit Fo = 5 m 3 / menit To = 303,15 K
ρ = 995,68 kg/m 3 cp = K kg kJ  181 , 4 4
V = 40 m 3 Kp = 3,002.10 -5        kJ menit K τ  = 5 (menit) θ = 1 menit Fo = 8 m 3 / menit To = 303,15 K ρ = 995,68 kg/m 3 cp = K kg kJ  181 , 4 5
V = 50 m 3 Kp = 3,002.10 -5
       kJ menit K τ  = 6,25 (menit) θ = 1 menit Fo = 8 m 3 / menit To = 303,15 K ρ = 995,68 kg/m 3 cp = K kg kJ  181 , 4

## Hasil
Nilai  parameter  pengendali  PID  (kc,  τi,  τd) untuk  sistem  yang  disetel  dengan  metode  optimasi ditunjukkan  pada  Tabel  2.  Gambar  4  s.d  8 menunjukkan  profil  output  T  dan  kalor  Q  pada perubahan set point maupun perubahan gangguan To untuk kelima sistem yang ditinjau.
Pengujian pada Sistem  no. 1 (Simulasi no  1) dilakukan  dengan  perubahan  set  point  berupa  fungsi step  dari  suhu  345,15  K  menjadi  362,15  K  pada waktu  t  =  1  menit  dilanjutkan  berubah  menjadi 347,15  pada  waktu  t  =  75  menit  dan  berubah  lagi menjadi  357,15  K  pada  t  =  120  menit  serta  adanya
gangguan  sebagai  perubahan  suhu  air  masuk  (To) berupa  fungsi  step  dari  suhu  303,15  K  menjadi 318,15  K  pada  t  =  30  menit  lalu  berubah  menjadi 293,15  K  pada  t  =  95  menit.  Profil  sistem  lintas tertutup dapat dilihat pada Gambar 4.
Kemudian  pada  Sistem  2  (Gambar  5)  juga dilakukan  pengujian  pada  sistem  dengan  perubahan set  point  berupa  fungsi  step  dengan  menaikkan  suhu dari suhu 338,15 K menjadi 355,15 K pada waktu t = 1  menit  dilanjutkan  menurunkan  suhu    menjadi 340,15  pada  waktu  t  =  75  menit  dan  mengubahnya lagi  menjadi  350,15  K  pada  t  =  120  menit,  serta adanya  gangguan  sebagai  perubahan  suhu  air  masuk
(To) berupa fungsi  step  dari  suhu  303,15  K  berubah menjadi 318,15 K pada t = 30 menit  dan dilanjutkan berubah menjadi 293,15 K pada t = 95 menit.
Tabel  2.  Parameter  Pengendali  yang  Disetel  dengan Metode Optimasi
Simulasi No kc τ i Τ d 1 199910 6,9668 0,6518 2 173480 8,0601 0,6572 3 195720 11,0548 0,6654 4 124940 6,9668 0,6518 5 247020 8,6177 0,6593
0 20 40 60 80 100 120 140 160 180 200 340 350 360 370 Te m p e ra tu re T (K ) 0 20 40 60 80 100 120 140 160 180 200 290 300 310 320 D is tr u b a n c e To ( K ) 0 20 40 60 80 100 120 140 160 180 200 0 5 x 10 6 Waktu (menit) K a lo r Q ( k J /m e n it ) Set Point Output
Gambar 4. Profil  output  dan  kalor  hasil  perubahan  set point dan gangguan dengan metode optimasi – Sistem 1 .
0 20 40 60 80 100 120 140 160 180 200 330 340 350 360 Te m p e ra tu re T (K ) 0 20 40 60 80 100 120 140 160 180 200 290 300 310 320 D is tr u b a n c e To ( K ) 0 20 40 60 80 100 120 140 160 180 200 0 5 10 x 10 6 Waktu (menit) K a lo r Q ( k J /m e n it ) Set Point Output
Gambar 5. Profil  output  dan  kalor  hasil  perubahan  set point dan gangguan dengan metode optimasi – Sistem 2
0 20 40 60 80 100 120 140 160 180 200 340 350 360 370 Te m p e ra tu re T (K ) 0 20 40 60 80 100 120 140 160 180 200 290 300 310 320 D is tu rb a n c e To ( K ) 0 20 40 60 80 100 120 140 160 180 200 0 5 10 x 10 6 Waktu (menit) K a lo r Q ( k J /m e n it ) Set Point Output
Gambar 6. Profil  output  dan  kalor  hasil  perubahan  set point dan gangguan dengan metode optimasi – Sistem 3
0 20 40 60 80 100 120 140 160 180 200 320 330 340 350 Te m p e ra tu re T (K ) 0 20 40 60 80 100 120 140 160 180 200 290 300 310 320 D is tu rb a n c e To ( K ) 0 20 40 60 80 100 120 140 160 180 200 0 5 10 x 10 6 Waktu (menit) K a lo r Q ( k J /m e n it ) Set Point Output
Gambar 7. Profil  output  dan  kalor  hasil  perubahan  set point dan gangguan dengan metode optimasi – Sistem 4
0 20 40 60 80 100 120 140 160 180 200 320 330 340 350 Te m p e ra tu re T (K ) 0 20 40 60 80 100 120 140 160 180 200 290 300 310 320 D is tu rb a n c e To ( K ) 0 20 40 60 80 100 120 140 160 180 200 0 5 10 x 10 6 Waktu (menit) K a lo r Q ( k J /m e n it ) Set Point Output
Gambar 8. Profil  output  dan  kalor  hasil  perubahan  set point dan gangguan dengan metode optimasi – Sistem 5
Demikian  dilakukan  juga  pengujian  pada Sistem 3 s.d 5 seperti ditunjukkan pada Gambar 6 s.d 8.  Pada  kelima  sistem  tersebut,  kenaikkan  suhu  set point  akan  menyebabkan  jumlah  kalor  yang  disuplai ke  dalam  tangki  harus  ditambah  agar  suhu  keluar dapat mencapai seperti set point atau mencapai  suhu yang  diinginkan,  begitu  pula  sebaliknya.  Selain  itu adanya gangguan berupa perubahan suhu masuk (To) akan  menyebabkan  jumlah  kalor  yang  disuplai berubah.
Tabel  3  menampilkan  harga  SSE  (sum squared  of  error)  yang  diperoleh  pada  kelima  sistem yang ditinjau ketika disetel dengan metode  optimasi. Harga SSE pada tabel  tersebut  dibandingkan  dengan harga  yang  diperoleh  ketika  dicoba  dengan  metode penyetelan  yang  lain  (yaitu Direct  Synthesis dan Ziegler-Nichols).  Tampak  bahwa  harga  SSE  yang diperoleh  dengan  metode  optimasi  adalah  paling kecil  jika  dibandingkan  dengan  kedua  metode  yang lain  (yang  memang  belum  optimum).  Detail  tentang metode Direct  Synthesis dan  Ziegler  Nichols  dapat dilihat pada literature [1,4].
Tabel 3. Nilai SSE Simulasi No SSE Optimasi Direct Synthesis Ziegler- Nichol 1 1070,9 1565,9 1308,7 2 1061,7 1555,4 1329,0 3 1049,3 1372,0 1586,7 4 1070,9 1565,9 1309,4 5 1058,3 1555,6 1340,6
Kesimpulan
Metode  optimasi  telah  digunakan  dalam penelitian ini untuk melakukan penyetelan pengendali PID  pada  sistem  tangki  pemanas  berpengaduk. Optimasi  menggunakan  metode  least  square  non linear merupakan metode efektif untuk menghasilkan nilai  parameter  pengendali  PID  terbaik.  Hal  ini terbukti  melalui  perbandingan  dengan  metode  – metode  lain  bahwa  nilai  SSE  dari  metode  Optimasi paling  kecil.  Dari  hasil  penelitian  tersebut  terbukti bahwa  dengan  melakukan  simulasi  proses  orde  satu sistem  lintas  tertutup  dengan  pengendalian,  hasil penyetelan  metode  Optimasi  didapat  output  yang paling mendekati perubahan set point.

## Ucapan Terimakasih
Penulis  mengucapkan  terimakasih  kepada Universitas  Surabaya  (Ubaya)  dan  mahasiswa: Ronald  Paulus  dan  Bob  Sadewo  yang  telah membantu melakukan penelitian ini.

## Daftar Pustaka
[1]  D.  Coughanowr,  S.E  LeBlanc,  Process  Systems Analysis and Control, 3rd ed, McGraw Hill, New York, 2009.
[2]  T.F  Edgar,  D.M  Himmelblau,  Optimization  of Chemical  Processes,  McGraw  Hill,  New  York, 1989.
[3] J. A. Romero, R. Sanchis, P. Balaguer, Journal of Process Control, 21, (2011) 840–851.
[4] D.E Seborg, T.F Edgar, D.A Mellichamp, Process Dynamics and Control,John Wiley & Sons, New York, 1989.
[5]  G.  Syrcos,  I.  K  Kookos,  Chemical  Engineering and Processing, 44, (2005) 41–49.
[6]  S.  Tavakoli,  M.  Tavakoli.    Department  of Automatic  Control  and  Systems  Engineering, The University of Sheffield, 2003.
[7]  R.  Agustriyanto,    A.  Fatmawati,      Prosiding Seminar  Nasional  Teknik  Kimia  Kejuangan, Yogyakarta, 2014.