PERANCANGAN PENGENDALI MODEL PREDICTIVE CONTROL TANPA CONSTRAINT PADA PROTON EXCHANGE MEMBRANE FUEL CELL

PERANCANGAN PENGENDALI MODEL PREDICTIVE CONTROL TANPA CONSTRAINT PADA PROTON EXCHANGE MEMBRANE FUEL CELL Feri Yusivar, Aries Subiantoro, Suria Departemen Teknik Elektro, Fakultas Teknik Universitas Indonesia Kampus Baru UI, Depok 16424, Indonesia E-mail:[email protected], [email protected], [email protected]

ABSTRAK Pada jurnal ini dirancang pengendali Model Predictive Control tanpa constraint pada sistem Proton Exchange Membrane Fuel Cell. Sebelumnya, dilakukan identifikasi sistem PEMFC dengan metode Least Square berdasarkan data masukan dan keluaran sistem nonlinier untuk mendapatkan model linier. Masukan sistem adalah flow H2, flow O2, dan arus. Kemudian dilakukan uji keabsahan model hasil identifikasi dengan membandingkan respon tegangan keluaran model nonlinier dengan model hasil indentifikasi. Selanjutnya, pengendali MPC tanpa constraint dirancang pada model identifikasi PEMFC dengan nilai parameter matriks faktor bobot Q dan R yang bervariasi. Dari hasil simulasi diperoleh respon pengendalian yang cukup baik meskipun hanya terbatas pada nilai matriks Q= 10 I3 dan R=5000 I6. Hal ini dapat dilihat dari respon yang dapat mengikuti sinyal acuan.. Kata Kunci: Least square, MPC tanpa constraint, model identifikasi, PEMFC, matriks faktor bobot Q dan R ABSTRACT In this journal, a Model Predictive Control without constraint was designed in Proton Exchange Membrane Fuel Cell. Before designing, PEMFC system identification had been done using Least Square method based on input and output data of the nonlinear system to obtain linear model. System inputs are H2 flow, O2 flow, and current. Then, validation test was done by comparing voltage response of nonlinear and identification model. After that, a MPC controller without constraint was designed in identification model. The weighting matrices Q and R were varied to observe their effect to control result. From the simulation result, a quite good control result was obtained even though it’s only limited to Q= 10 I3 and R=5000 I6. The quite good result can be observed from the tracking of setpoint. Keywords: Least square, MPC controller without constraint, identification model, PEMFC, weighting matrices Q and R

I.

PENDAHULUAN Fuel cell sebagai sumber energi alternatif yang terbarukan memiliki prospek yang cerah dalam penerapan ke depan. Fuel cell menggunakan hidrogen sebagai bahan bakar untuk menghasilkan energi dalam bentuk energi listrik dan mengeluarkan emisi air dan panas saja. Dalam penerapannya, fuel cell memiliki beberapa masalah. Salah satunya yaitu tegangan hasil keluaran yang tidak stabil. Maka, perlu adanya skema pengendalian untuk menjaga tegangan yang dihasilkan tetap stabil pada nilai yang diinginkan. Beberapa penelitian yang telah dilakukan mencoba mengendalikan tegangan keluaran fuel cell yang telah melewati inverter terlebih dahulu[8]. Dengan kata lain, tegangan yang dikendalikan adalah tegangan inverter. Pada percobaan ini, perancangan pengendali dilakukan pada aliran masukan H2 dan O2 fuel cell. Pengendali model predictive control (MPC) merupakan salah satu jenis pengendali lanjut yang sudah cukup banyak diterapkan pada aplikasi industri seperti industri pengolahan minyak, industri proses kimia, dan sebagainya. Hal ini disebabkan oleh kesederhanaan algoritmanya dan penggunaan

model non-parametrik yang membutuhkan sedikit informasi a-priori tentang proses dibandingkan model parametrik. Perancangan pengendali menggunakan metode MPC untuk menghasilkan respon tegangan yang lebih halus karena MPC mampu memprediksi perubahan sinyal kendali berikutnya. Namun, perancangan pengendali MPC memiliki kelemahan utama yaitu sangat bergantung pada model sistem. Apabila model sistem yang digunakan kurang baik, maka hasil pengendalian juga menjadi kurang baik. Dengan demikian, perlu dilakukan pemodelan yang baik terlebih dahulu sebelum perancangan pengendali MPC dilakukan. II.

PEMODELAN DINAMIK PEMFC Fuel cell merupakan suatu perangkat konversi energi yang dapat mengubah energi kimia menjadi energi listrik. Fuel cell terdiri dari beberapa bagian antara lain: a. Anoda (fuel electrode), bagian negatif dari fuel cell, berfungsi untuk menghantarkan elektron yang dilepaskan molekul hidrogen sehingga dapat digunakan oleh rangkaian eksternal. Anoda memiliki kanal yang di-etch padanya

b.

c.

d.

yang menyebarkan gas hidrogen secara merata pada permukaan katalis. Katoda (oxygen electrode), bagian positif dari fuel cell, memiliki kanal yang di-etch padanya yang mendistribusikan oksigen ke permukaan katalis. Katoda juga menghantarkan elektron kembali dari rangkaian eksternal menuju katalis, di mana mereka dapat berekombinasi dengan ion hidrogen dan oksigen untuk membentuk air. Elektrolit berupa proton exchange membrane. Ini merupakan bahan khusus, yang terlihat seperti bungkus plastik biasa, hanya menghantarkan ion bermuatan positif dan memblokir elektron. Oleh karena itu, membran ini disebut proton exchange membrane. Untuk suatu PEMFC, membran harus dihidrasi agar dapat berfungsi dan tetap stabil. Katalis, merupakan bahan khusus yang memfasilitasi reaksi antara oksigen dan hidrogen. Umumnya terbuat dari platinum nanoparticles yang dilapis sangat tipis di atas kertas karbon ataupun kain. Katalis ini bersifat kasar dan berpori sehingga permukaan maksimum dari platina dapat diekspos kepada hidrogen ataupun oksigen. Bagian katalis yang dilapisi platina menghadap PEM.

Pada awalnya, hidrogen (H2) dialirkan ke permukaan anoda. Bersamaan dengan itu, oksigen yang berasal dari udara bebas dialirkan ke permukaan katoda. Kemudian elektroda dihubungkan dengan beban dari luar sehingga terjadi proses kimiawi. Hidrogen yang mengalir ke dalam fuel cell akan menyentuh permukaan anoda sehingga H2 bereaksi secara kimiawi berupa reaksi reduksi menghasilkan ion hidrogen (H+) dan elektron (e-). Pada katoda terjadi reaksi oksidasi dan menghasilkan air. Ion hidrogen akan bergerak dari permukaan anoda menuju katoda melalui elektrolit dan elektron bergerak ke beban lalu menuju katoda. Cara kerja fuel cell tersebut dapat dilihat pada gambar 1.

H2 →

2 H + + 2e− ,

O 2 + 4 H + + 4e− → 2 H 2O, 2 H 2 + O 2 → 2 H 2O + electricity + heat. Tegangan yang dihasilkan oleh fuel cell tidak sepenuhnya ideal. Hal ini disebabkan oleh adanya arus yang mengalir. Seiring meningkatnya kerapatan arus, maka tegangan yang dihasilkan semakin turun dari tegangan idealnya. Faktor yang mempengaruhinya adalah rugi aktivasi, rugi ohmik, dan rugi perpindahan massa. Ketiga faktor tersebut dapat dilihat pada gambar 2.

Gambar 2. Kurva Karakteristik Fuel Cell

Tekanan parsial dari hidrogen, oksigen dan air pada sisi katoda didefinisikan sebagai variabel keadaan (state variable) sistem. Pada hukum gas ideal dinyatakan bahwa tekanan parsial setiap gas adalah proporsional terhadap banyaknya gas dalam cell, yang mana dari ketiga kontribusi tersebut relevan dan tergantung dari flow rate masukan, konsumsi gas dan flow rate keluaran gas [1].

Gambar 3 Ilustrasi Aliran Gas Dengan demikian, diperoleh persamaan keadaan yang dinyatakan dalam persamaan (1) dpH 2 RT = ( H 2−in − H 2 −used − H 2− out ) dt VA dpO2 RT = ( O2−in − O2−used − O2− out ) dt VC Gambar 1. Struktur Umum Fuel Cell

Energi listrik yang dihasilkan fuel cell berasal dari reaksi kimia berikut ini.

dpH 2 OC RT = ( H 2OC −in + H 2OC − produce − H 2OC −out ) dt VC (1) di mana H2-in ,O2-in dan H2OC-in berturut-turut adalah flow rate masukan hidrogen, oksigen dan air pada

katoda, sedangkan H2-out ,O2out dan H2OC-out adalah flow rate keluaran hidrogen, oksigen, dan air. H2used ,O2used dan H2OC-produced adalah gas-gas yang digunakan dan dihasilkan Pada sistem elektrokimia, dinyatakan bahwa antara gas yang digunakan dan yang dihasilkan dipengaruhi arus keluaran I, dan dinyatakan dengan persamaan (2): H2-used = 2O2used = 2*Kr*I = 2 KrAcI

(2)

dimana : Kr : N/(4F) Ac : luasan cell yang aktif I : kerapatan arus cell Demikian juga keluaran flow rate didefinisikan dan dinyatakan dalam persamaan (3) : H 2 − out = ( Anodain − 2 KrI ) FH 2 O2 − out = (Cathin − K r I ) FO2

(3)

H 2 Oout = (Cathin + 2 K r I ) FH 2Oc

di mana : Anodain= H2-in Cathin = N2-in + O2-in FH2,FO2 dan FH2O adalah tekanan masingmasing gas di dalam fuel cell, dan dinyatakan dalam persamaan (4) : FH 2 = FO2 =

pH 2 Pop pO 2 P op

FH 2O c =

(4)

pH 2 Oc Pop

Persamaan (3) dan (4) disubstitusi ke persamaan (1) maka diperoleh persamaan (5) yang merupakan persamaan keadaan orde-3 sebagai fungsi dari variabel keadaan dan fungsi masukan : dpH2 RT ⎛ pH2 ⎞ = ⎜ H2−in −2Kr Aci −(Anodain −2Kr Ai ⎟ c ) dt VA ⎜⎝ Pop ⎟⎠ dpO2 RT ⎛ pO2 ⎞ = ⎜O2−in −2Kr Aci −(Cathodain − Kr Ai ⎟ c ) dt VC ⎜⎝ Pop ⎟⎠

(5)

dpH2OC RT ⎛ pH2Oc ⎞ = ⎜ H2OC−in + 2Kr Aci −(Cathodain + 2Kr Ai ⎟ c ) dt VC ⎜⎝ Pop ⎟⎠

Dengan menggunakan persamaan (5) dan mensubstitusi Anodain dan Katodain maka apat diperoleh persamaan (6) seperti di bawah ini:

dpH 2 RT ⎛ pH 2 ⎞ = ⎜ H 2 − in − 2K r A c i − ( H 2 − in − 2 K r Ac i ) ⎟ dt VA ⎜⎝ Pop ⎟⎠ dpO2 RT ⎛ pO2 ⎞ = ⎜ O2 − in − K r A c i − (O2 − in + H 2 OC − in − K r Ac i ) ⎟ dt VC ⎜⎝ Pop ⎟⎠ dpH 2 OC RT ⎛ pH 2 Oc = ⎜ H 2 OC − in + 2K r A c i − (O2 − in + H 2 OC − in + 2 K r Ac i ) dt VC ⎜⎝ Pop

⎞ ⎟ ⎟ ⎠

(6) Dari persamaan perubahan tekananan parsial gas-gas yang digunakan, dapat dituliskan persamaan variabel keadaan dan keluaran model nonlinier sistem dinamik fuel cell PEM, yang dinyatakan dalam persamaan (5) dan persamaan (6) : ⎡ RT ⎢ & ⎡ x1 ⎤ ⎢ VA ⎢ x& ⎥ = ⎢ ⎢ 2⎥ ⎢ ⎢⎣ x&3 ⎥⎦ ⎢ ⎢ ⎢ ⎣

⎛ x 1− 1 ⎜⎜ P op ⎝ 0 0

⎡ RT ⎛ x ⎞⎤ ⎡ ⎤ ⎢ ⎜ -2K r A c (1 − 1 ) ⎟⎟ ⎥ ⎞⎤ ⎢ ⎥ VA ⎜⎝ Pop ⎠ ⎥ ⎢ 0 ⎟⎟ ⎥ ⎢ ⎥ ⎢ ⎥ ⎠⎥ ⎢ ⎛ ⎞⎥ ⎛ ⎞ ⎥ ⎢ RT ⎜ 1 − x2 ⎟ ⎥ u2 + ⎢ RT ⎜ - K r A c (1 − x2 ) ⎟ ⎥ u 3 u + 1 ⎢ V ⎜ ⎥ ⎢ VC ⎜⎝ Pop ⎟⎠ ⎥ Pop ⎟⎠ ⎥ ⎢ C ⎝ ⎥ ⎥ ⎢ ⎥ ⎢ ⎥ ⎥ ⎢ − RT x ⎥ ⎛ x3 ⎞ ⎥ RT ⎢ ⎥ 3 ⎥ )⎟ ⎜ 2K A (1 − ⎢ ⎦ ⎟⎥ ⎢ VC ⎜ r c P ⎣ VC Pop ⎦ op ⎠ ⎦ ⎝ ⎣

(7) Tabel 1: Parameter-parameter yang digunakan pada PEMFC[1] Parameter Simbol Nilai Satuan Konstanta gas

R

8,3144

J/mole*k

Konstanta Faraday Tekanan Kerja

F

96.439

C/mole

Pop / Pstd

100000

Pa

Luas aktif Cell

Ac

136,7

cm2

Volume Anoda

Va

6,495

cm2

Volume Katoda

Vc

12,96

cm2

Banyaknya Cell

N

35

Temperatur kerja

T

353

K

0.8

V

Tegangan tanpa beban Rasio Aliran Hidrogen-Oksigen Kerapatan Arus aktivasi Kerapatan Arus konsentrasi Konstanta Aktivasi Konstanta Konsentrasi Tahanan spesifik sehubungan dengan rugi resistif

0

E

1,1168 Io

1.4e-5

A

Il

7.7

A

a

0.16

V

b

0.063

V

Rloss

0.34

ohm

Pada PEMFC, persamaan tegangan keluaran cell adalah [1]:

V = N ( E0 +

RT ⎧ pH2 ( pO2 / Pstd )0.5 ⎫ ln ⎨ ⎬ − L) 2F ⎩ pH2O ⎭

(8)

dimana : V : tegangan keluaran stack. N : jumlah cell dalam stack E0 : tegangan rangkaian terbuka dari cell T : temperatur kerja L : rugi tegangan pH2, pO2, pH2O: tekanan parsial dari setiap gas dalam cell. R : konstanta gas (8.3144 J/mole oK) F : konstanta Faraday ( 96439 C/mole) Pstd : tekanan standar ( 101325 Pa) L dinyatakan dengan Rugi tegangan persamaan[1] : i + in i + in (9) ) - b ln(1) L = (i + in )r + a ln( io il dimana : i : keluaran kerapatan arus in : kerapatan arus internal terhadap rugi-rugi arus internal io : pertukaran kerapatan arus sehubungan rugirugi aktivasi il : batas kerapatan arus sehubungan rugi-rugi konsentrasi r : luas tahanan spesifik sehubungan rugi-rugi tahanan a,b : konstanta. dimana :

x = [ pH 2

pO2

pH 2 Oc ]

u = [ H 2 − in

O2 −in

i]

T

T

y =V Dari parameter-parameter diatas dilakukan simulasi model state yang didapat dengan output tegangan Nerst. Hasil dari simulasi ini dapat dilihat pada gambar 4.

III.

IDENTIFIKASI SISTEM

Identifikasi adalah penentuan karakteristik dinamik proses secara eksperimental berdasarkan sinyal-sinyal yang terukur ke dalam suatu model matematik. Model proses ditentukan berdasarkan data masukan dan keluaran sistem dengan menggunakan metode Kuadrat Terkecil (Least Square). Inti dari metode ini adalah bahwa kecocokan antara model dengan sistem yang akan diidentifikasi diperoleh dengan meminimumkan selisih kuadrat antara keluaran model dengan keluaran sistem yang diidentifikasi untuk semua N data pengamatan [5]. Selisih kuadrat antara keluaran model dan keluaran sistem dapat dinyatakan dalam persamaan fungsi kriteria sebagai berikut. N

N

i =1

i =1

J LS = ∑ ε i 2 = ∑ ( y (i ) − yˆ (i ) )

2

(10)

di mana

JLS

= fungsi kriteria

εi

= kesalahan prediksi data ke-i

y(i) = data keluaran ke-i yˆ (i) = prediksi keluaran ke-i Fungsi kriteria pada persamaan (10) disebut juga loss function. Nilai loss function yang semakin kecil menunjukkan bahwa model hasil identifikasi semakin mendekati sistem yang sebenarnya. Keluaran model untuk satu langkah prediksi ke depan dari model dinamik orde-n dengan N data pengamatan adalah sebagai berikut ⎡ a1 ⎤ ⎢ ⎥ u (0) L − y (− n) L u ( − n) ⎤ ⎢ M ⎥ ⎡ yˆ (1) ⎤ ⎡ − y (0) ⎢ˆ ⎥ ⎢ ⎥ ⎢a ⎥ − − − − (1) (1 ) (1) (1 ) y y n u u n L L (2) y ⎢ ⎥=⎢ ⎥⎢ n⎥ ⎢ M ⎥ ⎢ ⎥ ⎢ b1 ⎥ M M M M ⎢ ⎥ ⎢ ⎥⎢ ⎥ u ( N − n) M yˆ ( N ) ⎦ ⎣ − y ( N − 1) L − y ( N − n) u ( N − 1) ⎣4 1 24 3 14444444444244444444443⎦ ⎢ ⎥ bn ⎥⎦ ⎢⎣{ Ρ yˆ θˆ

(11) Agar persamaan (11) dapat diminimasi, maka persamaan (11) harus dinyatakan dalam bentuk T J LS = y − Ρθˆ y − Ρθˆ (12) T T T T T T ˆ ˆ ˆ ˆ = y y − θ Ρ y − y Ρθ + θ Ρ Ρθ

(

)(

)

Selanjutnya, dengan membuat turunan pertama J LS (θ ) terhadap θˆ menjadi nol : ∂J LS (θ )

Gambar 4. Hasil Simulasi Model Dinamik PEMFC

∂θ

= −2Ρ y + 2Ρ Ρθˆ = 0 T

T

(13)

θ =θˆ

maka

didapatkan rumus untuk parameter estimasi θˆ sebagai berikut

( )

θˆ = ΡT Ρ

−1

Ρ y T

menghitung (14)

PERANCANGAN PENGENDALI MODEL PREDICTIVE CONTROL

IV.

Model predictive control termasuk dalam kategori konsep perancangan pengendali berbasis model proses di mana model proses digunakan secara eksplisit untuk merancang pengendali dengan meminimasi suatu fungsi kriteria.

B ⎡ ⎤ ⎡ xˆ(k + 1| k ) ⎤ ⎡ A ⎤ ⎢ ⎥ ⎢ ⎥ M ⎢ ⎥ M ⎢ ⎥ M Hu − 1 ⎢ ⎥ ⎢ Hu ⎥ i ⎢ ∑ i =0 A B ⎥⎥ ⎢ xˆ(k + Hu | k ) ⎥ ⎢ A ⎥ ⎢ ⎢ ⎥ = ⎢ Hu +1 ⎥ x(k ) + ⎢ Hu i ⎥ u (k − 1) ⎥ ⎢ xˆ(k + Hu + 1| k ) ⎥ ⎢ A ⎢ ∑ i =0 A B ⎥ ⎢ ⎥ ⎢ M ⎥ ⎢ ⎥ M M ⎢ ⎥ ⎢ ⎥ ⎢ Hp −1 ⎥ i ⎢⎣ xˆ(k + Hp | k ) ⎥⎦ ⎢⎣ AHp ⎥⎦ ⎢∑ A B ⎥⎦ 1 424 3 i =0 ⎣14 4244 3 Ψ Γ 14444442444444 3 Lampau

B L 0 nxl ⎡ ⎤ ⎢ AB + B ⎥ L 0 nxl ⎢ ⎥ ⎢ ⎥ M O M Δuˆ(k ) ⎤ ⎢ Hu −1 i ⎥⎡ ⎥ B ⎥⎢ + ⎢∑ i =0 A B L M ⎥ ⎢ Hu i ⎥⎢ ⎢ ∑ i =0 A B L AB + B ⎥ ⎢⎣ Δuˆ(k + Hu − 1) ⎥⎦ ⎢ ⎥ M M M ⎢ ⎥ ⎢ Hp −1 i ⎥ Hp − Hu i A B A B L ∑ i =0 ⎥ ⎣⎢ ∑ i = 0 144444 42444444 3⎦ Θ 1444444444 24444444443

Formulasi Model Predictive Control Perhitungan sinyal kendali pada MPC dilakukan dengan meminimumkan suatu fungsi kriteria. Fungsi kriteria yang digunakan dalam algoritma MPC berbentuk kuadratik seperti berikut [2]: Hp Hu −1 ) V (k ) = ∑ || y (k + i | k ) − r (k + i | k ) ||Q2 (i ) + ∑ || Δuˆ (k + i | k ) ||2R (i ) i =1

i =0

(15) dengan : ) y (k + i | k )

= keluaran terprediksi untuk i-

langkah kedepan saat waktu k r (k + i | k ) trajectory)

= nilai trayektori acuan (reference

Δuˆ ( k + i | k ) = perubahan nilai sinyal kendali terprediksi untuk i-langkah kedepan saat waktu k

Q(i) dan R(i)

= faktor bobot

Hp

= prediction horizon

Hu

= control horizon

Pemilihan penggunaan Δuˆ ( k + i | k ) yang pada fungsi kriteria bertujuan untuk meminimumkan perubahan sinyal kendali yang masuk ke plant. Model Proses Model proses yang digunakan dalam perancangan pengendali MPC umumnya berupa model ruang keadaan diskrit linier: x(k + 1) = Ax(k ) + Bu (k ) (16)

y (k ) = C x(k )

(17)

Dengan x(k), u(k), dan y(k) adalah variabel keadaan, masukan sistem, dan keluaran sistem, sedangkan A, B, dan C adalah matriks keadaan, matriks masukan, dan matriks keluaran. Persamaan prediksi sepanjang prediction horizon untuk variabel keadaan dapat dinyatakan sebagai

Prediksi

(18) Sedangkan persamaan prediksi untuk keluaran sistem dapat dinyatakan sebagai ⎡ C 0 mxn L 0 mxn ⎤ ⎡ yˆ (k + 1 | k ) ⎤ ⎢ ⎡ xˆ (k + 1 | k ) ⎤ C L 0 mxn ⎥⎥ ⎢ ⎢ ⎥ ⎢0 mxn ⎥ = M M ⎢ ⎥ ⎢ M ⎥ M O M ⎥⎢ ⎢ yˆ (k + Hp | k )⎥ ⎢ ⎥ ⎢⎣ xˆ (k + Hp | k )⎥⎦ ⎣ ⎦ 0 mxn L C ⎦ mxn ⎣0 4 1 444 24444 3 Cy

(19)

Algoritma Pengendali MPC tanpa Constraint 1. Mendefinisikan jumlah masukan sistem (m), jumlah keluaran sistem (p), rentang horizon (Hp=Hu), model proses sistem, titik kerja linier, waktu pencuplikan 2. Membuat matriks dan Q = c1 I ( p × Hp )

R = c 2 I ( m × Hu ) 3.

Membuat matriks Cyψ, CyΓ, Cyθ

4.

Menghitung matriks H = (Cyθ )T QCyθ + R dan

5. 6.

invers H Menghitung matriks E(k) Menghitung matriks G (k ) = 2(Cyθ )T QE ( k )

7.

Mencari

8.

du optimal 1 −1 ΔU ( k ) opt = H G 2 Menghitung sinyal u (k ) = Δu (k ) + u ( k − 1)

dengan

kendali

Nilai u(k) disimpan untuk dipakai sebagai u(k-1) pada k selanjutnya nilai keadaan 10. Menghitung x (k + 1) = Ax ( k ) + Bu (k )

9.

11. Nilai x(k+1) disimpan untuk dipakai sebagai x(k) pada k selanjutnya

12. Menghitung nilai y(k ) = Cx (k ) + Du (k )

keluaran

Hasil pengendalian MPC sangat dipengaruhi oleh kualitas model sistem diskrit yang dikendalikan. Dengan demikian, perlu dilakukan identifikasi terhadap sistem PEMFC sebelum perancangan pengendali MPC dapat dilakukan. Tekanan parsial dari hidrogen, oksigen dan air pada sisi katoda didefinisikan sebagai variabel keadaan (state variable) sistem. Pada model PEMFC pada persamaan (7) dan (8) masukan sistem adalah flow H2, flow O2, dan arus. Keluaran sistem berupa tegangan. Kemudian, dari persamaan model dinamik yang telah diperoleh pada persamaan (7) dan (8) dilakukan identifikasi untuk mendapatkan model linier. Titik kerja yang diambil untuk identifikasi adalah u1ss = 0.22 mol/detik, u2ss = 0.197 mol/detik dan u3ss = 0.549 A/cm2[8]. Dari data masukan dan keluaran sistem nonlinier, diturunkan model ruang keadaan sistem PEMFC dengan rumusan ⎡b11 ⎡ x1 (k + 1) ⎤ ⎡ a11 a12 a13 ⎤ ⎥ ⎡ x1 (k ) ⎤ ⎢ ⎢ ⎥ ⎢ ⎢ x2 (k + 1) ⎥ = ⎢ a21 a22 a23 ⎥ ⎢ x (k ) ⎥ + ⎢b21 ⎢ x3 (k + 1) ⎥ ⎢ a31 a32 a33 ⎥ ⎢ 2 ⎥ ⎢b31 ⎥ ⎢⎣ x3 (k ) ⎥⎦ ⎢ ⎢ ⎥ ⎢ ⎢⎣ y (k ) ⎥⎦ ⎢⎣ a41 a42 a43 ⎥⎦ ⎢⎣b41

b12 b22 b32 b42

b13 ⎤ ⎡ k1 ⎤ ⎥ ⎡u1 (k ) ⎤ ⎢ ⎥ b23 ⎥ ⎢ ⎥ ⎢ k2 ⎥ u2 ( k ) ⎥ + ⎢ k3 ⎥ b33 ⎥ ⎢ ⎥ ⎢⎣u3 (k ) ⎥⎦ ⎢ ⎥ b43 ⎥⎦ ⎢⎣ k3 ⎥⎦

(20) Kemudian dilakukan identifikasi model dengan menggunakan langkah yang telah dijabarkan pada bab IV. Dari proses identifikasi diperoleh persamaan ruang keadaan model linier −1.3549e − 4 0.0053⎤ ⎡ x1 (k ) ⎤ ⎡ x1 (k + 1) ⎤ ⎡ 0.9997 ⎢ ⎥ ⎢ ⎢ ⎥ + = − − ( 1) 4.5208 5 0.9999 0.0012 ⎥⎥ ⎢ x2 (k ) ⎥ x k e 2 ⎢ ⎥ ⎢ ⎢⎣ x3 (k + 1) ⎥⎦ ⎣⎢ −5.375e − 6 −2.9788e − 6 1.0001 ⎦⎥ ⎢⎣ x3 (k ) ⎥⎦

⎡ 42.084 −0.4882 −1.015 ⎤ ⎡u1 (k ) ⎤ ⎡ 0.6664⎤ ⎢ ⎥ + ⎢⎢ 0.1143 21.6643 −0.2663⎥⎥ ⎢u2 (k ) ⎥ + ⎢⎢ 0.1502⎥⎥ ⎢⎣ 0.0132 −0.1024 0.56 ⎥⎦ ⎢⎣u3 (k ) ⎥⎦ ⎢⎣ 0.0178⎥⎦ ⎡ x1 (k ) ⎤ ⎢ ⎥ y (k ) = [ 0.0028 −0.0012 −0.0439] ⎢ x2 (k ) ⎥ ⎢⎣ x3 (k ) ⎥⎦ ⎡u1 (k ) ⎤ ⎢ ⎥ + [ 0.6344 −0.8121 −0.5587 ] ⎢u2 (k ) ⎥ + 26.4224 ⎢⎣u3 (k ) ⎥⎦ (21) Model hasil identifikasi terditri dari 3 masukan. Namun, masukan sistem pada model yang sebenarnya hanyalah flow H2 dan flow O2, sedangkan masukan ketiga yaitu kerapatan arus (i), merupakan umpan balik dari tegangan keluaran sistem yang telah dilewatkan pada beban dan dibagi dengan luas area sel aktif. Oleh karena itu, pengendali MPC yang dirancang hanya akan mengeluarkan dua buah sinyal kendali yaitu sinyal kendali flow H2 (u1) dan sinyal kendali flow O2 (u2). Kerapatan arus diubah menjadi fungsi gangguan df.

Dengan demikian, bentuk persamaan dari model identifikasi berubah menjadi: −1.3549e − 4 0.0053⎤ ⎡ x1 (k ) ⎤ ⎡ x1 (k + 1) ⎤ ⎡ 0.9997 ⎢ ⎥ ⎢ ⎢ ⎥ ( 1) 4.5208 5 0.9999 0.0012 ⎥⎥ ⎢ x2 (k ) ⎥ + = − − x k e ⎢ 2 ⎥ ⎢ ⎢⎣ x3 (k + 1) ⎥⎦ ⎢⎣ −5.375e − 6 −2.9788e − 6 1.0001 ⎥⎦ ⎢⎣ x3 (k ) ⎥⎦ ⎡ 42.084 −0.4882 ⎤ ⎡ −1.015 ⎤ ⎡ 0.6664 ⎤ ⎡ u (k ) ⎤ + ⎢⎢ 0.1143 21.6643 ⎥⎥ ⎢ 1 ⎥ + ⎢⎢ −0.2663⎥⎥ df (k ) + ⎢⎢ 0.1502 ⎥⎥ u (k ) ⎢⎣ 0.0132 −0.1024 ⎥⎦ ⎣ 2 ⎦ ⎢⎣ 0.56 ⎥⎦ ⎢⎣ 0.0178 ⎥⎦

⎡ x1 (k ) ⎤ ⎢ ⎥ y (k ) = [ 0.0028 −0.0012 −0.0439] ⎢ x2 (k ) ⎥ ⎢⎣ x3 (k ) ⎥⎦ ⎡ u (k ) ⎤ + [ 0.6344 −0.8121] ⎢ 1 ⎥ − 0.5587.df ( k ) + 26.4224 ⎣u2 ( k ) ⎦

(22) Untuk menguji kelayakan model hasil identifikasi untuk dijadikan model dasar perancangan pengendali, perlu dilakukan uji observability dan controllability. Uji observability bertujuan untuk mengetahui apakah sistem tersebut benar-benar dapat diobservasi dan untuk mengetahui apakah state-state yang diobservasi tersebut dapat mewakili keadaan sistem yang sebenarnya. Sedangkan uji controllability bertujuan untuk menentukan apakah sistem yang diwakili oleh model tersebut dapat dikendalikan oleh sebuah pengendali. Dengan mengasumsikan model ruang keadaan seperti pada persamaan (16) dan (17), uji observability dapat dilakukan dengan membentuk matriks observability seperti yang ditunjukkan oleh persamaan ⎡C T M AT C T MLM ( AT )n −1 C T ⎤ (23) ⎢⎣ ⎥⎦ Dengan n adalah jumlah state yang dimiliki sistem. Suatu sistem dinyatakan observable jika matriks observability memiliki rank sebanyak n (jumlah state). Untuk melakukan uji controllability dari suatu sistem, langkah yang harus dilakukan adalah membentuk matriks controllability (24) [ B M AB MKM An −1 B ] Setelah dilakukan pengujian observability dan controllability diperoleh bahwa model hasil identifikasi bersifat observable dan controllable sempurna. Langkah berikutnya yaitu menurunkan persamaan sinyal kendali berdasarkan struktur model x(k + 1) = Ax(k ) + Bu (k ) + Bd df (k ) + K x (25) y (k ) = C x(k ) + Du (k ) + Dd df ( k ) + K y Model linier yang digunakan mengandung matriks D (directfeedthrough) pada persamaan keluarannya, maka panjang horizon untuk prediksi dan sinyal kendali adalah sama (Hp=Hu) [2]. Directfeedthrough adalah kondisi sistem di mana keluaran sistem langsung dimasukkan/diumpanbalik ke dalam pengendali. Pada penurunan sinyal

kendali untuk model hasil identifikasi, gangguan yang diasumsikan adalah selisih keluaran tegangan saat ini dengan tegangan sebelumnya. Pada sistem sebenarnya, gangguan berupa perubahan arus akibat perubahan beban. Namun, untuk simulasi ini, beban diasumsikan tetap. Karena tegangan berbanding lurus dengan arus, maka asumsi ini masih dapat digunakan. Nilai Hp dan Hu yang digunakan adalah 3. Setelah dilakukan penjabaran persamaan sinyal kendali, persamaan disusun dalam bentuk matriks keluaran ⎡ y (k + 1 | k ) ⎤ ⎡ Δu (k | k ) ⎤ ⎢ y (k + 2 | k ) ⎥ = C ϕ * x(k ) + C Γ * u (k − 1) + C θ * ⎢ Δu ( k + 1| k ) ⎥ y y y ⎢ ⎥ ⎢ ⎥ ⎣⎢ y (k + 3 | k ) ⎦⎥ ⎣⎢ Δu (k + 2 | k ) ⎦⎥ ⎡ d (k | k ) ⎤ + C y Dw * ⎢⎢ d ( k + 1| k ) ⎥⎥ + C yα + β ⎣⎢ d (k + 2 | k ) ⎦⎥

(26)

di mana: ⎡ CA ⎤ ⎡ ⎤ CB + D ⎢ 2⎥ ⎢ ⎥ C y ϕ = ⎢CA ⎥ , C y Γ = ⎢ CAB + CB + D ⎥, ⎢ 3⎥ ⎢ ⎥ 2 ⎣CA B + CAB + CB + D ⎦ ⎣ CA ⎦ ⎡ CB + D D 0 ⎤ ⎢ ⎥ C yθ = ⎢ CAB + CB + D CB + D D ⎥ ⎢CA2 B + CAB + CB + D CAB + CB + D CB + D ⎥ ⎣ ⎦

⎡ CBd ⎢ C y Dw = ⎢ CABd ⎢ 2 ⎣CA Bd

⎤ ⎥ ⎥ ⎥ CABd CBd + Dd ⎦ ⎡K y ⎤ ⎡ ⎤ CK x ⎢ ⎥ ⎢ ⎥ (27) C yα = ⎢ CAK x + CK x ⎥ , β = ⎢K y ⎥ ⎢ ⎥ ⎢ 2 ⎥ ⎣CA K x + CAK x + CK x ⎦ ⎣⎢ K y ⎦⎥ Matriks penjejakan kesalahan dirumuskan dengan persamaan Dd CBd

0 Dd

E(k ) = T (k ) − C Y Ψ x(k ) − C Y Γu (k − 1) − C Y D w d (k ) − C Y α − β

(28) V.

HASIL SIMULASI Tujuan perancangan pengendali MPC tanpa constraint adalah menjaga tegangan keluaran terhadap perubahan trayektori acuan yang diberikan. Pengendali MPC digunakan untuk mengendalikan flow yang masuk ke dalam PEMFC baik berupa H2 maupun O2. Masukan pengendali MPC berupa nilai state model, tegangan setpoint, dan tegangan keluaran model pada waktu sebelumnya (k-1). Keluarannya adalah sinyal kendali H2 dan O2. Selisih tegangan pada sistem PEMFC dianggap sebagai gangguan. Masukan model hasil identifikasi berupa sinyal kendali hasil keluaran pengendali MPC dan arus. Arus diperoleh dari pembagian tegangan terhadap beban konstan pada 0.3 Ω dan luas area aktif sel. Blok diagram dari sistem PEMFC yang dikendalikan dengan pengendali MPC tanpa constraint dapat dilihat pada gambar 5.

Gambar 5. Blok Diagram Sistem PEMFC yang dikendalikan dengan pengendali MPC Pada simulasi pengendalian MPC tanpa constraint, dilakukan variasi terhadap nilai parameter matriks faktor bobot perubahan sinyal kendali R, matriks faktor bobot kesalahan prediksi Q, dan prediction horizon Hp beserta control horizon Hu. Pada simulasi pertama, dilakukan percobaan dengan nilai R bervariasi, nilai Q tetap pada 10 IHp, dan nilai Hp & Hu tetap pada 3. Percobaan dilakukan pada variasi nilai R sebesar 1000, 5000, 10000, 50000, 100000, 1000000. Simulasi dijalankan dan diperoleh hasil pengendalian dengan perbandingan respon transien pada gambar 6 dan 7.

Gambar 6. Perbandingan Steady State Error pada Nilai R Bervariasi

Gambar 7. Perbandingan Respon Transien pada Nilai R Bervariasi

Hasil simulasi terbaik diperoleh pada nilai parameter matriks R terbaik terdapat pada R=5000 I6 dengan nilai settling time, %overshoot, dan steady state error yang masih cukup baik yakni Ts= 184.7 detik, % OS=76.04%, dan steady state error sebesar 0.25. Pada R=1000 I6, nilai settling time tidak ada dan error steady state semakin lama dapat mencapai tak hingga karena sistem tidak pernah mencapai kestabilan. Dengan demikian, pengendali MPC tanpa constraint ini tidak dapat mengendalikan sistem dengan baik di semua nilai R. Grafik respon sistem untuk R=5000 I6 dan Q=10I3 dapat dilihat pada gambar 8.

Gambar 8. Respon Pengendalian pada R 5000 dan Q 10 Sinyal kendali dari respon pengendalian pada gambar 8 dapat dilihat pada gambar 9.

MPC tanpa constraint yang dirancang masih kurang robust. Pengaruh nilai matriks R pada pengendali MPC tanpa constraint menyebabkan nilai perubahan sinyal kendali yang berubah-ubah. Semakin besar nilai R yang diberikan, maka sinyal kendali yang diberikan ke sistem akan semakin ditekan. Pada nilai R yang semakin kecil, sinyal kendali semakin dilepas sehingga sinyal kendali akan diberikan lebih cepat dibandingkan dengan nilai R yang lebih besar serta memiliki amplitudo yang lebih besar. Pada simulasi kedua, dilakukan percobaan dengan nilai Q bervariasi, nilai R tetap pada 5000 I6, dan nilai Hp & Hu tetap pada 3. Percobaan dilakukan pada variasi nilai Q sebesar 1, 5, 10, 25, 50, 100. Dari hasil simulasi diperoleh bahwa Q yang bervariasi akan mempengaruhi perubahan sinyal kendali menjadi semakin besar dan tidak ditekan. Pada hasil simulasi dengan Q=25 I3, Q=50 I3, Q=100 I3, respon sistem tidak lagi dapat mengikuti setpoint dengan baik serta tidak pernah mencapai nilai steady state. Pada nilai Q lainnya, respon yang diperoleh masih cukup baik. Respon pengendalian optimal diperoleh pada nilai Q sebesar 10 I3. Pada simulasi ketiga, dilakukan percobaan dengan nilai Q dan R tetap pada 10I3 dan 5000I6, serta nilai Hp & Hu bervariasi pada nilai 6, 9, dan 12. Hasil respon pengendalian dan sinyal kendali dengan Hp dan Hu yang bervariasi dapat dilihat pada gambar 10, 11, dan 12.

Gambar 9. Sinyal Kendali pada R 5000 dan Q 10 Hasil penjejakan sinyal acuan sudah cukup baik. Namun, perbandingan nilai parameter respon sistem menunjukkan % overshoot yang buruk karena dapat merusak sistem akibat perubahan respon yang besar, respon pengendalian mengalami osilasi serta masih memiliki steady state error sehingga masih terdapat banyak kekurangan apabila pengendali ini diterapkan langsung pada alat. Hal ini disebabkan oleh nilai sinyal kendali dan respon tegangan yang dapat mencapai nilai negatif, Hp dan Hu yang bernilai kecil, model identifikasi yang tidak sepenuhnya dapat merepresentasikan sistem nonlinier yang sebenarmya, dan kekurangan yang dimiliki oleh pengendali MPC tanpa constraint. Oleh karena itu, pengendali

Gambar 10. Hasil Pengendalian pada Hp & Hu 6

Gambar 11. Hasil Pengendalian pada Hp & Hu 9

pemodelan yang lebih baik daripada metode kuadrat terkecil.

Gambar 12. Hasil Pengendalian pada Hp& Hu 12 Dari hasil simulasi, diperoleh bahwa respon pengendalian mengalami osilasi yang semakin berkurang seiring dengan peningkatan nilai Hp dan Hu seperti yang terlihat pada gambar 10, 11, dan 12. Hal ini disebabkan oleh nilai prediction horizon yang semakin besar. Nilai prediction horizon yang semakin besar menyebabkan perhitungan nilai prediksi perubahan sinyal kendali yang semakin banyak sehingga perubahan sinyal kendali yang terjadi lebih baik dan lebih halus (ditunjukkan dengan jumlah osilasi yang semakin sedikit). Respon pengendalian dan nilai sinyal kendali yang diperoleh masih memiliki osilasi. Hal ini diakibatkan oleh model hasil identifikasi yang masih memiliki perbedaan dengan model nonlinier, tidak adanya batasan dalam perancangan pengendali MPC, dan perancangan pengendali dilakukan pada model hasil identifikasi sehingga hasil pengendalian tidak cukup baik dan pengendali masih kurang robust. Oleh karena itu, diperlukan adanya teknik pemodelan yang lebih baik dan perancangan pengendali MPC nonlinier dengan constraint pada model nonlinier sehingga hasil pengendalian yang diperoleh lebih baik dan dapat diterapkan langsung pada plant PEMFC. VI.

KESIMPULAN Perancangan pengendali MPC sangat dipengaruhi oleh model yang digunakan karena pengendali MPC berbasis model proses. Untuk memperoleh hasil pengendalian yang baik, model yang digunakan untuk merancang pengendali harus baik pula. Hasil pengendalian MPC tanpa constraint tidak robust karena pengendali MPC tanpa constrait tidak memperhitungkan adanya batasan sinyal kendali sehingga nilai sinyal kendali bisa mencapai nilai yang tinggi dan dapat membahayakan sistem. Hasil pengendalian MPC tanpa constraint terbaik diperoleh pada parameter Q= 10 I3 dan R = 5000 I6 dengan nilai %OS sebesar 76,29%, settling time sebesar 209,7 detik, dan steady state error sebesar 0,1. Untuk memperoleh hasil perancangan pengendali yang lebih baik diperlukan pengendali MPC dengan constraint ataupun MPC nonlinier dan teknik

PUSTAKA [1] Aryani, Dharma. Identifikasi Sistem PEMFC (Proton Exchange Membrane Fuel Cell) dengan Metode Kuadrat Terkecil. Seminar Fakultas Teknik Program Studi Kontrol Industri Universitas Indonesia. Depok. Desember 2008. [2] Aryani, Dharma. Perancangan Pengendali Model Predictive Control Constrained pada Sistem Proton Exchange Membrane (PEMFC). Thesis Fakultas Teknik Program Studi Kontrol Industri Universitas Indonesia. Depok. Juni 2009. [3] Kim, Y.H et al. A Fuel Cell System with ZSource Inverters and Ultracapacitors. School of Electrical and Electronics Engineering, ChungAn University. Korea. [4] Nice, Karim, and Jonathan Strickland. How Fuel Cells Work. 18 September 2000. HowStuffWorks.com. 15 December 2009. [5] Subiantoro, Aries. Program Insentif Lecture Note Fakultas Teknik Universitas Indonesia Gasal/Genap 2006/2007. Depok. Agustus 2007. [6] Wingelaar, P.J.H, Duarte, J.L, Hendrix, M.A.M. Dynamic Characteristic of Fuel Cells. Eindhoven University of Technology. Netherlands. 2005. [7] Yusivar, Feri. Laporan Akhir Pengembangan Sistem Energi Baru dan Terbarukan Berbasis pada Teknologi Fuel Cell Program Insentif Riset Terapan Bidang Fokus: Sumber Energi Baru dan Departemen Teknik Elektro Terbarukan. Fakultas Teknik Universitas Indonesia. Depok. November 2008. [8] Na, Woonki & Gou, Bei and Diong, Bill. Nonlinear Control of PEM Fuel Cell by Exact Linearization. IEEE 2005.