TESIS
KENDALI MODEL PREDIKTIF TERDISTRIBUSI BERBASIS TEORI PERMAINAN DINAMIS KOOPERATIF
DISTRIBUTED MODEL PREDICTIVE CONTROL BASED ON COOPERATIVE DYNAMIC GAME THEORY
SUTRISNO 10/306185/PPA/3233
PROGRAM STUDI S2 MATEMATIKA FAKULTAS MATEMATIKA DAN ILMU PENGETAHUAN ALAM UNIVERSITAS GADJAH MADA YOGYAKARTA
2012
TESIS
KENDALI MODEL PREDIKTIF TERDISTRIBUSI BERBASIS TEORI PERMAINAN DINAMIS KOOPERATIF
DISTRIBUTED MODEL PREDICTIVE CONTROL BASED ON COOPERATIVE DYNAMIC GAME THEORY
Diajukan untuk memenuhi salah satu syarat memperoleh derajat Master of Science Matematika
SUTRISNO 10/306185/PPA/3233
PROGRAM STUDI S2 MATEMATIKA FAKULTAS MATEMATIKA DAN ILMU PENGETAHUAN ALAM UNIVERSITAS GADJAH MADA YOGYAKARTA
2012
i
HALAMAN PENGESAHAN TESIS KENDALI MODEL PREDIKTIF TERDISTRIBUSI BERBASIS TEORI PERMAINAN DINAMIS KOOPERATIF Telah dipersiapkan dan disusun oleh
SUTRISNO 10/306185/PPA/3233
Telah dipertahankan di depan Tim Penguji pada tanggal 27 Juli 2012
Susunan Tim Penguji
Dr. Salmah, M.Si. Pembimbing I
Dr. Fajar Adi Kusumo, M.Si. Ketua Tim Penguji
Dr. Indah Emilia W., M.Si. Pembimbing II
Dr. Ari Suparwanto, M.Si. Penguji
Tesis ini telah diterima sebagai salah satu persyaratan untuk memperoleh gelar Master of Science Tanggal Agustus 2012
Prof. Dr.rer.nat. Widodo, MS. Ketua Program Studi Monodisiplin S2/S3 Matematika ii
PERNYATAAN
Dengan ini saya menyatakan bahwa dalam Tesis ini tidak terdapat karya yang pernah diajukan untuk Tesis maupun untuk memperoleh gelar Master di suatu Perguruan Tinggi, dan sepanjang pengetahuan saya juga tidak terdapat karya atau pendapat yang pernah ditulis atau diterbitkan oleh orang lain, kecuali yang secara tertulis diacu dalam naskah ini dan disebutkan dalam daftar pustaka.
Yogyakarta, 27 Juli 2012
Sutrisno
iii
Kupersembahkan karya ini untuk : Ibu dan Ayah Segenap keluarga tercinta
iv
PRAKATA
Puji syukur penulis panjatkan atas rahmat yang dilimpahkan oleh Allah SWT, sehingga penulis dapat menyelesaikan tesis ini dengan judul “Kendali Model Prediktif Terdistribusi Berbasis Teori Permainan Dinamis Kooperatif”. Tesis ini disusun sebagai salah satu syarat untuk memperoleh gelar Master of Science di Program Studi S2 Matematika, Fakultas Matematika dan Ilmu Pengetahuan Alam, Universitas Gadjah Mada. Penulis menyampaikan terima kasih yang sebanyak-banyaknya kepada pihak-pihak yang telah berperan dalam penyusunan tesis ini, yaitu kepada : 1. Dr. Chairil Anwar selaku Dekan Fakultas MIPA UGM. 2. Dr. Lina Aryati, M.S selaku Ketua Jurusan Matematika Fakultas MIPA UGM. 3. Prof. Dr. Widodo, M.S selaku Ketua Pengelola Program Studi S2 Matematika Fakultas MIPA UGM. 4. Dr. Salmah, M.Si selaku Dosen Pembimbing I yang telah membimbing dan mengarahkan penulis dalam penyusunan tesis ini. 5. Dr. Indah Emilia W., M.Si selaku Dosen Pembimbing II yang telah membimbing dan mengarahkan penulis dalam penyusunan tesis ini. 6. Dr. Fajar Adi Kusumo, M.Si dan Dr. Ari Suparwanto, M.Si yang telah menjadi penguji untuk tesis ini. 7. Segenap Dosen Jurusan Matematika FMIPA UGM yang telah membagikan ilmunya selama masa studi. 8. Segenap
karyawan/karyawati
dan
civitas
akademika
Jurusan
Matematika FMIPA UGM yang telah mendukung perkuliahan penulis selama masa studi. 9. Segenap keluarga tercinta yang menjadi motivasi utama bagi penulis. 10. Rekan-rekan mahasiswa S2 Matematika angkatan 2009, 2010 dan 2011 yang memberikan dukungan dan motivasi selama masa studi.
v
11. Semua pihak yang telah memberikan bantuan dan dukungan sehingga penulis dapat menyelesaikan studi S2 Matematika di Fakultas MIPA UGM.
Semoga karya ini dapat bermanfaat untuk pembaca dan untuk kemajuan ilmu Matematika, khususnya teori kendali. Kritik dan saran dari pembaca akan penulis terima dengan baik sebagai bahan evaluasi.
Yogyakarta, 27 Juli 2012
Penulis
vi
DAFTAR ISI
HALAMAN JUDUL ........................................................................................... i HALAMAN PENGESAHAN ............................................................................. ii HALAMAN PERNYATAAN ............................................................................ iii HALAMAN PERSEMBAHAN ......................................................................... iv PRAKATA ......................................................................................................... v DAFTAR ISI .................................................................................................... vii DAFTAR LAMBANG DAN SINGKATAN ....................................................... x DAFTAR TABEL .............................................................................................. xi DAFTAR GAMBAR ........................................................................................ xii DAFTAR LAMPIRAN .................................................................................... xiii INTISARI ........................................................................................................ xiv ABSTRACT ..................................................................................................... xv
BAB I. PENDAHULUAN .................................................................................. 1 1.1 Latar Belakang dan Masalah ............................................................. 1 1.2 Tujuan dan Manfaat Penelitian .......................................................... 3 1.3 Tinjauan Pustaka ............................................................................... 4 1.4 Metode Penelitian ............................................................................. 6 1.5 Sistematika Penulisan ....................................................................... 7
BAB II. LANDASAN TEORI ............................................................................ 8 2.1 Derivatif dan Faktorisasi Matriks ....................................................... 8 2.2 Kekontinuan Fungsi Bernilai Skalar ................................................... 9 2.3 Optimisasi Fungsi Konveks ............................................................. 14 2.3.1 Himpunan dan Fungsi Konveks ......................................... 14 2.3.2 Fungsi Berbentuk Kuadratik .............................................. 18 2.3.3 Optimisasi Fungsi Konveks Tanpa Kendala ....................... 20 2.3.4 Pemrograman Kuadratik Berkendala .................................. 23 vii
2.4 Sistem Diskrit LTI ........................................................................... 29 2.4.1 Kestabilan Lyapunov Sistem Diskrit LTI ........................... 30 2.4.2 Keterkendalian Sistem Diskrit LTI .................................... 35 2.4.3 Keterstabilan Sistem Diskrit LTI ....................................... 36 2.5 Teknik Kendali MPC ...................................................................... 37 2.5.1 Ide Dasar Teknik Kendali MPC ......................................... 37 2.5.2 Strategi Teknik Kendali MPC ............................................ 38 2.5.3 Model Matematika dari Plant ............................................ 39 2.5.4 Fungsi Objektif ................................................................. 43 2.5.5 Hukum Kendali ................................................................. 44 2.5.6 Solusi Teknik Kendali MPC Tanpa Kendala ...................... 45 2.5.7 Solusi MPC Berkendala .................................................... 49 2.6 Teori Permainan Dinamis Kooperatif ............................................... 55 2.6.1 Solusi Pareto ..................................................................... 56 2.6.2 Solusi Nash-Bargaining ..................................................... 64
BAB III. KENDALI MODEL PREDIKTIF TERDISTRIBUSI BERBASIS TEORI PERMAINAN DINAMIS KOOPERATIF ......... 69 3.1 Model Systemwide Plant .................................................................. 69 3.1.1 Model Terdesentralisasi .................................................... 70 3.1.2 Model Interaksi ................................................................. 70 3.1.3 Model Composite .............................................................. 70 3.2 Model Prediksi Systemwide Plant .................................................... 71 3.3 Fungsi Objektif Systemwide Plant .................................................... 74 3.4 Feasible-Cooperation MPC ............................................................. 76 3.5 Skema Nash-Bargaining-MPC Terdistribusi .................................... 81 3.5.1 Solusi Nash-Bargaining MPC Terdistribusi ........................ 82 3.5.2 Model Negosiasi Nash-Bargaining FC-MPC ..................... 83 3.6 Sifat-sifat FC-MPC Sebagai Permainan ........................................... 84 3.6.1 Fisibilitas Algoritma MPC Terdistribusi ............................ 85 3.6.2 Kestabilan Sistem Lingkar Tertutup MPC Terdistribusi ..... 87
viii
BAB IV. APLIKASI KENDALI MODEL PREDIKTIF TERDISTRIBUSI PADA KANAL IRIGASI .................................................................. 93 4.1 Model Matematika Kanal Irigasi ...................................................... 93 4.1.1 Model Terdesentralisasi Kanal Irigasi ................................ 95 4.1.2 Model Interaksi Kanal Irigasi ............................................ 96 4.1.3 Model Composite Kanal Irigasi ......................................... 97 4.2 Deskripsi Systemwide Kanal Irigasi ............................................... 104 4.3 Fungsi Objektif Kanal Irigasi ......................................................... 105 4.4 Aplikasi FC-MPC Pada Kanal Irigasi ............................................ 106 4.4.1 Masalah Optimisasi FC-MPC Pada Kanal Irigasi ............ 106 4.4.2 Hasil Simulasi FC-MPC Pada Kanal Irigasi ..................... 107 4.5 Aplikasi Nash-Bargaining Pada Kanal Irigasi ................................ 111 4.5.1 Model Negosiasi Nash-Bargaining MPC Pada Kanal Irigasi ........................................................... 111 4.5.2 Hasil Simulasi Model Negosiasi Nash-Bargaining MPC Pada Kanal Irigasi ........................................................... 111 4.6 Perbandingan Performansi FC-MPC dan NB-MPC Pada Kanal Irigasi .......................................................................... 114
BAB V. PENUTUP ........................................................................................ 116 5.1 Kesimpulan ................................................................................... 116 5.2 Saran ............................................................................................. 117
DAFTAR PUSTAKA LAMPIRAN
ix
DAFTAR LAMBANG DAN SINGKATAN
: Himpunan semua bilangan real n
: Ruang vektor real berdimensi n : Himpunan semua bilangan asli
N ( x0 )
: Persekitaran titik x0 dengan jari-jari
N d ( x0 )
: Persekitaran terhapus dari titik x0 dengan jari-jari
S
: Klosur himpunan S
S
: Himpunan semua titik batas himpunan S : Himpunan kosong
f ( x0 )
: Gradien fungsi f di titik x0 : Norma Euclid
z
: Modulus dari bilangan kompleks z
Hp
: Horison prediksi
Hu
: Horison kendali
Hw
: Horison window
Hc
: Periode kendali
xˆ ( k )
: State prediksi dari x( k )
uˆ ( k )
: Input prediksi dari u (k )
M
: Himpunan bilangan-bilangan integer 1, 2,3,..., M . : Akhir suatu bukti
LTI
: Linear Time Invariant
MPC
: Model Predictive Control
FC-MPC
: Feasible Cooperation – MPC
NB-MPC
: Nash-Bargaining – MPC
x
DAFTAR TABEL
Tabel 4.2.1 Luas permukaan reach kanal irigasi .............................................. 104 Tabel 4.6.1 Perbandingan cost FC-MPC dan NB-MPC kanal irigasi ................ 115
xi
DAFTAR GAMBAR
Gambar 2.1 (a) Himpunan konveks K1 , dan (b) Himpunan tidak konveks K2 .... 14 Gambar 2.2 Ilustrasi teorema separasi untuk
2
............................................... 15
Gambar 2.3 Contoh grafik fungsi konveks satu variabel .................................... 17 Gambar 2.4.1 Ilustrasi kestabilan sistem single state (a) xe (b) xe
0 stabil asimtotik, dan (c) xe
0 stabil,
0 tidak stabil ................ 31
Gambar 2.5.1 Strategi teknik kendali MPC ....................................................... 38 Gambar 2.5.2 Diagram blok strategi teknik kendali MPC .................................. 39 Gambar 2.5.2 Grafik input dan output contoh (2.5.1) ........................................ 54 Gambar 2.6.1 Pareto frontier contoh (2.6.1) ...................................................... 64 Gambar 2.6.2 (a) Konsep Nash-Bargaining untuk 2 pemain (b) Solusi Nash-Bargaining N(S,d) ............................................. 65 Gambar 2.6.3 Titik disagreement dan solusi Nash-Bargaining contoh (2.6.2) .... 68 Gambar 4.1.1 Kanal irigasi 5 reach ................................................................... 93 Gambar 4.1.2 Ilustrasi masing-masing reach ..................................................... 94 Gambar 4.4.1 Hasil simulasi FC-MPC dengan wi sama ................................... 108 Gambar 4.4.2 Hasil simulasi FC-MPC dengan wi proporsional terhadap beban kendali (volume reach) .................................... 110 Gambar 4.5.1 Hasil simulasi NB-MPC dengan wi sama .................................. 112 Gambar 4.5.2 Hasil simulasi NB-MPC dengan wi proporsional terhadap beban kenali (volume reach) ...................................... 113
xii
DAFTAR LAMPIRAN Lampiran 1. Program MATLAB Contoh (2.5.1) .............................................. 121 Lampiran 2. Program MATLAB Contoh (2.6.1) .............................................. 122 Lampiran 3. Program MATLAB Contoh (2.6.2) .............................................. 123 Lampiran 4. Flow chart algoritma active-set .................................................... 124 Lampiran 5. Flow chart algoritma FC-MPC .................................................... 125 Lampiran 6. Flow chart NB-MPC ................................................................... 126 Lampiran 7. Program MATLAB simulasi FC-MPC kanal irigasi .................... 128 Lampiran 8. Program MATLAB simulasi NB-MPC kanal irigasi .................... 137
xiii
INTISARI Kendali Model Prediktif Terdistribusi Berbasis Teori Permainan Dinamis Kooperatif oleh Sutrisno 10/306185/PPA/3233 Masalah kendali pada systemwide plant merupakan masalah kendali pada sistem yang memuat beberapa subsistem yang saling berinteraksi satu sama lain. Kendali Model Prediktif (Model Predictive Control / MPC) merupakan teknik kendali optimal yang menentukan kendali optimal dengan membuat state dan input prediksi sepanjang horison, kemudian mengaplikasikan elemen pertama pada barisan kendali optimal yang dihitung menggunakan teknik optimisasi. Teknik kendali MPC klasik (hanya untuk sebuah sistem) dapat diaplikasikan pada systemwide menggunakan pendekatan desentralisasi atau sentralisasi, akan tetapi pendekatan desentralisasi mengabaikan interaksi sehingga dapat mengakibatkan performansi yang buruk, sedangkan pendekatan sentralisasi tidak fleksibel dan memuat perhitungan komputasi yang relatif besar. Untuk mengatasi hal ini, teknik kendali MPC diaplikasikan pada systemwide menggunakan pendekatan terdistribusi. Teknik kendali MPC terdistribusi dibentuk berdasarkan teori permainan dinamis kooperatif dimana subsistem dianggap sebagai pemain, kemudian diselesaikan dengan menentukan solusi Pareto (Feasible-CooperationMPC), yang dilanjutkan dengan menentukan solusi Nash-Bargaining (NashBargaining-MPC), yaitu solusi Pareto yang ditentukan berdasarkan titik disagreement para pemain. Solusi Pareto dan solusi Nash-Bargaining ditentukan menggunakan pembobotan fungsi objektif. Dapat dibuktikan bahwa MPC terdistribusi menghasilkan titik equilibrium 0 yang stabil asimtotik. Teknik kendali MPC terdistribusi diaplikasikan pada kanal irigasi 5 subsistem yang memuat cabang. Aplikasi ini disimulasikan menggunakan program MATLAB untuk menyelidiki performansi teknik kendali ini. Berdasarkan hasil simulasi ini, solusi Nash-Bargaining-MPC dengan pembobotan proporsional terhadap beban kendali menghasilkan total cost yang paling kecil. Kata kunci : MPC terdistribusi, kendali systemwide, permainan dinamis.
xiv
ABSTRACT Distributed Model Predictive Control Based On Cooperative Dynamic Game Theory By Sutrisno 10/306185/PPA/3233 Systemwide control problem is a problem that contains several subsystems and involving interaction each other. Model Predictive Control (MPC) is an optimal control technique that calculate the optimal control action by makes the state and input prediction along horizon, then apply the first element of the optimal control sequence that calculated using optimization technique. The classic MPC (only for single system) can be applied for systemwide using decentralized or centralized approach, but the decentralized approach ignoring the interaction, so it can made the bad performance, whereas the centralized approach is not flexible and involving the large scale computation. To handle these problems, MPC can be applied for each subsystem using distributed approach. Distributed MPC is based on cooperative dynamic game theory, where the subsystems is reputed as players, then this problem solved by finding the Pareto solution (Feasible-Cooperation-MPC), that continued by finding the Nash-Bargaining solution (Nash-Bargaining-MPC), namely the Pareto solution that determined according to the disagreement point of players. The Pareto and Nash-Bargaining solution was determined using weighting of objective functions. We can prove that distributed MPC gives equilibrium point 0 is asymptotical stable. Distributed MPC was applied to control an irrigation canal that has 5 subsystems and involves a branch. This study case was simulated using MATLAB to illustrate the performance of this control technique. From the results, the Nash-BargainingMPC with the proportional weighting to the control load gives the fewest total cost. Keywords: Distributed MPC, systemwide control, dynamic game.
xv
BAB 1 PENDAHULUAN
Seiring dengan berkembangnya teknologi, suatu sistem dibangun dengan ukuran yang relatif besar (Large Scale System). Sistem yang berskala besar biasanya terdiri dari subsistem-subsistem yang memuat interaksi satu sama lain yang disebut sebagai systemwide. Interaksi yang terjadi dapat berupa aliran energi, material, informasi, dan lain-lain.
1.1. Latar Belakang dan Permasalahan Untuk melakukan kendali pada systemwide, diperlukan teknik kendali yang cocok dan seoptimal mungkin. Teknik kendali model prediktif (Model Predictive Control/MPC) adalah teknik kendali optimal yang menentukan kendali optimal dengan membuat state dan input prediksi sepanjang horison, kemudian mengaplikasikan elemen pertama pada barisan kendali optimal yang dihitung menggunakan teknik optimisasi. Teknik kendali MPC merupakan teknik kendali untuk sistem diskrit yang dapat diaplikasikan pada sebuah sistem saja. Teknik kendali MPC dapat diaplikasikan pada systemwide dengan cara menerapkan MPC pada masing-masing subsistem secara independen. Teknik kendali ini disebut teknik kendali MPC terdesentralisasi (Decentralized MPC). Teknik kendali tersebut mengabaikan interaksi antar subsistem. Apabila interaksi antar subsistem lemah, teknik kendali MPC terdesentralisasi dapat menghasilkan performansi yang baik. Tetapi apabila interaksinya kuat, mengabaikan interaksi dapat mengurangi performansi sistem dan dapat berakibat fatal pada sistem secara keseluruhan. Sebagai contoh, apabila interaksi pada sistem pembangkit listrik diabaikan, maka sistem akan mati total ketika terjadi kerusakan pada salah satu subsistem. Untuk mengatasi interaksi antar subsistem, dikembangkan teknik kendali MPC tersentralisasi (Centralized MPC), yaitu satu buah MPC diaplikasikan untuk mengendalikan seluruh systemwide dengan menyertakan interaksi. Teknik kendali MPC tersentralisasi sangat sukar untuk diterapkan karena tidak fleksibel, yaitu
1
2
ketika salah satu subsistem membutuhkan perbaikan, maka subsistem yang lain juga harus dimatikan. Teknik kendali ini juga sukar dihitung secara komputasional karena menghasilkan optimisasi dengan dimensi yang relatif besar, sehingga jarang sekali diterapkan. Untuk mengatasi masalah yang tidak dapat ditangani oleh teknik kendali MPC terdesentralisasi dan MPC tersentralisasi, dikembangkan teknik kendali MPC terdistribusi (Distributed MPC). Ide dasar dari MPC terdistribusi diperoleh dari teori permainan dinamis kooperatif. Teori permainan dinamis kooperatif merupakan masalah permainan yang memuat beberapa pemain yang bekerjasama untuk mengoptimalkan fungsi objektif semua pemain. Sedangkan masalah kendali systemwide memuat beberapa subsistem yang akan meminimumkan fungsi objektif semua subsistem. Sehingga masalah kendali systemwide dapat dipandang sebagai permainan dinamis kooperatif dimana pemain-pemainnya adalah subsistem. Teknik kendali MPC terdistribusi menganggap subsistem sebagai pemain, menerapkan MPC pada masing-masing subsistem, kemudian menyelesaikannya secara terdistribusi, yaitu masing-masing MPC secara bersama-sama menentukan kendali yang mengoptimalkan semua fungsi objektif. Kendali optimal diperoleh dari masalah optimisasi yang dihasilkan oleh MPC terdistribusi. Pertama, kendali optimal dihitung dengan meminimumkan kombinasi linier konveks fungsi objektif semua subsistem yang disebut sebagai Feasible-Cooperation MPC (FC-MPC) yang akan memberikan solusi Pareto. Karena solusi Pareto tidak tunggal, selanjutnya dihitung solusi Nash-Bargaining yaitu solusi Pareto yang ditentukan berdasarkan titik disagreement masing-masing subsistem. Solusi Nash-Bargaining pada MPC terdistribusi selanjutnya disebut sebagai Nash-Bargaining MPC (NBMPC) yang ditentukan menggunakan model negosiasi Nash-Bargaining. Penelitian ini akan mambahas perumusan teknik kendali MPC terdistribusi dan aplikasinya pada kanal irigasi. Sehingga, masalah-masalah yang akan dibahas pada penelitian ini adalah sebagai berikut. 1) Memodelkan persamaan state space dari suatu systemwide yang memuat interaksi antar subsistem.
3
2) Memodelkan masalah MPC terdistribusi untuk systemwide. 3) Menyusun optimisasi FC-MPC dan NB-MPC berdasarkan teori permainan dinamis kooperatif. 4) Menyelidiki sifat-sifat MPC terdistribusi yang meliputi kestabilan dan fisibilitas (kondisi suatu iterasi dapat dilakukan atau tidak). 5) Mengaplikasikan MPC terdistribusi pada kanal irigasi kemudian mensimulasikannya menggunakan program MATLAB dan menyelidiki bagaimana MPC terdistribusi bekerja dan mengamati performansinya.
1.2 Tujuan dan Manfaat Penelitian Tujuan dari penelitian ini adalah sebagai berikut. 1) Menyajikan model berbentuk persamaan state space dari suatu systemwide. 2) Mengaplikasikan teknik kendali MPC terdistribusi untuk systemwide. 3) Memperoleh
bentuk
optimisasi
FC-MPC
dan
NB-MPC
untuk
menyelesaikan MPC terdistribusi. 4) Menentukan sifat-sifat kestabilan dan fisibilitas dari FC-MPC dan NBMPC. 5) Mendapatkan model systemwide dari sistem kanal irigasi yang terdiri dari 5 subsistem dan membuat program MATLAB untuk mensimulasikan hasil aplikasi MPC terdistribusi dan membandingkan performansi dari FC-MPC dan NB-MPC.
Manfaat yang dapat diperoleh dari penelitian ini adalah sebagai berikut. 1) Memberikan sumbangan terhadap pengembangan teori kendali diskrit, terutama pengembangan teknik kendali MPC terdistribusi. 2) Memberikan gambaran kepada pembaca tentang pemilihan teknik kendali diskrit yang sesuai untuk suatu systemwide. 3) Program MATLAB yang dibuat diharapkan dapat membantu untuk mensimulasikan aplikasi FC-MPC dan NB-MPC pada suatu systemwide.
4
4) Hasil simulasi yang dilakukan diharapkan dapat memberikan gambaran tentang performansi dari teknik kendali MPC terdistribusi. 5) Memberikan alternatif teknik kendali MPC terdistribusi yang mungkin dapat diaplikasikan untuk sistem-sistem lain.
1.3 Tinjauan Pustaka Suatu sistem diskrit memuat sistem persamaan diferensi yang solusinya dapat ditentukan secara rekursif (Elaydi, 2006; Ogata, 1995). Teknik kendali seperti MPC klasik dapat digunakan untuk mendesain kendali pada sebuah sistem diskrit. Teknik kendali MPC adalah teknik kendali optimal yang menentukan kendali optimal dengan membuat state dan input prediksi sepanjang horison, kemudian mengaplikasikan elemen pertama pada barisan kendali optimal yang dihitung menggunakan teknik optimisasi (Camacho dan Bordons, 1999). Untuk kasus tanpa kendala, masalah optimisasi pada MPC dapat diselesaikan dengan menentukan gradien dari fungsi objektif, sedangkan untuk kasus berkendala dapat diselesaikan menggunakan pemrograman kuadratik (Maciejowski, 2001). Untuk menentukan kendali optimalnya, teknik kendali MPC klasik menggunakan beberapa konsep aljabar seperti derivatif dan faktorisasi matriks yang diambil dari Howard Anton (1997) dan Larson (2009). Masalah permainan dinamis kooperatif memuat beberapa pemain yang bekerjasama
untuk
menentukan
strategi
optimal
sedemikian
sehingga
meminimumkan fungsi objektif semua pemain. Solusi optimal lengkap untuk masalah ini belum tentu ada apabila fungsi objektif saling konflik, yaitu ketika suatu strategi menjadi lebih baik untuk satu pemain tetapi mengakibatkan hasil pemain lain menjadi lebih buruk. Masalah ini dapat diselesaikan dengan menentukan solusi Pareto yang dapat diperoleh dengan meminimumkan kombinasi linier konveks fungsi objektif semua pemain (Jacob Engwerda, 2005). Karena solusi Pareto tidak tunggal, selanjutnya dapat ditentukan solusi NashBargaining yang merupakan solusi Pareto yang ditentukan berdasarkan titik disagreement para pemain (Jacob Engwerda, 2005). Titik disagreement para
5
pemain dapat ditentukan menggunakan solusi minimax masing-masing pemain (Narahari, 2009; Myerson, 1997). Masalah kendali systemwide merupakan masalah kendali yang memuat beberapa subsistem yang saling berinteraksi. Dengan menganggap subsitem sebagai pemain, masalah ini dapat dipandang sebagai masalah permainan dinamis kooperatif yang dapat diselesaikan menggunakan teknik kendali MPC terdistribusi. Teknik kendali MPC terdistribusi menentukan solusi Pareto sebagai kendali optimalnya (Venkat, 2006). Karena solusi Pareto tidak tunggal, Valencia et al., (2011), menyelesaikan masalah MPC terdistribusi dengan menentukan solusi Nash-Bargaining sebagai kendali optimalnya. Teknik kendali MPC terdistribusi dapat diaplikasikan pada suatu systemwide waktu diskrit. Valencia et al., (2011), mengaplikasikan MPC terdistribusi pada tangki reaktor yang terdiri dari tiga subsistem.
Selain itu, Doan et al., (2010), mengaplikasikan MPC
terdistribusi pada kanal irigasi 4 subsistem secara serial (tidak bercabang). Teknik kendali MPC terdistribusi memuat masalah optimisasi fungsi konveks. Karena fungsi objektif yang digunakan merupakan fungsi konveks, secara umum, peminimalnya merupakan peminimal global (Bazaraa et al., 2006). Masalah optimisasi fungsi konveks dapat diselesaikan menggunakan beberapa metode. Menurut Steven Boyd dan Vandenberghe (2004), masalah optimisasi ini dapat diselesaikan secara efisien menggunakan metode interior point seperti algoritma active-set, linear search algorithms, atau gradient-based algorithms. Algoritma active-set menyelesaikan masalah optimisasi kuadratik secara numerik. Metode ini dinilai paling efisien karena waktu iterasinya secara umum lebih cepat (Fletcher, 2000; Nocedal dan Wright, 1999).
Dengan memanfaatkan fungsi
‘fmincon’ pada MATLAB yang menggunakan algoritma active-set, metode ini dipilih untuk menyelesaikan masalah optimisasi pada MPC terdistribusi. Teori permainan yang diberikan oleh Jacob Engwerda (2006) merupakan permainan kontinu. Untuk dapat diterapkan pada teknik kendali diskrit, penulis memberikan teori permainan diskrit dengan cara memodifikasi persamaan dinamik dan fungsi-fungsi objektifnya untuk waktu diskrit. Selain itu, penelitian ini membahas teknik kendali MPC terdistribusi, yaitu FC-MPC yang
6
menghasilkan solusi Pareto kemudian dilanjutkan NB-MPC yang menghasilkan solusi Nash-Bargaining, dimana pada kedua teknik kendali tersebut, tidak ada keterangan bagaimana cara menentukan nilai pembobotan fungsi objektif. Menurut teori permainan statis, nilai pembobotan fungsi objektif dapat ditentukan sesuai dengan kontribusi para pemain (Myerson, 1997). Pada penelitian ini, nilai pembobotan fungsi objektif ditentukan berdasarkan beban kendali masing-masing subsistem. Selanjutnya, teori ini diaplikasikan pada sistem kanal irigasi yang nilainilai parameternya diambil dari Overloop (2006). Doan (2010) memberikan model kanal irigasi 4 subsistem secara serial, sedangkan pada penelitian ini, penulis memodifikasi susunan kanal irigasi, yaitu menambah subsistem menjadi 5 subsistem dan menyusunnya secara bercabang. Untuk MPC biasa, tutorial simulasi MPC pada MATLAB dapat diperoleh di Wang (2009). Sedangkan untuk MPC terdistribusi, penulis membuat program MATLAB untuk mensimulasikan aplikasi teknik kendali MPC tedistribusi pada kanal irigasi 5 subsistem.
1.4 Metodologi Penelitian Metode penelitian yang digunakan adalah studi literatur dan simulasi. Penulis mempelajari teori-teori yang diperlukan untuk menyusun teknik kendali MPC terdistribusi meliputi teori matriks, teori optimisasi, dan teknik kendali MPC biasa. Selain itu, penulis juga memberikan teori permainan dinamis kooperatif diskrit yang diperlukan untuk menyusun teknik kendali MPC terdistribusi. Penelitian selanjutnya adalah menyusun teknik kendali MPC terdistribusi, yaitu FC-MPC yang menghasilkan solusi Pareto, kemudian dilanjutkan dengan NBMPC yang menghasilkan solusi Nash-Bargaining. Selanjutnya, teknik kendali MPC terdistribusi diaplikasikan pada kanal irigasi dan melakukan simulasi menggunakan program MATLAB. Teknik kendali MPC terdistribusi dibentuk dengan langkah-langkah sebagai berikut. Langkah pertama adalah memodelkan systemwide menggunakan model composite, yaitu model yang memuat persamaan state space masingmasing subsistem digabungkan dengan model interaksi antar subsistem. Selanjutnya, model composite digunakan untuk menyusun prediksi state dan
7
prediksi masukan. Teori permainan dinamis kooperatif diaplikasikan ke systemwide dengan menganggap subsistem sebagai pemain untuk menyusun optimisasi FC-MPC. Selanjutnya, disusun algoritma secara iteratif untuk menentukan kendali optimal berdasarkan optimisasi FC-MPC yang akan menghasilkan solusi Pareto. Selanjutnya disusun model negosiasi NB-MPC yang digunakan untuk menentukan solusi Nash-Bargaining. Masalah optimisasi yang muncul baik pada FC-MPC maupun pada NB-MPC dihitung menggunakan algoritma active-set yang termuat pada fungsi ‘fmincon’ pada MATLAB. Terakhir, teknik kendali MPC terdistribusi diaplikasikan pada kanal irigasi 5 subsistem yang memuat cabang didalamnya dan disimulasikan menggunakan program MATLAB.
1.5 Sistematika Penulisan Penulisan tesis ini terdiri dari 5 Bab, yaitu Bab I Pendahuluan, Bab II Landasan Teori, Bab III Kendali Model Prediktif Terdistribusi Berbasis Teori Permainan Dinamis Kooperatif, Bab IV Aplikasi MPC Terdistribusi Pada Kanal Irigasi, dan Bab V Penutup. Bab I memuat latar belakang dan permasalahan, tujuan dan manfaat penelitian, tinjauan pustaka, metode penelitian, dan sistematika penulisan. Bab II memuat landasan teori, yaitu teori matriks, kekontinuan fungsi bernilai skalar, himpunan dan fungsi konveks, pemrograman kuadratik, sistem diskrit, teknik kendali MPC, dan teori permainan dinamis kooperatif. Bab III berisi pembahasan untuk merumuskan teknik kendali MPC terdistribusi. Bab IV berisi cara membuat model composite dari sistem kanal irigasi sebagai systemwide yang terdiri dari 5 subsistem yang memuat interaksi, kemudian mengaplikasikan MPC terdistribusi untuk mengatur ketinggian air pada kanal irigasi dengan meminimalkan fungsi objektif yang diberikan, serta mensimulasikannya
menggunakan program
MATLAB.
Bab
kesimpulan dan saran yang dapat diambil dari hasil penelitian ini.
V
memuat
BAB II LANDASAN TEORI Untuk menyusun suatu teknik kendali, dibutuhkan teori-teori yang mendasarinya. Teori-teori yang digunakan untuk menyusun teknik kendali MPC terdistribusi meliputi derivatif dan faktorisasi matriks, kekontinuan fungsi bernilai skalar, optimisasi fungsi konveks, sistem diskrit, teknik kendali MPC, dan teori permainan dinamis kooperatif diskrit.
2.1 Derivatif dan Faktorisasi Matriks Teori matriks sangat berperan dalam teori kendali, terutama jika menggunakan metode analisis state space. Derivatif matriks akan digunakan untuk menentukan gradien dari fungsi objektif pada teknik kendali MPC. Sedangkan faktorisasi matriks akan digunakan untuk membuktikan kekonveksan fungsi objektif pada permainan dinamis kooperatif diskrit.
n
Teorema 2.1.1 Diberikan x Jika A
n n
(ii) Jika A
n n
(iii) Jika A
n n
(iv) Jika A
n m
(i)
maka maka
x x
,y
m
Ax
AT ,
xT Ax
simetris, maka , maka
x
xT Ay
, dan A adalah suatu matriks konstan,
xT A AT , xT Ax
x
Ay dan
2 xT A,
y
xT Ay
AT x.
Teorema (2.1.2) berikut menjelaskan tentang diagonalisasi ortogonal dari suatu matriks simetris. Diagonalisasi ini akan digunakan untuk menentukan akar kuadrat dari suatu matriks simetris semidefinit positif (atau definit positif).
8
9
Teorema 2.1.2 Jika matriks bernilai real M n
n
merupakan matriks simetris,
maka M dapat didiagonalisasi secara ortogonal menjadi
PT
M
1
0
0
2
0 0
0
0
n
(2.1.1)
P
dengan P adalah suatu matriks ortogonal i.e. PT P
PPT
I dan
1
, 2 ,...,
n
adalah nilai eigen dari M. Lebih lanjut, jika M semidefinit positif atau definit positif, maka dapat dibentuk matriks 1 1
M2
PT
0
0 0
0 2
0
0
(2.1.2)
P n
yang disebut sebagai akar kuadrat dari M.
Akar kuadrat dari suatu matriks simetris semidefinit positif (atau definit positif) dapat digunakan untuk membentuk faktorisasi matriks. Faktorisasi matriks ini akan digunakan untuk membuktikan kekonveksan fungsi objektif pada teori permainan dinamis kooperatif diskrit. Teorema (2.1.3) berikut memberikan faktorisasi matriks simetris semidefinit positif (atau definit positif).
Teorema 2.1.3 Diberikan matriks bernilai real M n n . Jika M n
n
simetris
semidefinit positif (definit positif), maka M dapat difaktorisasi menjadi M
M2
dengan M adalah suatu matriks simetris semidefinit positif (definit positif).
2.2 Kekontinuan Fungsi Bernilai Skalar Konsep kekontinuan fungsi bernilai skalar diperlukan untuk membuktikan kekontinuan
fungsi
Lyapunov
pada
subbab
(2.4).
Diberikan
fungsi
10
n
f : Df f ( x)
.
Nilai
f
di
x ( x1 , x2 ,..., xn ) D f
titik
adalah
f ( x1 , x2 ,..., xn ) . Definisi berikut menjelaskan persekitaran dari suatu titik.
n
Definisi 2.2.1 Persekitaran titik x0
dengan jari-jari
0 didefinisikan
sebagai N ( x0 )
n
x
| x
x0 n
Persekitaran terhapus titik x0 N d ( x0 )
.
dengan jari-jari
0 didefinisikan sebagai
N ( x0 ) { x0 }.
Definisi berikut memberikan beberapa konsep pada suatu himpunan, antara lain titik limit, titik batas, klosur, titik interior, himpunan terbuka, dan himpunan tertutup. Definisi-definisi ini akan digunakan untuk mendefinisikan konsep limit fungsi.
n
Definisi 2.2.2 Diberikan himpunan S n
(a) Titik x0
,
disebut titik limit dari S jika untuk sebarang bilangan real
0 berlaku N d ( x0 ) n
(b) Titik x0
S
.
disebut titik batas dari S jika untuk sebarang bilangan real
0 berlaku N ( x0 )
S
dan N ( x0 )
Sc
.
Himpunan semua titik batas dari S disebut batas dari S, dinotasikan sebagai S . (c) Klosur dari S, dinotasikan sebagai S , adalah gabungan S dengan semua titik batasnya, yaitu S (d) Titik x0
n
S
S.
disebut titik interior dari S jika terdapat
N ( x0 )
S.
0 sehingga
11
(e) Himpunan S dikatakan terbuka jika setiap anggotanya merupakan titik interior. Himpunan S dikatakan tertutup jika S c terbuka.
Definisi selanjutnya menjelaskan konsep limit fungsi. Konsep ini akan digunakan untuk mendefinisikan kekontinuan fungsi.
Definisi 2.2.3. Fungsi f : D f
n
dikatakan berlimit L untuk x
x0 ,
ditulis
lim f ( x)
x
L,
x0
jika x0 adalah titik limit dari Df dan untuk sebarang
0 terdapat
0
sedemikian sehingga
f ( x) L untuk semua x
N d ( x0 ) .
Df
Teorema (2.2.4) berikut menjelaskan konsep nilai limit dari penjumlahan dan perkalian dua fungsi. Teorema (2.4.4) dapat digunakan untuk membuktikan kekontinuan suatu fungsi.
Teorema 2.2.4 Diberikan fungsi f : D
n
dan g : D
n
, x0
adalah titik limit dari D, dan
lim f ( x)
x
L1 , lim g ( x) L2 ,
x0
x
x0
maka (a) lim ( f x
x0
g )( x)
(b) lim ( fg )( x) x
x0
L1 L2 ,
L1L2 .
Selanjutnya, diberikan definisi kekontinuan suatu fungsi bernilai skalar. Konsep kekontinuan fungsi akan digunakan untuk membuktikan Teorema kestabilan sistem menggunakan fungsi Lyapunov pada subbab (2.4).
12
n
Definisi 2.2.5 Diberikan fungsi f : D f
.
(a) Fungsi f dikatakan kontinu di titik x0
lim f ( x)
x
D f jika
f ( x0 ) .
x0
Secara lengkap : Fungsi f dikatakan kontinu di titik x0 sebarang
0 , terdapat
x Df
D f jika untuk
0 sedemikian sehingga
N ( x0 )
f ( x)
f ( x0 )
.
(b) Fungsi f dikatakan kontinu di Df jika f kontinu di setiap x
Df .
Teorema (2.2.6) berikut memberikan sifat kekontinuan fungsi bernilai skalar yang berbentuk kuadratik. Sifat ini akan digunakan untuk membuktikan Teorema kestabilan Lyapunov.
Teorema 2.2.6 Diberikan fungsi f : n
n
. Jika f(x) berbentuk kuadratik, i.e.
n
f ( x)
aij xi x j , aij
, maka f(x) adalah fungsi kontinu.
i 1 j 1
Definisi berikut menjelaskan konsep himpunan terhubung dan himpunan tidak terhubung, yang akan digunakan untuk membuktikan teorema nilai tengah suatu fungsi bernilai skalar.
n
Definisi 2.2.7 Suatu himpunan S dinyatakan sebagai S
A
, B
A
, A
dikatakan terhubung jika S tidak dapat
B dimana
B
,
A
B
.
Jika S dapat dinyatakan seperti di atas, maka S dikatakan takterhubung.
Definisi berikut menjelaskan konsep suatu region di ini akan digunakan untuk membuktikan teorema nilai tengah.
n
. Konsep region
13
Definisi 2.2.8 Suatu region S di
n
adalah gabungan dari himpunan terbuka
terhubung dengan beberapa, semua, atau tidak ada dari batasnya.
Teorema (2.2.9) disebut teorema nilai tengah. Teorema ini diperlukan untuk membuktikan teorema kestabilan sistem menggunakan fungsi Lyapunov pada subbab (2.4).
n
Teorema 2.2.9 (Teorema Nilai Tengah) Diberikan fungsi f : D f kontinu pada region S
f ( x1 ) u maka terdapat x3
Df . Jika x1 , x2
f ( x2 ) , S sehingga f ( x3 )
Bukti : Diketahui f ( x1 ) u dinyatakan sebagai S
R
S dan
R
f ( x2 ) . Andaikan tidak ada x3, maka S dapat
T dengan
x S | f ( x) u dan T
Untuk sebarang x0
u.
x S | f ( x) u .
, kekontinuan dari f mengakibatkan terdapat
sedemikian sehingga f ( x) u untuk semua x x0
T . Sehingga R
T
N ( x0 )
0
S . Ini berarti bahwa
. Dengan cara analog, diperoleh R
T
.
Sehingga S takterhubung. Terjadi kontradiksi dengan yang diketahui bahwa S adalah suatu region, yaitu S terhubung. Terbukti bahwa terdapat x3
f ( x3 )
S sehingga
u. Selanjutnya diberikan teorema Bolzano-Weirstrass pada
n
. Teorema ini
akan digunakan untuk membuktikan teorema kestabilan sistem menggunakan fungsi Lyapunov pada subbab (2.4). Teorema ini diambil dari R.G. Bartle (1964) yang memuat bukti dari teorema ini.
Teorema 2.2.10 (Teorema Bolzano-Weirstrass) Semua barisan terbatas di memiliki subbarisan yang konvergen.
n
14
2.3 Optimisasi Fungsi Konveks Teori kendali optimal memuat masalah optimisasi yang dapat diselesaikan menggunakan metode-metode optimisasi. Pada tesis ini, fungsi objektif yang digunakan pada teknik kendali MPC dan teori permainan dinamis kooperatif adalah fungsi konveks. Sehingga pada pembahasan ini, akan diberikan teori optimisasi untuk fungsi koneks.
2.3.1 Himpunan dan Fungsi Konveks Definisi dan teorema berikut menjelaskan tentang himpunan dan fungsi konveks yang diambil dari S. Boyd (2009) dan K.V. Mital (1983).
n
Definisi 2.3.1 Himpunan K
x1 , x2
K dan sebarang x1 (1
Bentuk
) x2
x1 (1
dikatakan konveks jika untuk sebarang
dengan 0
1 , berlaku
K.
(2.3.1)
) x2 pada definisi (2.3.1) disebut kombinasi linier
konveks dari titik x1 dan x2 . Secara geometris, himpunan K
2
konveks
garis
jika
untuk
sebarang
dua
titik
di
K,
segmen
dikatakan yang
menghubungkannya termuat di K. Contoh himpunan konveks dan himpunan tidak konveks di ruang berdimensi dua dapat diilustrasikan pada Gambar (2.1).
K2 K1
(a)
(b)
Gambar 2.1 (a) Himpunan konveks K1, dan (b) Himpunan tidak konveks K2 Pada Gambar (2.1), himpunan yang digambarkan oleh (a) merupakan himpunan konveks karena segmen garis yang menghubungkan dua titik sebarang
15
di K1 termuat di K1. Sedangkan K2 tidak konveks karena terdapat segmen garis yang menghubungkan dua titik di K2 tidak termuat di K2. Berdasarkan definisi (2.3.1), himpunan
n
,n
merupakan himpunan konveks.
Teorema (2.3.2) berikut disebut teorema separasi. Teorema ini akan digunakan untuk membuktikan lemma (2.6.4) pada subbab (2.6). Untuk ruang berdimensi dua, teorema separasi dapat diilustrasikan oleh Gambar (2.2).
X
Gambar 2.2 Ilustrasi teorema separasi untuk
2
.
Teorema 2.3.2 Diketahui himpunan tak kosong, tertutup, terbatas, dan konveks
X
n
n
. Jika diberikan x0
sehingga xo
sedemikian sehingga untuk setiap x
X , maka terdapat p
X berlaku p T x
n
,p
0
p T x0 .
Bukti : Misalkan w merupakan titik di X yang terdekat dengan xo, sehingga untuk sebarang x
X
berlaku
w x0
x x0
dan
y
x (1
termuat di X. Sehingga
w x0
y x0
w x0 w x0
2
( w x0 )T (w x0 )
x (1
)w x0
x (1
) w x0
2
2
( x w)T ( x w) (w x0 )T ( w x0 )
2 ( w x0 )T ( x w) 2
Untuk
( x w)T ( x
w)
2 (w x0 )T ( x w).
0 , kedua ruas dibagi dengan
, diperoleh
)w, 0
1
16
( x w)T ( x w)
2(w x0 )T ( x w)
( x w)T ( x w)
2 pT ( x w),
dengan p
( w x0 ) merupakan vektor konstan karena x0 tertentu di
tertentu di X. Selanjutnya, untuk
0
n
dan w
0 , diperoleh
2 p( x w)
T
p ( x w) 0
pT x
pT w 0 pT x
Tetapi karena p T ( w
pT w. x0 )
pT p
0 maka p T w
p T x0 dan p T x
pT w
p T x0 .
Selanjutnya, apabila x0 mendekati w melalui titik-titik sepanjang garis yang menghubungkan x0 dan w, maka untuk x0 pT x
w , berlaku
p T x0 ,
untuk setiap x
X.
Konsep kekonveksan juga dapat diaplikasikan pada fungsi. Definisi berikut menjelaskan kekonveksan dan kekonkafan suatu fungsi. Definisi ini diambil dari K.V. Mital (1983).
Definisi 2.3.3 Diberikan fungsi f :
n
dan himpunan konveks K
Fungsi f ( x) dikatakan konveks di K jika untuk sebarang dua titik x1 , x2 sebarang
f
dengan 0
x1 (1
) x2
n
.
K dan
1 berlaku
f ( x1 ) (1
) f ( x2 ).
(2.3.2)
Fungsi f dikatakan konveks tegas jika tanda pertidaksamaan pada (2.3.2) adalah tegas. Fungsi f disebut konkaf (konkaf tegas) jika fungsi
f konveks (konveks
tegas).
Secara geometris, fungsi konveks untuk satu variabel dapat di ilustrasikan pada Gambar (2.3) berikut.
17
f(x)
f(x1)+(1
)f(x2)
f( x2+(1
)x2)
(x2, f(x2))
(x1, f(x1))
x1
x
x2
Gambar 2.3 Contoh grafik fungsi konveks satu variabel
Teorema (2.3.4) berikut menjelaskan tentang kekonveksan suatu fungsi yang dibentuk dari kombinasi linier beberapa fungsi konveks. Kekonveksan fungsi yang berbentuk seperti ini akan digunakan pada teori permainan dinamis kooperatif diskrit.
n
Teorema 2.3.4 Diberikan fungsi f i : i
dan sebarang
0 dan i 1,2,..., m . Jika himpunan fungsi
konveks di K
n
, maka fungsi
f
1 1
f
f1 , f 2 ,..., f m f
2 2
m m
dengan
i
masing-masing
juga konveks di K.
m
Untuk kasus khusus Teorema (2.3.4) dimana
m
1 , maka kombinasi
i 1 m
linier konveks dari
f1 , f 2 ,..., f m
yaitu
f juga konveks. Teorema (2.3.5)
i i i 1
berikut menjelaskan kekonveksan suatu fungsi yang memuat logaritma.
Teorema 2.3.5 Jika diberikan fungsi f : himpunan konveks dan f ( x) fungsi konveks di K.
n
0 untuk semua x
konkaf di K K , maka
n
dengan K
log f ( x) adalah
18
2.3.2 Fungsi Berbentuk Kuadratik Pada teknik kendali MPC dan teori permainan dinamis kooperatif diskrit, fungsi objektif yang akan digunakan adalah fungsi berbentuk kuadratik. Selanjutnya, diberikan definisi fungsi yang berbentuk kuadratik yang diambil dari K.V. Mital (1983).
n
Definisi 2.3.6 Diberikan fungsi f :
. Fungsi berbentuk kuadratik
f ( x)
didefinisikan sebagai f ( x)
c11 x12
c22 x22
c1n x1 xn
dengan x
cnn xn2
c23 x2 x3
( x1 , x2 ,..., xn )T
Misalkan aij
n
a ji
c12 x1 x2
c13 x1 x3
cn 1 xn 1 xn
dan cij
1 cij , i 2
(2.3.3)
, i, j 1,2,..., n .
j . Kemudian aij
a ji disubstisusikan ke
fungsi kuadratik (2.3.3) menggantikan cij , dan aii disubsitusikan menggantikan
cii , maka fungsi kuadratik (2.3.3) dapat ditulis sebagai f ( x)
a11 x12
a22 x22
2a1n x1 xn n
ann xn2
2a23 x2 x3
2a12 x1 x2 2a( n
2a13 x1 x3
x x
1) n n 1 n
n
aij xi x j i 1 j 1
x1
T
a11
a12
a1n
x1
x2
a21
a22
a2 n
x2
xn
an1
an 2
ann
xn
T
x Ax. Karena aij
a ji maka matriks A adalah matriks simetris. Definisi berikut
menjelaskan kedefinitan fungsi yang berbentuk kuadratik.
19
n
Definisi 2.3.7 Diberikan fungsi f :
xT Ax dan A adalah
dengan f ( x)
matriks real simetris berukuran n n . Fungsi f dikatakan (i)
Definit positif jika xT Ax
0 untuk semua x
(ii) Semi definit positif jika xT Ax sedikit sebuah vektor x (iii) Definit negatif jika xT Ax
0
0 untuk semua x
0 sehingga xT Ax
0 untuk semua x
(iv) Semi definit negatif jika xT Ax paling sedikit sebuah vektor x
0 dan terdapat paling
0 0
0 untuk semua x 0 sehingga xT Ax
0 dan terdapat
0
Teorema (2.3.8) berikut dapat digunakan untuk menyelidiki kedefinitan suatu fungsi berbentuk kuadratik menggunakan nilai eigen. Teorema ini diambil dari K.V. Mital (1983).
Teorema 2.3.8 Misalkan nilai-nilai eigen dari matriks real simetris A berukuran n n adalah
i
, i 1, 2,..., n . Maka fungsi berbentuk kuadratik f ( x)
xT Ax
adalah (i) Definit positif jika dan hanya jika
0 untuk semua i,
i
(ii) Semi definit positif jika dan hanya jika
i
0 dan
j
0 untuk paling
sedikit satu j, (iii) Definit negatif jika dan hanya jika
i
0 untuk semua i,
(iv) Semi definit negatif jika dan hanya jika
i
0 dan
j
0 untuk paling
sedikit satu j.
Teorema (2.3.9) berikut menjelaskan kekonveksan dari suatu fungsi berbentuk kuadratik semidefinit positif (atau definit positif). Sifat kekonveksan fungsi yang berbentuk kuadratik akan digunakan pada teknik kendali MPC dan teori permainan dinamis kooperatif diskrit.
20
n
Teorema 2.3.9 Diketahui x kuadratik f ( x)
n
dan fungsi f :
. Jika fungsi berbentuk
xT Ax adalah semidefinit positif (definit positif), maka f ( x) n
konveks (konveks tegas) di
.
Bukti : Hanya dibuktikan untuk sifat konveks. Diketahui xT Ax n
sebarang dua titik x1 , x2
f ( x1 ) (1
) f ( x2 )
0
f ( x1 ) (1
x1 (1
) x2 dengan 0
) x2T Ax2
x1T Ax1 (1
) x2T Ax2
x1 (1
x1T Ax1 (1
) x2T Ax2
x1T
x1T Ax1 (1
) x2T Ax2
x1T Ax1
xT Ax
2
1
) x1T Ax1
(1
(1
) x1T Ax1
x2T Ax2
(1
)( x1
) f ( x2 )
dan
x1
x2
f ( x) 0 n
yang artinya f(x) konveks di
(1
T
(1
A
x1 (1
) x2T
)2 x2T Ax2
) x2T Ax2
x2 )T A( x1
) x2
x1T Ax1 2 (1
(1
(1
1
1 berlaku
f ( x)
x1T Ax1 (1
2
Karena
dan x
0. Untuk
) x2
Ax1 (1 ) x1T Ax2
2 (1
) Ax2
(1
) 2 x2T Ax2
) x1T Ax2
).2 x1T Ax2
2 x1T Ax2
x2 ) 0. adalah
sebarang
maka didapat
vektor
f ( x)
di
f ( x1 ) (1
n
,
dan
) f ( x2 )
.
2.3.3 Optimisasi Fungsi Konveks Tanpa Kendala Konsep utama dalam teori optimisasi adalah nilai minimum/maksimum dari suatu fungsi. Definisi berikut memberikan pengertian dari nilai minimum suatu fungsi. Masalah optimisasi yang akan dibahas adalah masalah minimisasi, karena masalah maksimisasi dapat dibawa ke masalah minimisasi dengan cara mengalikan fungsi objektifnya dengan –1.
Definisi 2.3.10 Fungsi f : titik x0
n
n
jika untuk semua x
dikatakan memiliki nilai minimum global di n
berlaku f ( x0 )
f ( x ). .
21
Definisi berikut memberikan pengertian peminimum lokal dari suatu fungsi. Konsep minimum lokal akan digunakan pada teorema syarat cukup suatu fungsi untuk mencapai peminimumnya.
Definisi 2.3.11 Fungsi f : titik
f ( x0 )
x0
n
jika
n
terdapat
dikatakan memiliki nilai minimum lokal di persekitaran
N ( x0 )
sedemikian
sehingga
f ( x ) untuk semua x N ( x0 ) .
Teorema (2.3.12) berikut memberikan syarat cukup suatu fungsi untuk mencapai peminimum dan dapat digunakan untuk menentukan titik minimum suatu masalah minimisasi. Teorema ini diambil dari K.V. Mital (1983).
Teorema 2.3.12 Diberikan fungsi f : di titik x0 . Jika
f ( x0 )
n
diferensiabel sampai tingkat dua
0 dan hessian H ( x0 ) definit positif, maka x0 adalah
peminimum lokal tegas.
Karena pada tesis ini, fungsi yang dibahas adalah fungsi konveks, selanjutnya diberikan konsep minimum lokal dan minimum global pada fungsi konveks. Teorema (2.3.13) berikut menjelaskan tentang sifat-sifat peminimum lokal pada fungsi konveks.
Teorema 2.3.13 Diberikan fungsi f : K
n
, himpunan K konveks, dan
fungsi f(x) konveks di K. Jika f ( x) memiliki titik peminimum lokal, maka titik peminimum tersebut tersebut adalah titik peminimum global. Jika f ( x) memiliki titik peminimum lokal di beberapa titik, maka titik peminimum global terjadi di setiap kombinasi linier konveks dari titik-titik tersebut.
22
Berdasarkan Teorema (2.3.13), peminimum lokal suatu fungsi konveks adalah peminimum global. Sehingga, dari Teorema (2.3.13), dapat diberikan syarat cukup suatu fungsi konveks mencapai minimum global.
n
Akibat 2.3.14 Diberikan fungsi f : K
, himpunan K konveks, fungsi
f ( x) konveks di K dan diferensiabel sampai tingkat dua di titik x0 K . Jika
f ( x0 )
0 dan hessian H ( x0 ) definit positif, maka x0 adalah peminimum
global.
Teorema (2.3.15) berikut memberikan sifat ekuivalensi dua buah masalah optimisasi. Sifat ini akan digunakan pada subbab (2.6) untuk menyelesaikan masalah optimisasi Nash-Bargaining.
Teorema 2.3.15 Diberikan fungsi fi : K i fungsi fi ( x) konveks di Ki, dan di
n
, himpunan K i konveks,
adalah suatu konstanta. Maka
M
x
arg max x K
di
f i ( x)
i 1
kendala : f i ( x)
di , i 1, 2,..., M
jika dan hanya jika M
x
di
arg max log x K
fi ( x)
i
i 1
kendala : f i ( x)
di , i 1, 2,..., M
untuk suatu bilangan real
Bukti : Karena f i ( x)
i
0, i 1,2,..., M .
di , i 1,2,..., M , maka
Diketahui M
x
arg max x K
di
f i ( x)
i 1
kendala : f i ( x)
di , i 1, 2,..., M
di
fi ( x)
0, i 1,2,..., M .
23
maka M
M
fi ( x )
di f1 ( x )
d1
i 1
fM (x )
dM
d1
f1 ( x )
1
dM
fM ( x )
M
log d1
f1 ( x )
1
dM
fM ( x )
M
d1
f1 ( x)
d1
f1 ( x)
log d1
M
dM 1
f1 ( x)
dM 1
f M ( x) f M ( x) dM
M
f M ( x)
M
di
log
i
fi ( x )
log
di
fi ( x)
i
i 1
i 1
untuk suatu bilangan real
fi ( x)
di
i 1
0, i 1, 2,..., M . Sehingga
i
M
x
di
arg max log x K
fi ( x)
i
.
i 1
kendala : f i ( x)
di , i 1, 2,..., M .
2.3.4 Pemrograman Kuadratik Berkendala Masalah
kendali
MPC
berkendala
dapat
diselesaikan
dengan
membentuknya menjadi masalah pemrograman kuadratik berkendala. Diberikan fungsi V : oleh V ( )
n
n
dengan V ( ),
1 2
T
T
berbentuk kuadratik yang didefinisikan
. Diasumsikan matriks
semidefinit positif (atau
definit positif), sehingga menurut Teorema (2.3.9), V ( ) adalah konveks (atau n
konveks tegas) di
. Masalah pemrograman kuadratik dapat dinyatakan dalam
bentuk umum
minn
V( )
kendala : H dengan n n,
n
,
1 2
T
T
(2.3.4)
,
,
(2.3.5)
h,
(2.3.6)
merupakan matriks real simetris definit positif berukuran
merupakan matriks real berukuran m n , dan H merupakan matriks real
berukuran p n , dan vektor-vektor
, h, dan
secara berturut-turut merupakan
M
24
vektor real berdimensi n, m dan p. Didefinisikan fungsi lagrangian yang memuat pengali lagrange
dan
1 2
L1 ( , , , t )
1 2
serta vektor slack t yaitu T
T
t
T
T
t
H H
Syarat perlu dan syarat cukup untuk
h t
h.
(2.3.7)
mencapai minimum global
diberikan oleh kondisi Karush-Kuhn-Tucker (KKT) (S. Boyd, 2004). Kondisi KKT yaitu harus terdapat vektor pengali lagrange t
0 dan
, dan vektor slack
0 sedemikian sehingga
L1 ( , , )
0
L1 ( , , )
0
L1 ( , , )
0
T
t
T
t H
HT
0,
0,
h 0,
0,
atau dapat ditulis kembali menjadi
HT
T
H
, h,
(2.3.9) (2.3.10)
t
tT
(2.3.8)
0
(2.3.11)
Kondisi KKT (2.3.8)-(2.3.11) dapat digunakan untuk menyelesaikan masalah pemrograman kuadratik (2.3.4)-(2.3.6) secara analitik. Karena teknik kendali MPC terdistribusi akan diselesaikan secara numerik, maka akan digunakan algoritma active-set untuk menyelesaikan pemrograman kuadratik berkendala secara numerik. Metode active-set mengasumsikan bahwa terdapat solusi fisibel, yaitu solusi yang memenuhi kendala (2.3.5)-(2.3.6). Suatu solusi fisibel pasti memenuhi kendala persamaan (2.3.6) dan himpunan bagian dari pertidaksamaan (2.3.5). Himpunan bagian ini dikatakan aktif (active-set), dimana kendala pertidaksamaan (2.3.5) menjadi kendala persamaan. Misalkan a menotasikan elemen-elemen baris
25
dari pertidaksamaan (2.3.5) yang merupakan anggota active set, sehingga a
a
. Misalkan
active-set
menyatakan solusi fisibel pada saat iterasi ke-r, metode
r
berjalan
memperbaiki
solusi
dengan
arah
yang
r
meminimumkan fungsi objektif (2.3.4) ketika memenuhi kendala persamaan H
h dan
a r
a
set). Jika solusi
fisibel, yaitu memenuhi
r
dilanjutkan dengan
, tanpa memperhatikan kendala pertidaksamaan (inactive-
r 1
r
pencarian dengan arah
. Jika solusi
, maka iterasi
r
tidak fisibel, maka dibuat garis
r
untuk menentukan lokasi titik yang menghilangkan
fisibilitas, yaitu salah satu titik yang tidak aktif menjadi aktif. Solusi pada titik ini digunakan untuk iterasi berikutnya, yaitu
r 1
r
, 0
r
r
1 dan
kendala aktif yang baru ditambahkan ke active-set. Selanjutnya adalah menentukan kapan suatu solusi mencapai optimum global, atau kapan perbaikan solusi dapat dilakukan. Untuk mengevaluasinya, digunakan kondisi KKT (2.3.8)-(2.3.11) untuk setiap solusi
r 1
. Karena kesemua
elemen t adalah positif, maka kondisi komplementer dari (2.3.11) yang mengakibatkan kesemua elemen dari inactive
adalah
nol.
Sehingga
yang berkorespondensi dengan kendala hanya
elemen-elemen
dari
yang
berkorespondensi dengan kendala aktif yang perlu di evaluasi. Jika kesemua bernilai nonnegative ( dilanjutkan. Diberikan
0 ), maka solusi global telah didapat, selain itu iterasi q
0 , karena
q
merupakan pengali lagrange yang
berkorespondensi dengan kendala pertidaksamaan ke-q, yang merupakan kendala aktif, nilai negatif tersebut menandakan bahwa nilai fungsi objektif dapat direduksi dengan membuat kendala menjadi inactive, yaitu berpindah ke solusi dimana
q
q
. Sehingga kendala ke-q dihilangkan dari active-set. Jika ada
lebih dari satu elemen
yang negatif, maka kendala yang dihilangkan adalah
kendala yang berkorespondensi dengan elemen
yang paling kecil. Sehingga
active set baru terpilih, dan kesemua proses diatas diulang dengan cara mengganti r
dengan
r 1
. Misalkan V ( r )
1 2
T r
T r
r
. Karena V
r 1
V
r
26
maka metode ini menjamin bahwa proses iterasi dapat dihentikan saat optimum global karena kekonveksan dari fungsi objektif. Selanjutnya, akan diberikan masalah optimisasi pada saat iterasi ke-r secara lebih detail. Fungsi objektif pada saat iterasi ke-r adalah
V
1 2
r
1 2 1 2 1 2 1 2
V
T
T
r
r
T
T
r
r
T
(2.3.12)
r
r
T
T
r
r
T r
r
2
T
T
r
r
1 2
r
T
T
T
r
T
r
r
r
1 2
T
T
T
T
T
.
r
Karena untuk setiap iterasi r nilai V
T r
r
(2.3.13)
merupakan konstanta, maka suku
tersebut dapat dihilangkan dari (2.2.13). Misalkan
r
r
, maka masalah
pemrograman kuadratik diatas dapat dinyatakan sebagai
minn
1 2
T
T
0 dan
kendala : H
(2.3.14)
,
r
0.
a
(2.3.15)
Masalah optimisasi (2.3.14)-(2.3.15) merupakan masalah optimisasi konveks kuadratik, tetapi hanya memuat kendala persamaan. Masalah optimisasi tersebut dapat diselesaikan dengan membentuk fungsi lagrange, kemudian menggunakan kondisi KKT. Karena masalah (2.3.14)-(2.3.15) hanya memuat kendala persamaan, maka hanya digunakan kondisi KKT (2.3.8)-(2.3.9). Fungsi lagrange masalah optimisasi (2.3.14)-(2.3.15) dengan pengali lagrange
L2 (
dan
,
, didefinisikan oleh
,
)
1 2
T
T r
H
a
.
27
Kondisi KKT untuk masalah diatas yaitu harus terdapat vektor pengali langrange dan
sedemikian sehingga HT
L2 (
,
,
)
0
L2 (
,
,
) 0
L2 (
,
,
) 0
H
T a
r
,
0, 0.
a
Persamaan-persamaan diatas dapat ditulis menjadi sistem persamaan linier
H a
HT 0
T
0
0
0
0
0
a
r
yang dapat diselesaikan menggunakan eliminasi Gaussian. Panjang langkah saat iterasi ke-r adalah a
ditentukan dengan cara sebagai berikut.
r
r
r
a r
a
r
a
r
a r a r a
a r
min 1,
Sehingga
a
,
a
r
0.
r
,
a
r
0 .
r
Berdasarkan uraian diatas, dapat disusun metode active-set yang diberikan oleh algoritma berikut. Diberikan r = 0, solusi fisibel awal Untuk r
0
dan active set
0
0,1, 2,... berjalan iterasi berikut.
(a) Dibentuk
r
r
(b) Selesaikan (2.3.14)-(2.3.15), yaitu dengan menyelesaikan
H a
HT 0
T
0
0
0
0
0
a
untuk memperoleh (c) Jika
r
r
r
(2.3.16)
.
0 , langsung ke langkah (f). Jika
r
0 ke langkah (d).
.
28
(d) Tentukan panjang langkah a r
min 1,
r
a
Jika
,
a
r
1,
r
{a} . Selain itu,
r 1
r
0
r
tambahkan
kendala r 1
r
ke-a
ke
active-set,
sehingga
.
(e) Update r 1
r
r
r
r
r 1, kembali ke langkah (a).
(f) Cek nilai
. Jika
q
*
didapat, yaitu
r
0 untuk semua q, maka solusi global telah
, stop iterasi. Selain itu, untuk
0 , kendala ke-q
q
yang paling kecil dihilangkan dari active-set, sehingga
r 1
r
\ {q}
kemudian kembali ke langkah (b).
Flowchart dari algoritma active-set dapat dilihat pada lampiran (4). Akan ditunjukkan bahwa algoritma active-set konvergen, yaitu mencapai titik minimum untuk iterasi yang berhingga. Jika solusi dari (2.3.16) adalah saat iterasi ini, yaitu active-set
r
, yaitu
r
r
r
0 maka titik
adalah peminimum global fungsi objektif (2.3.4) untuk peminimum global V ( ) . Jika solusi
r
bukan solusi dari
masalah asli (2.3.4)-(2.3.6), (yaitu terdapat pengali lagrange yang bernilai negatif),
maka
pengambilan
mengakibatkan, V
r 1
V
langkah r
r 1
r
r
r
untuk
k
0
yang artinya V ( ) bergerak turun terhadap
iterasi. Selanjutnya, karena iterasi-iterasi berikutnya memiliki nilai V yang lebih kecil dari peminimum global untuk active-set ini tidak akan pernah kembali ke active-set
r
r
, maka terjamin bahwa algoritma
. Karena banyaknya kemungkinan
dari active-set berhingga, maka algoritma tidak akan beriterasi selamanya. Jadi, algoritma ini menjamin bahwa iterasi berakhir pada solusi peminimum V ( ) .
29
2.4 Sistem Diskrit LTI Teknik kendali MPC terdistribusi merupakan teknik kendali untuk sistem diskrit. Untuk menyelidiki kestabilan sistem hasil kendali MPC terdisribusi, digunakan teori kestabilan Lyapunov. Pada tesis ini, sistem yang dibahas adalah sistem diskrit autonomous, dimana persamaan diferensinya tidak memuat variabelnya secara eksplisit. Pada pembahasan ini, variabel dari sistem adalah waktu, sehingga persamaan diferensi yang dibahas tidak bergantung waktu, yang disebut dengan sistem time invariant. Pada subbab ini, akan diberikan teori sistem diskrit khusus untuk sistem diskrit Linear Time Invariant (LTI), konsep kestabilan Lyapunov, sifat keterkendalian sistem, dan sifat keterstabilan sistem. Didefinisikan k menyatakan waktu instant. Model matematika dari sistem diskrit LTI dapat dinyatakan sebagai sistem persamaan state space x ( k 1) y (k )
Ax ( k )
Bu ( k ),
x (0)
(2.4.1a)
x0
(2.4.1b)
Cx ( k ) Du (k )
dengan x(k) menyatakan vektor state berdimensi n, u(k) menyatakan vektor masukan berdimensi p, dan y(k) menyatakan vektor keluaran berdimensi q. Matriks-matriks A, B, C, dan D secara berturut-turut merupakan matriks-matriks konstan berukuran n×n, n×p, q×n, dan q×p. Solusi dari persamaan diferensi (2.4.1a) untuk sebarang k
dapat
ditentukan menggunakan prosedur iteratif berikut. x (1)
Ax (0) Bu (0) Ax0 Bu (0),
x(2)
Ax(1) Bu (1) A Ax0 A2 x0
x(3)
Bu (0)
Bu (1)
ABu (0) Bu (1),
Ax(2)
Bu (2)
A3 x0
A2 Bu (0)
A A2 x0 ABu (1)
ABu (0) Bu (1)
Bu (2)
Bu (2),
k 1
x (k )
Ak x0
Ak i 0
i 1
Bu (i ).
(2.4.2a)
30
Solusi x(k) disubstitusikan ke persamaan (2.3.1b) diperoleh solusi k 1
y(k ) CAk x0
CAk
i 1
Bu (i )
Du(k )
(2.4.2b)
i 0
2.4.1 Kestabilan Lyapunov Sistem Diskrit LTI Untuk menganalisis kestabilan sistem, ada dua metode yang diberikan oleh Lyapunov, yaitu metode Lyapunov pertama (indirect method) dan metode Lyapunov kedua (direct method). Metode Lyapunov pertama menganalisis kestabilan sistem menggunakan solusi eksplisit dari sistem, sedangkan metode Lyapunov kedua tidak menggunakan solusi sistem, tetapi menggunakan fungsi Lyapunov. Pada tesis ini, akan diberikan metode Lyapunov pertama, yaitu menyelidiki kestabilan sistem menggunakan nilai eigen. Konsep ini juga digunakan untuk menyelidiki sifat keterstabilan sistem. Sedangkan metode Lyapunov kedua akan digunakan untuk menyelidiki kestabilan sistem hasil kendali MPC terdistribusi. Diberikan sistem diskrit LTI (2.4.1). Untuk sistem tanpa kendali, sistem (2.4.1) dapat dituliskan sebagai
x(k 1)
Ax(k ) , x(0)
x0 .
(2.4.3)
Definisi dan teorema berikut menjelaskan tentang titik ekuilibrium dan sifat kestabilan dalam arti Lyapunov yang diambil dari S. Elaydi (2005).
Definisi 2.4.1 Diberikan sistem (2.4.3). State xe
n
yang memenuhi Axe
disebut state ekuilibrium.
Didefinisikan state baru z (k ) z ( k 1)
x ( k 1) Ax( k )
n
dengan z (k )
xe xe
A z (k )
xe
xe
z ( k 1)
Az ( k )
Axe
xe
z ( k 1)
g z (k )
x(k ) xe , maka
xe
31
dengan g 0
0 jika Axe
xe . Sehingga dapat selalu diasumsikan bahwa state
ekuilibrium sistem (2.4.3) adalah xe
0 . Untuk selanjutnya, asumsi ini dipakai.
Definisi (2.4.2) berikut memberikan definisi kestabilan dari sistem diskrit linier LTI (2.4.3).
Definisi 2.4.2 Misalkan x(k) menyatakan solusi sistem (2.4.3) saat k i. State ekuilibrium xe real
0 dikatakan stabil jika untuk sebarang bilangan
0 terdapat bilangan real
x( k )
x0 ii. State ekuilibrium xe
lim x(k )
k
iii. Titik ekulibrium xe
,
( )
0 sehingga
untuk semua k
.
0 dikatakan stabil asimtotik jika xe stabil dan
0, 0 dikatakan tidak stabil jika tidak memenuhi (i)
Secara geometris, kestabilan dalam arti Lyapunov untuk sistem single state dapat diilustrasikan pada Gambar (2.4.1).
(a) xe
(a) (b) (c) Gambar 2.4.1 Ilustrasi kestabilan sistem single state 0 stabil, (b) xe 0 stabil asimtotik, dan (c) xe 0 tidak stabil
Misalkan z mendefinisikan modulus dari bilangan kompleks z yaitu z
a bi ,
a 2 b 2 . Untuk metode Lyapunov pertama, Teorema (2.4.3) berikut
dapat digunakan untuk menyelediki kestabilan sistem menggunakan nilai eigen dari matriks sistem, yang diambil dari Olsder (2004). Konsep kestabilan nilai eigen ini akan digunakan untuk mendefinisikan sifat keterstabilan sistem.
32
Teorema 2.4.3 Diberikan sistem (2.4.3) dan
1
, 2 ,..., k , k
n adalah nilai-nilai
eigen berbeda dari matriks A, maka (a) Titik ekuilibrium xe adalah stabil asimtotik jika dan hanya jika i
1,
i 1, 2,..., k .
(b) Titik ekuilibrium xe adalah stabil jika dengan
j
i
1, i 1, 2,..., k dan tepat satu j
1 serta banyaknya multisiplitas
j
sama dengan banyaknya
vektor eigen bebas linier yang bersesuaian dengan
j
.
(c) Untuk kasus lain, xe adalah tidak stabil.
Teorema (2.4.3) dapat digunakan untuk menyelidiki kestabilan sistem tanpa kendali (2.4.3). Selanjutnya akan diberikan cara untuk menyelidiki sistem dengan kendali (2.4.1) menggunakan teori kestabilan fungsi Lyapunov. Menurut teori mekanik klasik, suatu sistem “getaran” bersifat stabil jika energi totalnya menurun secara kontinu sampai state ekuilibriumnya tercapai. Metode Lyapunov kedua merupakan perumuman dari sifat tersebut. Pada kenyataannya, untuk sistem matematis murni, fungsi energi sulit untuk didefinisikan. Lyapunov mengenalkan fungsi Lyapunov yang merupakan fungsi energi tiruan. Untuk sistem kontinu, fungsi Lyapunov adalah fungsi kontinu bernilai skalar yang definit positif, derivatif parsial pertama (terhadap semua argumennya) kontinu di daerah sekitar titik asal, dan derivatif pertama terhadap waktu sepanjang trayektori adalah definit negatif (atau semidefinit negatif). Untuk sistem diskrit, fungsi Lyapunov tidak menggunakan derivatifnya, tetapi menggunakan variasi terhadap state sistem. Fungsi Lyapunov dapat bergantung terhadap waktu atau tidak. Pada tesis ini, khusus dibahas fungsi Lyapunov yang tidak bergantung terhadap waktu.
33
Metode Lyapunov kedua diberikan oleh definisi dan teorema berikut diambil dari S. Elaydi (2005), yang dapat digunakan untuk menyelidiki kestabilan suatu sistem menggunakan fungsi Lyapunov. Misalkan V :
n
menyatakan fungsi bernilai skalar. Variasi dari V
terhadap sistem (2.4.1) dapat didefinisikan oleh
V ( x) V x(k 1) Jika
V ( x)
V x(k ) .
(2.4.4)
0 , maka V tidak naik sepanjang solusi dari (2.4.1) yang diberikan
oleh (2.4.2).
Definisi 2.4.4 Fungsi V :
n
dikatakan fungsi Lyapunov untuk sistem (2.4.1)
jika 1. V kontinu 2.
V ( x)
Misalkan Br ekuilibrium xe
n
0 untuk semua x
x
n
: x
r
menyatakan bola terbuka dari state
0 dengan jari-jari r. Selanjutnya diberikan Teorema (2.4.5) yang
memberikan syarat cukup kestabilan suatu sistem menggunakan fungsi Lyapunov. Teorema ini akan digunakan untuk menyelidiki kestabilan sistem hasil kendali MPC terdistribusi pada Bab III.
Teorema 2.4.5 Misalkan V adalah fungsi Lyapunov untuk sistem (2.4.1), 1. Jika V definit positif maka state ekuilibrium xe 2. Jika xe
0 stabil dan
V ( x)
maka state ekuilibrium xe
0 adalah stabil.
0 untuk semua x
n
1
xe ,
0 adalah stabil asimtotik.
Bukti : Diberikan sistem (2.4.1) dengan titik ekuilibrium xe Lyapunov V. Pilih sebarang
dengan x
0 sehingga B
n 1
0 dan fungsi
. Karena (2.4.1a) kontinu,
34
yaitu memenuhi definisi (2.2.5), maka terdapat x
B
2
B 1 . Selanjutnya, diambil 0
maka V ( x ) ( )
0 sedemikian sehingga jika
2
min V ( x) |
x
1
2
. Didefinisikan
.
(2.4.5)
Menurut teorema nilai tengah (Teorema 2.2.9), terdapat 0
( ) untuk semua x dengan x
sehingga V ( x)
sedemikian
.
Selanjutnya akan dibuktikan dengan kontradiksi, bahwa
x0
B
x (k ) B untuk semua k
0.
(2.4.6)
Andaikan (2.4.6) tidak terjadi, artinya terdapat x0 positif m sehingga x ( m)
B 2 , yang berakibat x ( m 1)
B
V x(m 1)
( ).
Kontradiksi
V x(m 1)
V x(m)
... V ( x0 )
bilangan real
dengan
x(k )
Untuk sebarang x0
n
sehingga berlaku x (k )
( )
, karena xe
konvergen ke nol. Karena xe
0 sehingga
0 stabil asimtotik, yaitu n
0, untuk semua x0
B
bahwa
0 stabil.
Selanjutnya akan dibuktikan bahwa xe
lim x(k )
diketahui
untuk semua k
yang artinya, titik ekuilibrium xe
k
yang
B 1 , sehingga
( ) . Jadi, terbukti bahwa untuk sebarang
0 terdapat bilangan real
x0
B dan bilangan integer
.
0 stabil, maka dapat diasumsikan x0
untuk semua k
. Andaikan
x(k )
B tidak
0 stabil maka x(k ) terbatas, sehingga menurut
teorema Bolzano-Weirstrass (Teorema 2.2.10), terdapat subbarisan x(ni ) yang konvergen, katakanlah ke y terbuka dari y dengan 0 h( x)
V x ( k 1) V x(k )
n
. Misalkan E
E . Karena
V ( x)
B 1 menyatakan persekitaran
0 , maka dapat didefinisikan fungsi
yang terdefinisi dengan baik dan kontinu di E, dengan
h( x) 1 untuk semua x
E . Selanjutnya jika
h( y),1 , maka terdapat
0
35
sedemikian sehingga jika x E
berakibat h( x)
. Sehingga untuk ni yang
cukup besar, berlaku
V x(ni )
2
V x(ni 1)
V x(ni
2)
...
ni
V x0
yang artinya 0.
lim V x( ni )
ni
Tetapi karena V kontinu, sehingga menurut definisi (2.2.5), lim V x( ni ) ni
maka V ( y )
xe
0 , yang artinya y
0 . Sehingga terbukti bahwa lim x(k ) k
V ( y) ,
0 . Jadi,
0 stabil asimtotik.
Fungsi Lyapunov untuk suatu sistem tidak tunggal, sehingga untuk membentuk fungsi Lyapunov, biasanya dipilih bentuk yang paling sederhana, yaitu bentuk kuadratik n
n
V ( x)
qij xi x j
xT Qx
(2.4.7)
i 1 j 1
dengan x
( x1 ,..., xn )T dan Q adalah suatu matriks simetris berukuran n n .
Berdasarkan Teorema (2.2.6), fungsi berbentuk kuadratik (2.4.7) adalah kontinu, sehingga tinggal dicek apakah
V ( x)
0 terhadap state sistemnya.
2.4.2 Keterkendalian Sistem Diskrit LTI Berdasarkan solusi (2.4.2), dapat dianalisis sifat-sifat sistem selain kestabilan, antara lain keterkendalian. Definisi dan teorema berikut menjelaskan sifat keterkendalian sistem diskrit LTI yang diambil dari Olsder (2004).
Definisi 2.4.6 Sistem (2.4.1) dengan kondisi awal x (0) penuh jika untuk setiap state x0 , x1
n
, k
u(0), u(1),...u( j ) sedemikian sehingga x (k , x0 , u )
x0 dikatakan terkendali
, terdapat barisan kendali
x1 , dengan j
berhingga.
36
Berdasarkan definisi (2.4.6), sifat keterkendalian merupakan syarat perlu suatu sistem dapat dikendalikan. Teorema (2.4.7) berikut dapat digunakan untuk menentukan sifat keterkendalian sistem.
Teorema 2.4.7 Sistem (2.4.1), atau pasangan matriks ( A, B) , adalah terkendali penuh jika dan hanya jika rank B AB
An 1 B
n.
Sifat keterkendalian yang dimaksud pada Teorema (2.4.7) adalah keterkendalian penuh, artinya jika state sistem berdimensi n, maka state xi terkendali untuk semua i 1, 2,..., n . Untuk sistem yang tidak terkendali penuh, yaitu terdapat subruang tidak terkendali berdimensi m
n
dimana state
xi , i 1, 2,..., m tidak terkendali, maka perlu peninjauan lebih lanjut untuk memastikan suatu sistem dapat dikendalikan.
2.4.3 Keterstabilan Sistem Diskrit LTI Definisi dan teorema berikut menjelaskan sifat keterstabilan sistem yang diambil dari Olsder (2004).
Definisi 2.4.8 Sistem (2.4.1), atau pasangan matriks ( A, B) , dikatakan terstabilkan jika terdapat matriks F sedemikian sehingga semua nilai eigen dari
A BF dapat ditempatkan di open unit disc, yaitu daerah di bidang kompleks yang memiliki modulus kurang dari satu.
Berdasarkan definisi (2.4.8), maka syarat sistem yang tidak terkendali penuh supaya dapat dikendalikan adalah sistem terstabilkan, yaitu semua nilai eigen sistem dapat ditempatkan di open unit disc. (2.4.9) berikut dapat digunakan untuk menyelidiki sifat keterstabilan sistem.
37
( A, B) , bersifat
Teorema 2.4.9 Sistem (2.4.1), atau pasangan matriks terstabilkan jika dan hanya jika
rank
I
A B
untuk semua nilai eigen
n i
dari A dengan
i
1, i 1,2,..., k
n.
2.5 Teknik Kendali MPC Kendali model prediktif merupakan salah satu teknik kendali modern yang berpengaruh secara signifikan di bidang industri Teknik kendali MPC banyak diaplikasikan ke berbagai bidang, seperti bidang ekonomi, proses kimia, sistem mekanik, dan lain-lain. Konsep utama dari teknik kendali MPC adalah strategi receding
horizon
dan
prediksi,
kemudian
kendali
optimal
ditentukan
menggunakan teknik optimisasi.
2.5.1 Ide Dasar Teknik Kendali MPC Ide dasar dari teknik kendali MPC terdiri dari tiga konsep berikut. 1) Menggunakan model matematika yang bekerja pada plant untuk memprediksi keluaran pada waktu instant yang akan datang (horison) 2) Menghitung
barisan
kendali
prediksi
sepanjang
horison
yang
meminimumkan fungsi objektif, dalam hal ini adalah meminimumkan selisih antara prediksi keluaran dengan trayektori referensi 3) Menggunakan strategi receding horizon, dimana pada setiap waktu instant sepanjang horison, barisan kendali yang dihitung bergantung kepada barisan kendali yang telah dihitung pada setiap waktu instant sebelumnya, kemudian pada barisan kendali yang didapat, hanya sinyal kendali pertama yang diaplikasikan pada sistem.
Berdasarkan ide dasarnya, teknik kendali MPC dapat didefinisikan sebagai teknik kendali yang menggunakan model untuk memprediksi keluaran yang akan datang sepanjang horison dan menentukan barisan kendali prediksi yang mengoptimalkan fungsi objektif. Tujuan umum dari teknik kendali MPC yaitu
38
membawa prediksi keluaran yang akan datang sepanjang horison supaya sedekat mungkin dengan trayektori referensi. Pada teknik kendali MPC, tidak ada ketentuan berapa panjang horison yang digunakan, tetapi yang paling sering digunakan adalah 10 langkah. Semakin panjang horison prediksi, maka akan semakin banyak variabel yang digunakan pada optimisasinya, sehingga secara umum membutuhkan waktu penyelesaian yang semakin lama.
2.5.2 Strategi Teknik Kendali MPC Strategi teknik kendali MPC dapat dibagi menjadi 3 tahap utama berikut. Strategi ini dapat diilustrasikan pada Gambar (2.5.1). Pertama, keluaran yang akan datang sepanjang horison prediksi (Hp), namakan sebagai prediksi keluaran, diprediksikan pada setiap waktu instant k menggunakan model yang bekerja pada plant. Prediksi keluaran saat k dinotasikan sebagai yˆ k i | k , i 1, 2,..., H p yang nilainya bergantung pada masukan dan keluaran yang lalu sampai waktu instant-k. Barisan kendali prediksi yang akan datang saat k dinotasikan sebagai
uˆ k i | k , i
0,1,..., H p 1, dihitung setiap waktu instant, kemudian digunakan
untuk menghitung prediksi keluaran dan sinyal kendali untuk waktu instant selanjutnya.
Gambar 2.5.1 Strategi teknik kendali MPC
Kedua, barisan sinyal kendali yang akan datang dihitung dengan mengoptimalkan fungsi objektif. Fungsi objektif yang diambil biasanya berbentuk
39
fungsi kuadratik yang merupakan selisih antara prediksi keluaran dengan trayektori referensi ditambah suku yang menyatakan pembobotan kendali. Trayektori acuan
Masukan dan keluaran yang lalu MODEL Masukan yang akan datang Pengoptimal
Fungsi objektif
Prediksi keluaran
+
Error yang akan datang
kendala
Gambar 2.5.2 Diagram blok strategi teknik kendali MPC
Ketiga, sinyal kendali u k | k diaplikasikan ke sistem, kemudian dihitung nilai y k 1
yang merupakan nilai keluaran baru saat (k+1). Selanjutnya
dihitung sinyal kendali u k 1| k 1 dan seterusnya sepanjang periode kendali. Strategi MPC ini dapat digambarkan sebagai diagram blok (2.5.2).
Teknik kendali MPC terdiri dari 3 elemen yaitu model prediksi, fungsi objektif, dan hukum kendali. Model prediksi yang dipakai yaitu model matematika yang bekerja pada plant dalam bentuk persamaan state space. Fungsi objektif digunakan sebagai upaya meminimalkan usaha yang dilakukan dan meminimalkan selisih antara prediksi keluaran dan trayektori referensi.
2.5.3 Model Matematika dari Plant Model matematika yang digunakan untuk menentukan state dan keluaran prediksi sepanjang horison prediksi adalah persamaan diferensi dari plant. Misalkan k menyatakan waktu instant, persamaan diferensi pada plant diasumsikan berbentuk persamaan state space diskrit LTI
40
x (k )
Ax(k 1) Bu (k 1),
(2.5.1a)
y(k ) Cx(k ),
(2.5.1b)
dengan x menyatakan vektor state berdimensi n, u menyatakan vektor masukan berdimensi p, dan y menyatakan vektor keluaran berdimensi m. Matriks-matriks A, B, dan C secara berturut-turut adalah matrik real konstan berukuran
n n, n p, dan m n . Model (2.5.1) ekuivalen dengan sistem diskrit (2.4.1) dengan asumsi matriks D = 0, yang artinya keluaran hanya bergantung pada state, dan tidak dipengaruhi secara langsung oleh masukan.
Diasumsikan state x(k ) diketahui nilainya yang dihitung berdasarkan state nilai-nilai x(k 1) dan u(k 1) . Prediksi state xˆ(k
i | k ) untuk i 1, 2,..., H p
ditentukan menggunakan prosedur rekursif berikut.
xˆ (k 1 k )
Ax(k ) Buˆ (k k ),
xˆ (k
Ax(k 1 k )
2 k)
Buˆ ( k 1 k )
A Ax(k ) Buˆ(k k ) A2 x ( k )
xˆ (k 3 k )
Ax(k
ABuˆ (k k )
2 k)
A A2 x(k ) A3 x(k ) xˆ (k
H p k)
Ax(k H
Buˆ ( k
Buˆ (k 1 k ) Buˆ (k 1 k ),
2 k)
ABuˆ (k k ) Buˆ (k 1 k )
A2 Buˆ(k k )
ABuˆ (k 1 k )
H p 1 k ) Buˆ (k
A p x(k )
A
Hp 1
Buˆ(k Buˆ (k
2 k) 2 k ),
H p 1 k)
Buˆ (k k ) ... Buˆ (k
H p 1 k ).
Pada iterasi diatas, digunakan notasi prediksi masukan uˆ(k | k ) karena pada saat menentukan prediksi, nilai u(k | k ) belum diketahui. Misalkan Hu menyatakan horison kendali, yaitu panjang horison yang ditentukan untuk menghitung prediksi masukan. Selanjutnya, diasumsikan bahwa nilai masukan hanya berubah pada waktu k , k 1,..., k H u 1 dan bernilai konstan setelahnya,
41
sehingga untuk Hu
j H p 1 nilai uˆ (k
j | k ) uˆ (k
diambil Hu = 1. Selanjutnya didefinisikan
H u 1) . Tetapi biasanya
uˆ(k i | k ) menyatakan displacement uˆ (k i | k ) uˆ (k i | k ) uˆ (k i 1| k ) .
prediksi masukan yang dirumuskan oleh
Pada saat k, nilai u(k 1) sudah diketahui, sehingga prediksi masukan sepanjang horison kendali ditulis sebagai
uˆ (k k )
uˆ (k | k ) u (k 1),
uˆ(k 1 k )
uˆ (k 1| k )
uˆ(k
2 k)
uˆ (k
uˆ (k
Hu 1 k )
uˆ (k | k ) u (k 1), uˆ (k 1| k )
uˆ (k | k ) u (k 1),
H u 1| k ) ...
uˆ (k | k ) u (k 1).
2 | k)
uˆ (k
Prediksi masukan uˆ disubstitusikan ke prediksi state xˆ diperoleh
xˆ (k 1 k )
Ax(k ) B
xˆ (k
A2 x ( k )
2 k)
uˆ (k | k ) u (k 1) , uˆ(k k ) u (k 1)
AB
uˆ (k 1 k )
B
A2 x ( k )
uˆ (k k ) u(k 1)
AB uˆ (k k )
ABu (k 1) B uˆ (k 1 k )
B uˆ(k k ) Bu(k 1) A2 x(k ) ( A I ) B uˆ (k k ) B uˆ (k 1 k ) ( A I ) Bu (k 1), xˆ ( k
3 k)
A3 x( k ) AB B
A2 B
uˆ (k | k ) u (k 1)
uˆ (k 1 | k ) uˆ ( k
uˆ (k | k ) u (k 1) uˆ (k 1| k )
2 | k)
A3 x( k )
A2 B uˆ (k | k )
AB uˆ ( k 1| k ) B uˆ (k
xˆ(k
Hu k )
uˆ ( k | k ) u ( k 1)
A2 Bu ( k 1)
AB uˆ ( k | k )
ABu ( k 1)
2 | k ) B uˆ (k 1| k ) B uˆ (k | k )
A3 x( k )
A2
B uˆ( k
2 k)
A I B uˆ (k | k ) A2
AHu x(k )
AH u
B uˆ(k
Hu 1 k )
1
Bu ( k 1)
A I B uˆ (k 1| k )
A I Bu ( k 1),
... A I B uˆ (k k ) ... AH u
1
... A I Bu (k 1),
42
xˆ (k
A H u 1 x (k )
Hu 1 k)
AHu
A I B uˆ(k
xˆ (k
H p k)
H
A p x(k ) A
H p Hu
A
Hp 1
A
Hp 1
... A I B uˆ (k k ) ... Hu 1 k )
A Hu
... A I Bu(k 1),
... A I B uˆ (k k ) ...
... A I B uˆ (k
Hu 1 k )
... A I Bu (k 1).
Persamaan state prediksi sepanjang horison Hp dapat dituliskan dalam bentuk matriks menjadi B A
xˆ (k 1 k )
Hu 1
Ai B xˆ(k xˆ (k
AHu x (k ) AHu 1
Hu k ) Hu 1 k )
i 0 Hu
Ai B
u (k 1)
i 0
xˆ(k
H p k)
A
Hp Hp 1
Ai B i 0
B
0
AB
B
0
Ai B
B
Hu 1
uˆ (k k )
i 0
.
Hu i
AB
AB
B
i 0
Hp 1
H p Hu
Ai B i 0
Ai B i 0
Prediksi keluaran y dapat dituliskan sebagai
uˆ(k
Hu 1 k )
(2.5.2)
43
yˆ (k 1 k ) Cxˆ (k 1 k ), yˆ (k yˆ (k
2 k ) Cxˆ (k H p k ) Cxˆ (k
2 k ), H p k ),
atau dapat ditulis dalam bentuk matriks menjadi yˆ ( k 1 k ) yˆ ( k 2 k ) yˆ ( k
H p k)
xˆ ( k 1 k ) xˆ ( k 2 k )
C
0
0
0
C
0
0
0
C
xˆ ( k
(2.5.3)
.
H p k)
2.5.4 Fungsi Objektif Fungsi objektif untuk teknik kendali MPC terdiri dari dua suku. Suku pertama didefinisikan sebagai selisih antara prediksi keluaran dengan trayektori referensi
r (k i | k ), i
H w , H w 1,..., H p , dengan Hw menyatakan horison
window, yaitu waktu instant dimana selisih
yˆ k i | k
r k i
mulai
ditentukan, tetapi biasanya diambil Hw = 1. Suku kedua adalah pembobotan masukan yang mewakili besarnya usaha kendali yang dilakukan, sehingga suku ini dapat direpresentasikan sebagai fungsi energi tiruan. Suku pertama diasumsikan semidefinit positif, dan suku kedua diasumsikan definit positif, sehingga menurut Teorema (2.3.9), fungsi objektif ini stricly konveks. Hal ini dilakukan supaya eksistensi peminimum global terjamin berdasarkan Teorema (2.3.13). Selain itu, rumus fungsi energi (energi potensial dan energi kinetik) juga berbentuk kuadratik. Hal ini menjadi salah satu motivasi untuk menyusun fungsi objektif yang berbentuk kuadratik. Fungsi objektif yang digunakan diasumsikan berbentuk kuadratik yang didefinisikan oleh persamaan berikut. Hp
yˆ k
J (k )
i|k
r k
i
T
Q (i ) yˆ k
i Hw
r k
i
T
(2.5.4)
Hu 1
uˆ k i 0
i|k
i|k
T
R (i )
uˆ k
i|k .
44
Fungsi objektif (2.5.4) dapat dinyatakan menggunakan notasi bentuk kuadratik P
2
PT QP , yaitu
Q
J (k )
k
k
2
k
2
,
(2.5.5)
dengan yˆ ( k yˆ ( k
(k )
yˆ ( k
Q Hw 0
Hw k) Hw 1 k)
r (k r (k
(k )
,
H p k)
0 0
0
Q Hp
0
Hw 1 k)
r (k
0 Q Hw 1
uˆ ( k k )
Hw k) ,
(k )
uˆ ( k
H p k)
,
uˆ ( k 1 k ) Hu 1 k)
R 0
0
0
0
R 1
0
0
0
R Hu 1
.
Matriks pembobotan Q(i) untuk setiap i diasumsikan semidefinit positif. Sedangkan matriks R(i)
untuk setiap i diasumsikan definit positif. Sehingga
menurut Teorema (2.3.9), fungsi objektif (2.5.5) merupakan fungsi stricly konveks. Pada umumnya, suatu proses terbatas pada suatu rentang nilai tertentu. Kendala masalah kendali sistem diskrit secara umum dapat dinyatakan oleh
umin
u (k )
umax ,
(2.5.6)
umin
u (k ) umax ,
(2.5.7)
ymin
y (k )
(2.5.8)
ymax .
Kendala (2.5.6) menyatakan batasan displacement masukan. Kendala (2.5.7) menyatakan batasan nilai masukan yang dapat diaplikasikan ke sistem. Kendala (2.5.8) menyatakan batasan nilai keluaran sistem.
2.5.5 Hukum Kendali Hukum kendali digunakan untuk menentukan barisan kendali optimal dari teknik kendali MPC. Berdasarkan Teorema (2.3.9), fungsi objektif (2.5.4)
,
45
merupakan fungsi kuadratik stricly konveks. Sehingga pada pembahasan ini, hukum kendali yang digunakan yaitu menentukan sinyal kendali optimal yang meminimumkan fungsi objektif menggunakan teknik optimisasi fungsi konveks. Selain itu, hukum kendali yang dilakukan adalah menggunakan konsep receding horizon, dimana barisan kendali optimal dihitung berdasarkan nilai kendali pada waktu sebelumnya dan hanya elemen pertama dari barisan kendali optimal yang diperoleh yang diaplikasikan ke sistem untuk setiap waktu instant, dan berulang untuk waktu instant berikutnya.
2.5.6 Solusi Teknik Kendali MPC Tanpa Kendala Tinjau kembali persamaan prediksi keluaran (2.5.2). Persamaan (2.5.2) disubstitusikan ke prediksi keluaran (2.5.3) didapat B A yˆ( k 1 k ) yˆ( k 2 k )
C 0
0 C
AHu AHu
yˆ ( k
H p k)
0
0
Hu 1
Ai B
0 0
i 0 1
x(k )
Hu
Ai B
u ( k 1)
i 0
C A
Hp Hp 1
Ai B i 0
B
0
AB
B
0
Ai B
B
Hu 1
uˆ(k k )
i 0
,
Hu i
AB
AB B
i 0
Hp 1
H p Hu
Ai B i 0
Ai B i 0
uˆ(k
Hu 1 k )
46
B A yˆ(k 1 k ) yˆ (k 2 k )
C 0 0 C
0 0
A
Hu 1
Hu
AHu yˆ(k
H p k)
0
x(k )
1
C
0
C 0 0 C
i 0 Hu
Ai B 0
A
Ai B
0 0
i 0
C
0
u (k 1)
Hp Hp 1
Ai B i 0
B
0
AB C
0
0
Hu 1
0
C
0
i 0
B
0
Ai B
B
uˆ(k k ) .
Hu i
0
0
AB
C
AB
uˆ(k
B
Hu 1 k )
i 0
Hp 1
H p Hu
Ai B
Ai B
i 0
i 0
Didefinisikan matriks-matriks B A C 0 0 C 0
0
0 0
A Hu , AH u 1
C
Hu 1
C 0 0 C
Ai B i 0 Hu
Ai B 0
A
0 0
0
C
i 0
Hp Hp 1
Ai B i 0
,
47
B
0
AB C
0
0
Hu 1
0
C
0
i 0
B
0
Ai B
B .
Hu i
0
0
AB
C
AB
B
i 0
Hp 1
H p Hu i
Ai B
AB i 0
i 0
Maka persamaan prediksi keluaran (2.5.5) dapat dituliskan menjadi
(k )
x( k )
u( k 1)
( k ).
(2.5.9)
Didefinisikan
(k )
(k )
x( k )
u (k 1).
(2.5.10)
Persamaan (2.5.10) disebut dengan tracking error, yang menyatakan selisih antara trayektori referensi dengan sistem free response
(k )
menandakan bahwa
0 . Apabila
(k )
(k) = 0,
0 merupakan langkah yang tepat karena free
response sama dengan trayektori referensi. Persamaan (2.5.9) disubstitusikan ke fungsi objektif (2.5.3) diperoleh
J (k )
x (k )
u (k 1)
(k )
x( k )
(k )
(k )
(k )T ( k )T
Misalkan
2
T
T
T
(k ),
k
2
u(k 1) T k
2
(k )
(k ) 2 k
(k )
k T
T
2
2
k
2
2
( k )T ( k )T
k
(k )
k
T
k
(k ) (2.5.11)
k . T
dan karena
(k )T
(k )
konstanta ,
maka persamaan (2.5.11) menjadi
J (k )
( k )T
(k )T
(k ) konstanta.
(2.5.12)
48
(k ) optimal dapat ditentukan dari gradien J(k).
Berdasarkan teori optimisasi,
Persamaan (2.5.12) mencapai optimal jika (k )
J(
( k ))
2 2
(k )
(k )
J(
(k )) 0 .
0
(k ) 1 2
(k ) 1
(k )
1 2
(k ) 1 2
1
1
.
Sehingga diperoleh barisan kendali optimal
1 2
(k )opt
1
(2.5.13)
.
Vektor masukan pada persamaan sistem (2.5.1) memiliki p elemen. Sehingga, berdasarkan strategi receding horizon, hanya digunakan p baris pertama
(k )opt , yang dapat dituliskan sebagai
dari matriks
u(k )opt
1 p ,0 p ,...,0 p
(k )opt .
(2.5.14)
( H u 1) kali
Solusi kendali optimal ditulis dalam bentuk
u(k )opt bukan dalam bentuk
uˆ(k )opt , karena sudah dapat ditentukan masukan optimal sebenarnya, dan bukan u(k )opt merupakan solusi optimal yang
prediksi lagi. Untuk menjamin bahwa
meminimumkan J(k), perlu diselidiki lebih lanjut dengan menentukan matriks Hessian dari J(k), yaitu 2
J (k )2
2
2
T
.
Diasumsikan bahwa Q(i ) 0 untuk setiap i, sehingga
(2.5.15)
T
0 . Supaya
menjamin solusinya adalah peminimal tegas, matriks hessian harus definit positif, yang dipenuhi bila apabila
> 0, artinya R(i ) 0 untuk setiap i. Untuk kasus khusus
=0, untuk menjamin solusi minimum dan menjamin
invertibel
49
maka haruslah T
T
0 . Untuk kasus lainnya, yaitu jika
0 , maka haruslah
0.
2.5.7 Solusi MPC Berkendala Tinjau kembali kendala masalah MPC (2.5.6)-(2.5.8). Ketiga kendala tersebut dapat dituliskan menjadi bentuk umum kendala masalah MPC yang dapat dinyatakan oleh E
uˆ ( k k ),..., uˆ ( k
F uˆ ( k k ),..., uˆ ( k
G yˆ (k
H u 1 k ), 1 p
H u 1 k ), 1 p
H w k ),..., yˆ ( k
T
H p k ), 1 p
T
0
(2.5.16) (2.5.17)
0 T
0,
(2.5.18)
dengan E, F, dan G merupakan matriks-matriks konstan. Kendala (2.5.16) merepresentasikan batasan perubahan displacement aktuator, kendala (2.5.17) merepresentasikan batasan kemampuan aktuator, dan (2.5.18) merepresentasikan batasan nilai keluaran. Kendala (2.5.16)-(2.5.18) dapat ditulis dalam bentuk
E
F
G dengan
(k )
0
1p (k ) 1p
(k ) 1p (k )
(2.5.19)
0
(2.5.20)
0,
(2.5.21) T
uˆ ( k k ) ,..., uˆ ( k
dinyatakan dalam variabel dalam
H u 1 k )T
T
. Karena fungsi objektif (2.5.12)
(k ) , maka ketiga kendala diatas akan dinyatakan
(k ) .
Misalkan matriks F berbentuk F
F1 , F2 ,..., FH u , f , dimana masing-
masing Fi adalah matriks berukuran q m dan f adalah vektor berdimensi q, sehingga (2.5.20) dapat ditulis sebagai
50
u(k | k ) u (k 1| k ) F1 , F2 ,..., FH u , f
0 u (k
H u 1| k ) 1
Fu 1 (k | k )
F2u (k 1| k ) ... FH u u ( k
H u 1| k )
f
0
Hu
ˆ(k i 1 k ) Fu i
0.
f
(2.5.22)
i 1
i 1
Karena uˆ(k
uˆ (k
i 1 k ) u (k 1)
j k ) , maka (2.5.22) dapat
j 0
ditulis kembali menjadi Hu
i 1
uˆ (k
Fi u (k 1) i 1
j k)
f
uˆ (k
j k)
0
j 0 Hu
Hu
Fu 1) i (k i 1
i 1
Fi i 1
f
0
j 0
F1 uˆ (k k ) F2 uˆ (k k ) F2 uˆ (k 1 k ) ... Hu
FHu uˆ(k k ) ... FHu uˆ (k
Hu 1 k )
Fu 1) i (k
f
0
i 1 Hu
Hu
Fj uˆ(k k ) j 1
Fj uˆ (k 1 k ) ... j 2
(2.5.23)
Hu
FHu uˆ (k
Hu 1 k )
Fj u (k 1)
f
0.
j 1
Hu
Selanjutnya didefinisikan
Fj
i
dan
1
,...,
, maka
Hu
j i
(2.5.23) dapat ditulis kembali menjadi Hu
Hu
Fj j 1
Hu
Fj j 2
Fj
uˆ(k k ) uˆ(k 1 k )
j Hu
Hu
F j u (k 1) j 1
uˆ (k
Hu 1 k )
f
0
51
1
2
(k )
Hu
(k )
u (k 1)
u ( k 1)
1
f
f,
1
(2.5.24)
dimana ruas kanan (2.5.24) merupakan suatu vektor yang diketahui nilainya saat k. Dengan langkah-langkah yang sama, akan dirubah (2.5.21) dalam
(k ) .
Untuk merubahnya, digunakan persamaan (2.5.6), sehingga (2.5.21) dapat ditulis sebagai
x( k )
G
u(k 1)
(k )
1
Misalkan G
,g
0.
(2.5.25)
dengan g menyatakan kolom terakhir dari G, sehingga
pertidaksamaan (2.5.25) ditulis menjadi
x(k )
u(k 1)
(k ) Misalkan E
W
(k )
x( k )
u (k 1)
g
0
(2.5.26)
g.
(2.5.27)
W , w maka pertidaksamaan (2.5.19) dituliskan menjadi (k ) w 0
W
(k )
w.
(2.5.28)
Kendala pertidaksamaan (2.5.24), (2.5.27), dan (2.5.28) dituliskan dalam bentuk matriks menjadi u (k 1) f x(k ) u ( k 1) 1
(k ) W
g .
(2.5.29)
w
Selanjutnya, akan dikonstruksikan masalah MPC berkendala menjadi masalah pemrograman kuadratik berkendala. Fungsi obektif yang akan diminimumkan adalah fungsi objektif pada (2.5.12) yaitu
min
(k )
(k )T
(k )T
(k )
(2.5.30)
dengan kendala u (k 1) f x(k ) u ( k 1) 1
(k ) W
w
g .
(2.5.31)
52
Untuk kasus khusus dimana kendalanya berbentuk interval, yaitu
ul (k i ) uˆ (k i | k ) uh (k matriks
i), i
0,1, 2,..., H u 1 ,
maka
akan
terbentuk
yang relatif sederhana. Pembentukannya akan dibahas sebagai berikut.
Kendala interval tersebut dapat dituliskan menjadi
uˆ(k i | k ) uh (k i) 0 uˆ(k i | k ) ul (k i) 0.
0,1, 2,..., H u 1, kendala diatas dapat ditulis menjadi
Untuk i
uˆ(k | k ) uh (k ) 0 uˆ( k | k ) ul (k ) 0 uˆ(k 1| k ) uh (k 1)
0
uˆ( k 1| k ) ul (k 1) 0 uˆ(k H u 1| k ) uh (k H u 1) 0 uˆ(k H u 1| k ) ul (k H u 1) 0. i 1
Substitusi uˆ(k
uˆ (k
i 1 k ) u (k 1)
j k ) ke pertidaksamaan diatas
j 0
diperoleh u (k 1) u ( k 1) u ( k 1) u ( k 1) u ( k 1) u ( k 1)
uˆ ( k k ) uh ( k ) 0 uˆ( k k ) ul (k ) uˆ (k k ) uˆ( k k ) uˆ ( k k ) uˆ( k k )
0
uˆ( k 1 k ) uh ( k 1)
0
uˆ (k 1 k ) ul (k 1) 0 uˆ( k 1 k ) ... uˆ (k 1 k ) ...
uˆ ( k uˆ ( k
H p 1 k ) uh ( k H p 1 k ) ul ( k
Dengan memindahkan suku-suku konstant ke ruas kanan, diperoleh uˆ( k k )
u ( k 1) uh ( k )
uˆ ( k k ) u ( k 1) ul ( k ) uˆ ( k k ) uˆ ( k k )
uˆ( k 1 k )
u ( k 1) uh ( k 1)
uˆ ( k 1 k ) u ( k 1) ul ( k 1)
H u 1) 0 H u 1)
0.
53
uˆ ( k k )
uˆ ( k 1 k ) ...
uˆ ( k k )
uˆ ( k
uˆ ( k 1 k ) ...
H p 1 k)
uˆ ( k
u (k 1) uh (k
H u 1)
H p 1 k ) u ( k 1) ul ( k
H u 1).
Pertidaksamaan diatas dapat dituliskan dalam bentuk matriks menjadi 1 1 1
0 0 1
0 0 0
0 0 0
0 0 0
0 0 0
1
1
0
0
0
0
uˆ(k k ) uˆ( k 1 k ) uˆ(k 2 k ) uˆ( k 3 k )
1 1
1 1
1 1
1 1
1 1
1 1
uˆ(k H u 2 k ) uˆ( k H u 1 k )
uh ( k )
1 1
ul (k ) uh (k 1) ul ( k 1)
1 1 u (k 1)
.
uh ( k H u 1) ul ( k H u 1)
1 1
Pertidaksamaan diatas bersesuaian dengan pertidaksamaan (2.5.22) dimana 1 1 1
0 0 1
0 0 0
0 0 0
0 0 0
0 0 0
1
1
0
0
0
0 , f
1
1
1
1
1 1
1 1
1
1
1
1
uh ( k ) ul ( k ) uh ( k 1) ul ( k 1)
.
uh ( k H u 1) ul (k H u 1)
Berdasarkan uraian diatas, teknik kendali MPC berkendala menjadi masalah pemrograman kuadratik berkendala pertidaksamaan. Salah satu metode yang dapat digunakan untuk menyelesaikan pemrograman kuadratik berkendala pertidaksamaan secara numerik adalah metode active-set yang sudah diberikan di sub-subbab (2.4.4).
54
Contoh 2.5.1 Suatu sistem memiliki persamaan state space x(k 1)
x(k )
1,5
y (k )
x(k ),
u1 (k ) , x(0) u2 ( k )
1,5
0
dengan fungi objektif 4
T
J (k )
x ( k 1 i ) r 10 x ( k 1 i ) r
u ( k i )T 100 u (k i ).
i 0
Akan dikendalikan state x(k) untuk menuju set point r 10 dengan horison prediksi H p
5 , horison window H w
1 , dan horison kendali H u
5 . Simulasi
dilakukan dengan 12 langkah waktu instant. Kendali optimal ditentukan menggunakan persamaan (2.5.14) yang termuat dalam progam MATLAB pada Lampiran (1). Sistem hasil kendali dapat disimulasikan sebagai grafik inputoutput (2.5.2) berikut.
Input 1
3 2 1 0
2
4
6 Waktu Instant
8
10
12
2
4
6 Waktu Instant
8
10
12
2
4
6 Waktu Instant
8
10
12
Input 2
-1 -2 -3 0
Output
20 10 0 0
Gambar 2.5.2 Grafik input dan output Contoh (2.5.1)
Berdasarkan grafik input-output (2.5.2), dapat dilihat bahwa output mencapai set point pada waktu instan ke-8. Sebelum mencapai set point r = 10, input 1 dan input 2 bekerja memberikan kendali optimal yang dihitung oleh kontroler MPC sehingga state mencapat set point.
55
2.6 Teori Permainan Dinamis Kooperatif Diskrit Diketahui terdapat M pemain. Masalah permainan dinamis kooperatif memuat M pemain yang bekerjasama untuk menentukan strategi optimal sedemikian sehingga meminimumkan fungsi objektif semua pemain. Persamaan dinamik dari masalah permainan kooperatif dapat dinyatakan sebagai persamaan state space LTI diskrit, yaitu
x(k 1)
Ax(k ) B1u1 (k ) B2u2 (k ) ... BM uM ( k ), x(0) n
dengan x(k )
menyatakan vektor state, ui (k )
pi
x0 ,
(2.6.1)
menyatakan strategi
(vektor kendali) pemain ke-i, A dan Bi secara berturut-turut merupakan matriksmatriks real konstan berukuran n n dan n pi , i = 1, 2, …, M. Misalkan B
B1 , B2 ,..., BM
dan u
u1 , u2 ,..., u M
T
, persamaan diferensi
(2.6.1) dapat dituliskan sebagai
x(k 1)
Ax(k ) Bu(k ), x(0)
x0 ,
(2.6.2)
yang merupakan persamaan state space diskrit LTI. Solusi persamaan diferensi tersebut dapat ditentukan dengan prosedur iterasi seperti pada solusi sistem diskrit, yaitu k 1
x( k )
Ak
Ax0
i 1
Bu (0).
(2.6.3)
i 0
Fungsi cost untuk waktu proses sepanjang N yang akan diminimumkan oleh masing-masing pemain ke-i diasumsikan berbentuk kuadratik yang dinyatakan oleh J i xi ( k ), ui ( k
1N 2k
1
xi ( k )T Qi xi ( k ) ui ( k )T Ri ui ( k ) ,
(2.6.4)
0
dengan matriks-matriks Qi diasumsikan simetrik semidefinit positif dan Ri diasumsikan definit positif untuk semua i. Didefinisikan himpunan konveks pi
menyatakan ruang strategi pemain ke-i, u ( k )
i
dan u(k )
dengan
1
2
M
min J i xi ( k ), ui ( k ) , i 1, 2,..., M , ui
u1 (k ), u2 ( k ),..., uM ( k )
T
. Diberikan masalah minimisasi (2.6.5)
56
dengan kendala ui
i
. Berdasarkan Teorema (2.3.9), Ji stricly konveks di
i
i 1,2,..., M , sehingga masalah optimisasi (2.6.5) merupakan
untuk setiap
masalah optimisasi fungsi stricly konveks dengan M fungsi objektif.
Definisi 2.6.1 Strategi
dikatakan solusi optimal lengkap masalah
optimisasi (2.6.5) jika dan hanya jika untuk semua
berlaku J i ( )
J i ( ),
untuk setiap i 1,2,..., M .
Pada umumnya, solusi optimal lengkap tidak dapat ditentukan ketika terdapat konflik antar fungsi objektif, yaitu ketika suatu strategi menjadi optimal untuk satu pemain, tetapi mengakibatkan hasil pemain lain menjadi lebih buruk. Sehingga diberikan konsep Pareto efficient atau solusi optimal Pareto.
2.6.1 Solusi Pareto Konsep dasar dari Pareto efficient yaitu menentukan solusi yang sudah tidak dapat diminimumkan lagi oleh semua pemain secara simultan.
Definisi 2.6.2 Strategi ˆ
disebut Pareto efficient masalah optimisasi (2.6.5)
jika himpunan pertidaksamaan
Ji ( )
J i ( ˆ),
i 1,..., M ,
(2.6.6)
dan paling sedikit satu pertidaksamaan bersifat strict, tidak terjadi untuk manapun
ˆ.
dengan
J1 ( ˆ ),..., J M ( ˆ)
M
Nilai
objektif
yang
berkorespondensi
yaitu
disebut solusi Pareto. Himpunan yang beranggotakan
semua solusi Pareto disebut Pareto frontier.
M
( 1,
Didefinisikan
2
,...,
M
)|
i
0 dan
i
1 . Dibentuk
i 1
kombinasi
linier
konveks
dari
fungsi-fungsi
objektif
(2.6.5)
yaitu
57
M
J( )
i
J i ( ) . Lemma berikut menjelaskan bahwa peminimum dari
i
J i ( ) merupakan Pareto efficient.
i 1 M
J( ) i 1
. Jika ˆ
Lemma 2.6.3 Diberikan
sehingga
M
ˆ arg min
i
Ji ( ) ,
(2.6.7)
i 1
maka ˆ adalah Pareto efficient untuk masalah optimisasi (2.6.5).
M
dan ˆ arg min
Bukti : Diketahui
i
J i ( ) . Andaikan ˆ bukan
i 1
Pareto efficient, maka terdapat
Ji ( ) dan J j ( )
J i ( ˆ),
sehingga
i 1,..., M ,
J j ( ˆ ) untuk paling sedikit satu j. Akibatnya,
J( )
1 1
J ( ) ...
2 2
M
J ( ˆ)
JM ( )
1 1
M
2
J 2 ( ˆ ) ...
M
J M ( ˆ)
M i
Ji ( )
i
i 1
J i ( ˆ ).
i 1 M
Dengan kata lain ˆ bukan peminimum
i
J i ( ) . Hal ini kontradiksi dengan
i 1 M
yang diketahui bahwa ˆ peminimum
i
J i ( ) . Sehingga ˆ adalah Pareto
i 1
efficient.
Selanjutnya
diberikan
lemma
berikut
yang
menjamin
eksistensi
pembobotan kombinasi linier konveks apabila diberikan suatu Pareto efficient.
Lemma 2.6.4 Diasumsikan ruang strategi
i
konveks dan Ji konveks untuk setiap
i 1, 2,..., M . Jika ˆ merupakan Pareto efficient untuk masalah minimisasi (2.6.5), maka terdapat
sehingga untuk semua
berlaku
58
M
M i
Ji ( ˆ)
i
i 1
(2.6.8)
J i ( ).
i 1
Z :
M
z
M
, didefinisikan himpunan Z
Bukti : Untuk semua | zi
yaitu
J i ( ˆ) , i 1,..., M .
Ji ( )
(2.6.9)
Didefinisikan pula
Z:
Z .
Karena ˆ Pareto efficient, maka 0 Z . Karena Ji konveks, maka Z konveks dan
z
gabungannya
juga
Z , z Z , dan zi
J i ( ) (1 Ji
Sehingga
yaitu
Z
konveks.
Untuk
sebarang
[0,1] berlaku ) zi
(1
konveks,
z (1
)z
) Ji ( )
(1
Z
(1
Ji ( ˆ)
)
J i ( ˆ ).
Z.
Dengan
)
berdasarkan Teorema (2.3.2), terdapat p
mengambil
0 sehingga pT z
x0
pT x0
0, x0
Z,
pT z 0
untuk semua z Z . Berdasarkan persamaan (2.6.9), untuk setiap i 1,..., M , dapat dipilih bilangan zi yang relatif besar. Karena pT z
0 , berakibat elemen ke-i
dari p tidak bernilai negatif untuk setiap i 1,..., M . Diberikan sebarang z maka terdapat
zi
M
dan
J i ( ) J i ( ˆ)
i
,
0, i 1,.., M sehingga
i
, i 1,..., N .
Oleh karena itu, untuk semua
Z
dan untuk semua
(2.6.10) i
0 berlaku
M
pT z
pi J i ( ) J i ( ˆ )
i
i 1
berlaku
Akibatnya untuk semua M
pi J i ( ) J i ( ˆ)
i
0
i
0
i 1 M
M
pi J i ( ˆ )
pi J i ( ) i 1
i 1
0.
(2.6.11)
59
M
M
pi J i ( ˆ )
pi J i ( ) i 1
M
M
pi J i ( ˆ ).
pi J i ( ) i 1
Dengan mengambil
i
i 1
(2.6.12)
i 1
i
pi
:
, untuk semua
M
berlaku
pj j 1
M
pi M
i 1
M
pi
Ji ( )
M i 1
pj
Ji ( ˆ) ,
pj
j 1
j 1
M
M
atau
i
Ji ( )
i
i 1
J i ( ˆ ).
i 1
Berdasarkan lemma (2.6.3) dan lemma (2.6.4) dapat dibentuk akibat berikut yang menyatakan bahwa terdapat pemetaan bijektif antara
dan
Pareto frontier.
Akibat 2.6.5 Diasumsikan 2, …, M. Diberikan
i
konveks dan Ji stricly konveks untuk semua i = 1,
dan M
ˆ
arg min
i
Ji ( ) .
(2.6.13)
i 1
Jika ˆ
ˆ 2,
1
1
2
, maka terdapat pemetaan bijektif antara
dan
Pareto frontier.
Bukti : Karena
i
konveks dan Ji konveks tegas untuk semua i = 1, 2, …, M,
maka terdapat bijeksi antara semua solusi ˆ
. Selanjutnya akan ditunjukkan bahwa solusi ˆ
dengan (1)
dari masalah optimisasi (2.6.13)
ˆ
(2)
,
(1)
(2)
(1)
dan ˆ
( 2)
dengan
berkorespondensi dengan solusi Pareto yang berbeda.
60
Andaikan ada solusi Pareto J dengan J ( ˆ (1)
(2)
. Hal ini berakibat ˆ
( 2)
(1)
J( ˆ
)
( 2)
) untuk beberapa
merupakan solusi dari masalah optimisasi
M (1) i i
min
J( ) .
i 1
Tetapi berdasarkan asumsi, fungsi cost tersebut bersifat stricly konveks, sehingga
ˆ
memiliki solusi minimum tunggal. Akibatnya,
( 2)
ˆ
(1)
, sehingga terjadi
kontradiksi dengan yang diketahui bahwa strategi optimal yang berbeda (i)
berkorespondensi dengan bijektif antara
yang berbeda. Jadi, terbukti terdapat pemetaan
dan Pareto frontier.
Berdasarkan uraian dan lemma-lemma diatas, dapat dibentuk Teorema (2.6.6) berikut. Teorema (2.6.6) inilah yang memberikan solusi Pareto kooperatif dari suatu permainan dinamis kooperatif.
Teorema 2.6.6 Diberikan masalah optimisasi min
J i xi (k ), ui ( k )
Diasumsikan Qi
1N 2k
1
xi (k )T Qi xi ( k ) ui ( k )T Riui ( k ) .
(2.6.14)
0
0, Ri >0, J i (u) konveks untuk semua i 1,..., M . Himpunan
semua solusi Pareto Kooperatif diberikan oleh J1 u (
, J2 u (
,..., J M u (
,
,
yang berkorespondensi dengan Pareto efficient yang diperoleh dari masalah minimisasi M
u ( )
arg min u U
dengan kendala x(k 1)
i
Ax(k ) B1u1 (k ) ... BM uM (k ), x(0)
Bukti : Diberikan fungsi cost masing-masing pemain Ji
1N 2k
1
x ( k )T Qi x ( k ) ui ( k )T Ri ui ( k ) , i 1, 2,..., M , 0
(2.6.15)
J i ( xi ( k ), ui ( k )),
i 1
x0 .
61
k 1
dengan x(k )
Ak x0
i 1
Ak
Bu (i ), B
B1 , B2 ,..., BM , u
u1, u2 ,..., uM
T
.
i 0
Ambil sebarang dua input u(k) dan v(k), dengan u
T
u1 , u2 ,..., u M
dan v
T
v1 , v2 ,..., vM
.
Didefinisikan notasi k 1
xu (k )
Ak x0
Ak
i 1
Ak
i 1
Bv(i).
1
x
Bu (i ),
i 0 k 1
xv (k )
Ak x0 i 0
Akan dibuktikan Ji konveks. Ji
u
(1
1N 2k
)v
1N 2k
u (1
ui
0
1
)v
T
(k )
ui
T
) vi
(1
u (1
Ri
)v
ui T
xu ( k ) (1
0
Qi x
) xv ( k ) Qi
(1
T
) vi
Ri
ui
(k ) (1
) vi
xu ( k ) (1 (1
) xv ( k )
) vi
Berdasarkan Teorema (2.1.4), karena Qi simetrik semidefinit positif dan Ri simetrik definit positif, maka dapat difaktorisasi menjadi Qi 2 dan Ri 2 . Didefinisikan k 1
xu (k )
Qi Ak x0
Qi xu (k )
Ak
i 1
Ak
i 1
Bu (i)
i 0 k 1
xv (k ) Qi xv (k ) Qi Ak x0
Bv(i)
i 0
ui (k )
Riui (k )
vi (k )
Ri vi (k ),
maka Ji
u (1
)v
1N 2k
1 0
xu (k ) (1 ui
(1
T
) xv ( k ) QiQi )vi
T
Ri Ri
ui
xu ( k ) (1 (1
)vi
) xv ( k )
.
62
1N 2k 1N 2k
Qi xu (k ) (1
1
(1
Ri ui
0
)vi
T
T
Ji
u (1
Karena dan
)v
xu ( k ) ui
(1
1N 2k 1N 2k
(1
ui
)Qi xv (k )
) Ri vi ) xv ( k )
.
)vi
x
2
x T x yang memenuhi
y maka
1
xu (k ) (1
2
) xv (k )
(1
ui
)vi
2
0 1
xu (k )
2
) xv (k )
(1
ui
(1
2
)vi
.
0
(1 ) vi
x
(1
xu (k ) (1
Dengan menggunakan definisi norma Euclidean ketaksamaan segitiga x y
Qi xu (k ) (1
Ri ui
i
) xv ( k )
(1
ui
0
T
) Ri vi
xu (k ) (1
1
T
)Qi xv (k )
) xv ( k ) 2
ui
2
2
xu ( k )
2
(1
) xv ( k )
2
2
(1
) vi ,
maka Ji
u (1
)v
1N 2k 1N 2k
1
xu (k )
2
(1
) xv ( k )
(1
) xvT Qi xv
2
ui
2
2
(1
) vi
(1
)viT Ri vi
) xvT Qi xv
(1
0 1
xu T Qi xu
uiT Riui
0
1N 1 xu T Qi xu uiT Ri ui 2k 0 1N 1 T xu Qi xu uiT Riui 2k 0 J i (u ) (1 ) J i (v ),
1N 2k
1
(1
)viT Ri vi
0
(1
)
1N 2k
1
xvT Qi xv
viT Ri vi
0
yang artinya Ji konveks untuk setiap i = 1, 2, …, M. Berdasarkan lemma (2.6.4) dan akibat (2.6.5), maka himpunan semua solusi Pareto Kooperatif diberikan oleh J1 u (
, J2 u (
,..., J M u (
,
,
yang berkorespondensi dengan Pareto efficient yang diperoleh dari masalah minimisasi
63
M
u ( )
arg min u U
i
J i ( xi ( k ), ui ( k )),
i 1
dengan kendala x(k 1)
Ax(k ) B1u1 (k ) ... BM u M (k ), x(0)
x0 .
Berdasarkan Teorema (2.6.6) diatas, dapat ditentukan Pareto frontier dengan cara menyelesaikan masalah minimisasi (2.6.15) untuk beberapa nilai
.
Contoh (2.6.1) memberikan gambaran bagaimana permainan dinamis kooperatif diskrit diselesaikan sehingga diperoleh Pareto frontier.
Contoh 2.6.1 Diketahui terdapat 2 pemain dinyatakan oleh persamaan dinamik
x(k 1)
x(k ) u1 (k ) u2 (k ), x0 1.
Fungsi objektif yang akan diminimalkan masing-masing pemain dinyatakan oleh N 1
N 1
x 2 (k ) u2 2 (k ).
x 2 ( k ) u12 ( k ) , J 2 ( x, u2 )
J 1 ( x , u1 )
k 0
k 0
Dibentuk kombinasi linier konveks dari kedua fungsi objektif yaitu J ( x, u1 , u2 )
J1
(1
)J2
N 1
x 2 (k ) u12 (k )
) x 2 ( k ) u2 2 ( k )
(1
k 0 N 1
x 2 (k )
u12 (k ) (1
x 2 (k )
u1 (k ) u2 ( k )
)u 2 2 ( k )
k 0 N 1 k 0
T
0 0
(1
u1 (k ) . ) u2 ( k )
Masalah ini diselesaikan menggunakan waktu proses N = 2. Pareto efficient ditentukan menggunakan (2.6.15) kemudian disubstitusikan ke fungsi objektif J1 dan J2 untuk memperoleh solusi Pareto kooperatif untuk nilai berjalan dari 0,1 sampai 0,9 dengan langkah 0,01. Langkah-langkah ini termuat dalam program MATLAB pada Lampiran (2), sehingga menghasilkan grafik Pareto frontier yang diberikan oleh Gambar (2.6.1) berikut.
64
Gambar 2.6.1 Pareto frontier Contoh 2.6.1
Gambar (2.6.1) menampilkan grafik Pareto frontier Contoh (2.6.1). Sumbu datar menyatakan nilai fungsi objektif pemain ke-1 (J1), sedangkan sumbu tegak menyatakan nilai fungsi objektif pemain ke-2 (J2). Sehingga titik-titik pada Pareto frontier menyatakan solusi Pareto J1, J 2 .
2.6.2 Solusi Nash-Bargaining Berdasarkan Teorema (2.6.6), masalah minimisasi (2.6.5) mempunyai lebih dari satu solusi Pareto. Selanjutnya akan digunakan konsep Nash-Bargaining untuk menentukan solusi Pareto yang terbaik menurut para pemain. Untuk kondisi tanpa kerjasama (non-cooperative), maka masing-masing pemain dapat menentukan solusi minimax, yaitu solusi terbaik yang dapat diperoleh pada kondisi-kondisi terburuk yang dapat terjadi. Solusi ini disebut titik disagreement. Untuk pemain ke-i, titik disagreement dinotasikan sebagai di
J i (u d ) , dengan ud
ui d
arg min max J i (u(k )) ui
u
kendala : ui ( k ) u i (k ) dengan
u
(u1d , u2d ,..., uM d ) diperoleh dari masalah minimax
i
i
i
,
(2.6.16) i
,
(u1 , u2 ,..., ui 1, ui 1 ,..., uM )
dan
i
1
i 1
i 1
M
.
Gambar (2.6.2.a) menampilkan konsep Nash-Bargaining untuk 2 pemain. Daerah
65
S merupakan daerah fisibel untuk J yang kurang dari d, dimana d merupakan titik disagreement, dan kurva P adalah Pareto frontier.
(a)
(b)
Gambar 2.6.2 (a) Konsep Nash-Bargaining untuk 2 pemain (b) Solusi Nash-Bargaining N(S,d)
Gambar (2.6.2.b) menampilkan solusi Nash-Bargaining N(S,d). Secara geometris, solusi Nash-Bargaining adalah titik di tepi S, yaitu bagian dari Pareto frontier, yang menghasilkan persegi panjang terluas (A,N,B,d). Sehingga solusi Nash-Bargaining, N (S , d ) , ditentukan dengan cara memilih titik di S sehingga perkalian dari selisih Ji dan di bernilai maksimum, yaitu M
N (S , d ) arg max J S
kendala : J Untuk suatu
di
Ji
(2.6.17)
i 1
S, J
d.
, didefinisikan masalah optimisasi M
di
max log J S
Ji
i
(2.6.18)
i 1
kendala : J
S, J
d.
Berdasarkan Teorema (2.3.15), menyelesaikan masalah optimisasi (2.6.17) sama dengan menyelesaikan masalah optimisasi (2.6.18). Fungsi objektif (2.6.18) dapat ditulis secara ekuivalen menjadi
66
M
di
log
Ji
i
log d1
1
J1
d2
2
J2
... d M
JM
M
i 1
log d1 1
1
J1
log d1
log d2
J1
2
J2
log d 2
2
... log d M ...
J2
M
JM
log d M
M
JM
M i
log di
Ji .
i 1
Sehingga masalah optimisasi (2.6.18) dapat ditulis sebagai M
max J S
i
log di
Ji
(2.6.19)
i 1
kendala : J
S, J
d.
Baik pada solusi Pareto maupun solusi Nash-Bargaining, masalah optimisasinya memuat pembobotan fungsi objektif
yang nilainya
ditentukan secara bebas oleh para pemain. Untuk permainan statis noncooperative, Myerson (1997) memberikan cara untuk menentukan nilai yaitu
berdasarkan
kontribusi
masing-masing
pemain.
Myerson,
(1997)
mengilustrasikan cara ini dengan sebuah permainan “Divide the Dollars” dimana pemain 1 adalah seorang kepala keluarga dengan total 4 anggota keluarga, dan pemain 2 adalah seorang saja. Pemain 1 dan pemain 2 akan memperoleh sejumlah uang apabila strategi mereka berjumlah kurang dari atau sama dengan $100, dan $0 apabila strategi mereka berjumlah lebih dari $100. Myerson (1997) memberikan bobot pemain 1 adalah 1
2
5
1
4
5
, sedangkan pemain 2 diberi bobot
, hal ini dilakukan untuk menyamakan pendapatan perkapita mereka,
sehingga
1
2
1 . Berdasarkan ilustrasi ini, dapat diperumum untuk M
pemain, dan pembobotannya adalah
1
,
2
,...,
M
sedemikian sehingga
M i
1 . Konsep ini akan digunakan pada teknik kendali MPC terdistribusi
i 1
untuk menentukan pembobotan fungsi objektif masing-masing subsistem.
67
Masalah maksimisasi (2.6.19) digunakan untuk menentukan solusi NashBargaining. Contoh (2.6.2) berikut merupakan kelanjutan dari contoh (2.6.1). Pada contoh (2.6.1), diperoleh solusi Pareto sepanjang Pareto frontier. Contoh (2.6.2) berikut memberikan gambaran cara menentukan solusi Nash-Bargaining.
Contoh 2.6.2 Diberikan permainan yang memuat dua pemain pada contoh (2.6.1). Pemain ke-1 menentukan titik disagreement dengan masalah optimisasi
u1d
arg min max J1 (u (k )). u1
u2
Pemain ke-2 menentukan titik disagreement dengan masalah optimisasi
u2d Diperoleh
d1
arg min max J 2 (u (k )). u2
J 1 (u d )
u1
32 dan d 2
untuk permainan ini adalah d
J 2 (u d )
32 , sehingga titik disagreement
(32,32) .
Misalkan fungsi objektif pemain ke-1 diberikan bobot objektif pemain ke-2 diberikan bobot
2
1
3 dan fungsi 4
1 , maka masalah optimisasi Nash4
Bargaining diberikan oleh
3 1 max log 32 J1 log 32 J 2 J 4 4 kendala : J1 32 dan J 2 32.
Masalah optimisasi diatas menghasilkan solusi J1 18,3 dan J 2
20,0
yang merupakan solusi Nash-Bargaining dari permainan ini. Langkah-langkah ini termuat dalam program MATLAB pada Lampiran (3). Titik disagreement dan solusi Nash-Bargaining dapat diilustrasikan pada Gambar (2.6.3).
68
Gambar 2.6.3 Titik disagreement dan Solusi Nash-Bargaining Contoh (2.6.2)
Gambar (2.6.3) menampilkan Pareto frontier dan solusi Nash-Bargaining untuk masalah permainan contoh (2.6.1). Berdasarkan Gambar (2.6.3), dapat dilihat bahwa solusi Nash-Bargaining N (18,3;20) merupakan solusi Pareto yang yang menghasilkan persegi panjang terluas yang ditarik dari titik disagreement
d (32;32) .
BAB III KENDALI MODEL PREDIKTIF TERDISTRIBUSI BERBASIS TEORI PERMAINAN DINAMIS KOOPERATIF
Teknik kendali MPC terdistribusi merupakan sebuah metode untuk mendesain kendali pada systemwide yang berbasis MPC. Teknik kendali MPC diterapkan pada masing-masing subsistem dan kendali optimal ditentukan secara terdistribusi, yaitu besarnya kendali dihitung oleh semua subsistem secara bersama-sama melalui pertukaran informasi nilai kendali dengan meminimumkan fungsi objektif subsistem. Konsep ini dibentuk berdasarkan teori permainan dinamis kooperatif, dimana subsistem dianggap sebagai pemain yang akan meminimalkan fungsi objektifnya masing-masing dengan melakukan kerjasama (kooperatif). Teknik kendali MPC terdistribusi yang akan dibahas adalah tipe kendali state feedback yaitu kendali optimal ditentukan berdasarkan umpan balik state sistem.
3.1 Model Systemwide Plant Diberikan plant yang terdiri dari M subsistem. Didefinisikan
M
yang
menyatakan himpunan bilangan-bilangan integer 1, 2,3,..., M . Metode state space digunakan untuk menyatakan model matematika dari systemwide plant. Model systemwide yang akan dibahas terdiri dari tiga macam, yaitu model terdesentralisasi, model interaksi, dan model composite. Model terdesentralisasi memuat persamaan state space masing-masing subsistem secara independen, sedangkan model interaksi memuat persamaan state space yang menyatakan interaksi antar
subsistem.
Model terdesentralisasi dan model interaksi
digabungkan sehingga diperoleh model composite. Model composite inilah yang akan digunakan untuk menyusun teknik kendali MPC terdistribusi. Model composite dipilih karena memuat persamaan state space masing-masing subsistem sekaligus memuat interaksi antar subsistem.
69
70
3.1.1 Model Terdesentralisasi Model terdesentralisasi adalah model yang memperhatikan setiap subsistem secara independen, sehingga model ini tidak memuat interaksi antar subsistem (interaksi dianggap nol). Diberikan model terdesentralisasi untuk setiap subsistem i
yang dinyatakan oleh persamaan state space diskrit LTI yaitu
M
xii (k 1)
Aii xii (k ) Biiui (k ),
(3.1.1a)
yi (k ) Cii xii (k ), nii
dengan xii
mi
, ui
nii nii
Aii
dan output,
(3.1.1b) , dan yi nii mi
, Bii
zi
berturut-turut menyatakan state, input, zi nii
, Cii
merupakan matriks-matriks
konstan.
3.1.2 Model Interaksi Diberikan sebarang subsistem i terhadap subsistem j
M
M
. Efek dari interaksi subsistem i
, j i dinyatakan oleh persamaan state space LTI
diskrit
xij (k 1)
Aij xij (k ) Biju j ( k ),
(3.1.2a)
M
yi ( k )
(3.1.2b)
Cij xij ( k ), j 1
dengan
xij
subsistem j
nij
menyatakan state interaksi, u j
i , dan yi
memuat interaksi,
Aij
zi
mj
menyatakan input
menyatakan output subsistem ke-i yang sudah
nij nij
, Bij
nij m j
, Cij
zi nij
merupakan matriks-
matriks konstan.
3.1.3 Model Composite Model composite merupakan gabungan dari model terdesentralisasi dan model interaksi. Untuk setiap subsistem i, vektor state model terdesentralisasi xii ditambahkan dengan state yang ditimbulkan dari efek interaksi dengan semua subsistem yang lain. State model composite subsistem i
dinyatakan oleh
71
xi
xiT1 ,
, xiiT ,
T , xiM
T
ni
M
, ni
n . Model composite dinyatakan oleh
j 1 ij
persamaan state space diskrit LTI berikut xi ( k 1)
Ai xi ( k ) Biui ( k )
(3.1.3a)
Wij u j (k ), j i
yi (k ) Ci xi (k ),
(3.1.3b)
dengan
Ai1 Ai
Aii
, Bi
0
0
Bii , Wij 0
Bij , 0
AiM Ci
Ci1
Cii
CiM ,
merupakan matriks-matriks konstan. Diasumsikan state equilibrium dari model composite (3.1.3) adalah nol.
3.2 Model Prediksi Systemwide Plant Model yang digunakan untuk membentuk prediksi adalah model composite (3.1.3). Untuk sebarang subsistem i masukan pada waktu instant k berturut-turut
dinyatakan
j, j
oleh
M
, prediksi state dan prediksi
0 berdasarkan nilai awal pada saat k secara
xˆi k
j|k
ni
dan
uˆi k
j|k
mi
.
Berdasarkan pada MPC klasik, untuk menentukan state dan input prediksi diperlukan parameter panjang horison. Pada pembahasan ini, digunakan panjang horison prediksi berhingga Hp, horison window Hw = 1, dan horison kendali Hu = Hp. Prediksi state model composite untuk setiap subsitem i sepanjang horison prediksi Hp ditentukan dengan prosedur rekursif berikut.
xˆi ( k 1 k )
Ai xi ( k )
Bi uˆi ( k | k )
Wijuˆ j (k | k ), j i
xˆi ( k
2 k)
Ai xˆi ( k 1| k )
Bi uˆi ( k 1| k )
Wij uˆ j ( k 1| k ) j i
72
Biuˆi ( k | k )
Ai Ai xi (k )
Wijuˆ j (k | k ) j i
Bi uˆi (k 1 k )
Wij uˆ j (k 1 k ) j i
Ai 2 xi (k )
Ai Biuˆi (k | k ) Biuˆi (k 1| k )
ˆ AW i ij u j ( k | k )
Wij uˆ j (k 1| k ),
j i
xˆi (k 3 k )
j i
Ai xˆi (k
2 | k ) Biuˆi (k
Wijuˆ j (k
2 | k)
2 | k)
j i
Ai 2 xi (k ) Ai
Ai Biuˆi (k | k ) Biuˆi (k 1| k )
ˆ AW i ij u j ( k | k )
Wijuˆ j (k 1| k )
j i
j i
Biuˆi (k 2 k )
Wijuˆ j (k
2 k)
j i
Ai3 xi (k )
Ai 2 Biuˆi (k | k )
Ai 2Wij uˆ j (k | k )
ˆ AW 1| k ) i ij u j ( k
j i
xˆi (k
H p k)
Ai Biuˆi (k 1| k ) Biuˆi (k Wij uˆ j (k
j i
H
Ai p xi (k ) Ai
Ai
Hp 1
j i
Biuˆi (k | k )
Hp 1
Wijuˆ j (k | k )
Biuˆi (k Wijuˆ j (k
j i
H p 1| k )
H p 1| k ).
j i
Prediksi state sepanjang horizon Hp dapat ditulis sebagai
Ai
xˆi (k 1 k ) xˆi (k 2 k ) xˆi (k 3 k )
xˆi (k
Ai 2 Ai 3
H p k)
Ai
Ai
Hp 1
Bi
xi (k )
Hp
Bi Ai Bi Ai 2 Bi
0 Bi Ai Bi Ai
Hp 2
0 0 Bi Bi
Ai
Hp 3
0 0 0 0 Bi
2 | k)
Bi
uˆi (k k ) uˆi ( k 1 k ) uˆi (k 2 k ) uˆi (k
H p 1 k)
2 | k ),
73
Wij AW i ij
0 Wij
0 0
0 0
Ai 2Wij
AW i ij
Wij
0 0
j i
Ai
Hp 1
Wij
Ai
Hp 2
Wij
Ai
Hp 3
Wij
Wij
uˆ j (k k ) uˆ j (k 1 k ) uˆ j (k 2 k ) uˆ j (k
.
H p 1 k)
Misalkan
Ai
uˆi (k k ) ui
uˆi (k 1 k ) uˆi (k 2 k ) uˆi (k
Eij
Ai 2 Ai3
, fi
H p 1 k)
Ai
0 Wij
0 0
Ai 2Wij
AW i ij
Wij
Hp 1
Wij
Ai
Hp 2
Wij
Ai
0
0
0
Ai Bi
Bi
0
0
2
Ai Bi
Bi
0 , 0
Ai Bi Ai
Hp
Wij AW i ij
Ai
, Eii
Bi
Hp 1
Bi
Ai
Hp 2
Bi
Ai
Hp 3
Bi
Bi
0 0 0 , 0
Hp 3
Wij
Wij
maka state prediksi model composite sepanjang horizon Hp dapat ditulis menjadi
xˆi (k 1 k ) xˆi (k 2 k ) xˆi (k 3 k )
fi xi (k ) Eiiui (k )
Eiju j (k ).
(3.2.1)
j i
xˆi (k
H p k)
Teknik kendali MPC terdistribusi akan diselesaikan secara numerik, sehingga algoritma MPC terdistribusi yang akan dibentuk memuat iterasi. Misalkan q menyatakan iterasi algoritma ini. Sehingga state prediksi (3.2.1) dituliskan sebagai
74
xˆi (k 1 k ) xˆi (k 2 k ) xˆi (k 3 k )
Eij u jq 1 (k ).
fi xi (k ) Eiiui (k )
(3.2.2)
j i
xˆi (k
H p k)
3.3 Fungsi Objektif Systemwide Plant Stage cost, yaitu cost pada setiap waktu instan k dari masing-masing subsistem, didefinisikan oleh 1 xi (k )T Qi xi ( k ) ui ( k )T Riui ( k ) , 2
Li xi ( k ), ui ( k )
dengan Qi
0, dan Ri
(3.3.1)
0 adalah matriks pembobotan yang simetris. Fungsi cost
subsistem i untuk horison berhingga Hp didefinisikan oleh Hp 1 i
xˆi (k ), uˆi (k ); xi (k )
Li xˆi (k 1
j | k ), uˆi (k
j | k) .
(3.3.2)
j 0
Misalkan
i
diag Qi (0),..., Qi ( H p 1) ,
i
diag Ri (0), Ri (1),..., Ri ( H p 1) ,
state prediksi model composite (3.2.2) disubstitusikan ke fungsi objektif (3.3.2), sehingga diperoleh
i
xˆi ( k ), uˆi (k ); xi (k ) Hp 1
j 0
1 2
1 T xˆi (k 1 2
j | k )Qi xˆi (k 1
j | k ) uˆiT ( k
j | k ) Riuˆi (k
T
f i xi ( k )
Eii ui ( k )
Eij u j
q 1
(k )
i
j i
Eij u j q 1 ( k )
fi xi ( k ) Eiiui ( k ) j i
ui ( k )T
u (k )
i i
j | k)
75
1 2
xi (k )T fi T
ui (k )T Eii T
u j q 1 (k )T Eij T j i
f x (k )
i i i
i
Eii ui (k )
j
Eij u j q 1 (k )
ui (k )T
u (k )
i i
j i
1 xi (k )T fi T 2
f x (k )
i i i
i
Eiiui (k )
j
Eiju j q 1 (k )
j i
ui (k )T Eii T
f x (k )
i i i
i
Eii ui (k )
j
Eij u j q 1 (k )
j i
u j q 1 (k )T Eij T
f x (k )
i i i
i
Eii ui (k )
j
j i
j i T
ui ( k )
1 2
Eij u j q 1 (k )
u (k )
i i
xi (k )T fi T
f x (k ) xi (k )T fiT
i i i
i
xi (k )T f iT
Eii ui (k )
j
Eij u j q 1 (k )
j i
ui (k )T Eii T
f x (k ) ui (k )T EiiT
i i i
i
ui (k )T EiiT
Eiiui (k )
j
Eij u j q 1 (k )
j i
u j q 1 (k )T Eij T
j
u j q 1 (k )T Eij T
f i xi (k )
j i
u j q 1 (k )T Eij T j i
Eii ui (k )
j
j
Eii ui (k )
Eij u j q 1 (k )
j i
ui (k )T 1 ui ( k )T 2
u (k )
i i
u (k ) ui (k )T Eii T
i i
xi ( k )T f iT
i
i
Eii ui (k )
Eii ui (k ) ui (k )T EiiT
ui ( k )T Eii T
j
Eij u j q 1 (k )
j i
f x (k )
i i i
u j q 1 (k )T Eij T j i
1 xi ( k )T fi T 2
xi ( k )T fi T
f x (k )
i i i
u j q 1 (k )T Eij T j i
j
j i
j
Eij u j q 1 (k )
j i
j
u j q 1 ( k )T Eij T
f i xi ( k ) j i
j j i
Eij u j q 1 ( k ) .
76
Karena 1 xi ( k )T fi T 2
xi (k )T f iT
f x (k )
i i i
j
Eiju j q 1 (k )
j i
u j q 1 ( k )T Eij T
u j q 1 ( k )T Eij T
f x (k )
j i i
j i
j i
j
Eij u j q 1 ( k )
j i
konstanta
dan penambahan konstanta pada fungsi objektif suatu masalah optimisasi tidak mempengaruhi titik optimal, maka konstanta tersebut dapat dihilangkan dari fungsi objektif, sehingga (3.3.2) dapat ditulis menjadi i
1 ui ( k )T 2
xˆi ( k ), uˆi ( k ); xi ( k )
xi ( k )T fi T
u ( k ) ui ( k )T Eii T
i i
i
i
Eii ui ( k ) ui (k )T Eii T
ui (k )T Eii T
j
Eij u j q 1 (k )
j i
Eiiui ( k ) f x (k )
i i i
u j q 1 (k )T Eij T
j
Eiiui (k ) .
j i
(3.3.3) 3.4 Feasible-Cooperation MPC Diberikan kendali
ui (k
j)
mi
menyatakan himpunan kendali fisibel untuk aksi
i
untuk
Diasumsikan himpunan
i
j
0,1,..., H p 1
masing-masing
subsistem
i.
merupakan himpunan tidak kosong dan konveks yang M
memuat titik asal di interiornya ( 0
i ). Misalkan
i
i
. Untuk
i 1
menyederhanakan, didefinisikan
i
xˆi ( k ), uˆi (k ); xi ( k )
i
ui (k ) . Secara umum,
masalah MPC untuk systemwide didefinisikan oleh masalah optimisasi M
min u (k )
i
(u ( k ))
(3.4.1)
i 1
kendala : ui (k )
i
, i 1, 2,..., M .
Masalah optimisasi MPC (3.4.1) merupakan masalah optimisasi MPC tersentralisasi karena variabel optimisasinya adalah
u ( k ) . Sebagaimana
77
dijelaskan pada latar belakang penelitian ini, MPC tersentralisasi sangat sukar untuk diselesaikan terutama apabila systemwide memuat subsistem-subsistem yang relatif banyak. Sehingga, dilakukan solusi pendekatan dengan cara menyelesaikannya secara terdistribusi. Ada beberapa cara pendekatan untuk menyelesaikan masalah optimisasi MPC secara terdistribusi, salah satunya adalah dengan memberikan pembobotan pada masing-masing fungsi objektif, yaitu dengan membentuk fungsi objektif
i
yang mengukur pengaruh kontrol lokal
terhadap performansi systemwide secara keseluruhan. Ada banyak pilihan yang dapat dibuat, bentuk yang paling sederhana adalah membentuk kombinasi konveks tegas dari fungsi objektif semua subsistem. Misalkan
ui (k ), u i (k )
r
dengan
i
u (k ) ,
(3.4.2)
u1T (k ),..., u iT 1 ( k ), u iT 1 ( k ),..., u MT ( k )
u i (k )
u i (k )
r
dengan
i
...
1
i 1
i 1
...
T
untuk M
ui (k )
i
dan
. State prediksi dari
model composite sudah disubstitusikan ke fungsi objektif, sehingga i
xˆi , uˆi ; xi (k )
i
u1 , u2 ,..., uM ; xi (k ) .
(3.4.3)
M
Untuk setiap subsistem i, wr
0,
wr
1 , fungsi objektif yang baru
r 1
didefinisikan oleh M i
u (k )
wr
r
(3.4.4)
ui ( k ), u i ( k ) .
r 1
Untuk menentukan kendali optimal, dilakukan proses negosiasi yang memuat iterasi dan pertukaran nilai kendali antar subsistem yang dilakukan pada waktu
sampel.
menyatakan u qi (k ) u jq ( k ), j
Misalkan
indeks
uiq (k )
iterasi
saat
uiqT ( k ),..., uiqT ( k
T
,
( q 1,2,..., qmax ).
k
u1q*T (k ),..., uiq*1T (k ), uiq*1T ( k ),..., uMq*T ( k )
H p)
T
dengan
q
Misalkan
menyatakan nilai optimal dari
1,..., i 1, i 1,..., M . Trayektori masukan optimal uiq* ( k ) diperoleh
dari solusi masalah optimisasi FC-MPC yang didefinisikan oleh
78
1
(i) : Feasible-Cooperation MPC M
ui ( q ) (k ) arg min uˆi
kendala : uˆi (k
wr
r
uiq 1 (k ),..., uiq 11 (k ), ui (k ), uiq 11 (k ),..., uMq 1 (k ); xr ( k )
r 1
j | k)
i,
uˆi (k
j | k ) 0,
xr (k )
xˆr (k ), r
0
j
j
H p, M
H p 1,
.
Fungsi objektif dari optimisasi FC-MPC dapat disederhanakan sebagai M
wi
xˆi (k ), uˆi (k ); xi (k )
i
i 1
1 2
M
wi ui (k )T
u (k ) ui (k )T Eii T
i i
i
Eiiui (k )
i 1
xi (k )T fiT
i
T
ui (k ) Eii
Eiiui (k ) ui (k )T Eii T T j
Eiju j
q 1
f x (k )
i i i
(k )
j i
u j q 1 (k )T EijT
j
Eiiui (k )
j i
1 M wi ui (k )T i ui (k ) ui (k )T EiiT i Eii ui (k ) 2i1 1 M wi xi (k )T fi T i Eii ui (k ) ui (k )T EiiT i fi xi (k ) 2i1 1 2
(3.4.5)
M
ui (k )T Eii T
wi i 1
j
Eij u j q 1 (k )
u j q 1 (k )T Eij T
j i
j
Eii ui (k ) .
j i
Suku pertama (3.4.5) dapat ditulis sebagai :
1 2
M
wi ui (k )T
T T iui ( k ) wi ui ( k ) Eii
i Eii ui ( k )
i 1
1 ui (k )T 2
1 2
M
ui (k )T wi
i
i 1
M
wi
i
wi Eii T
i
Eii ui (k )
i 1
1 ui (k )T wi 2
M i
wi Eii T
i Eii
wj j i
j
w j E ji T
j
E ji ui (k ).
wi Eii T
i
Eii ui (k )
79
Suku kedua (3.4.5) dapat ditulis sebagai :
1 2
M
wi xi (k )T f i T
i
Eii ui (k ) wi xi (k )T f iT
i
Eiiui ( k )
i 1 M
1 2
2 wi xi ( k )T fi T
i
Eii ui (k )
i 1
M
wi xi (k )T fi T
i
Eii ui (k )
i 1 M
wi xi (k )T fi T
i
Eii ui (k )
i 1 T
M T
wi Eii f i
T
x (k )
w j E ji f j
i i
j
x j (k )
ui (k ).
j 1
Suku ketiga (3.4.6) dapat ditulis sebagai : 1 2
M
ui (k )T EiiT
wi i 1
j
Eij u j q 1 (k )
u j q 1 (k )T Eij T
j i
1 2
j
Eiiui (k )
j i
M
u j q 1 (k )T Eij T
wi 2 i 1
j
Eii ui (k )
j i
M
u j q 1 (k )T Eij T
wi i 1
j
Eii ui (k )
j i
M
T
M
wl Eli
T l
Elju j
q 1
(k )
ui (k ).
j 1, j i l 1
Berdasarkan penjabaran diatas, (3.4.5) menjadi M
wi
ui (k ), u i (k ); xi (k )
i
i 1 M
1 ui (k )T wi 2
wi EiiT
i
i Eii
wj j i
T
M T
wi Eii fi
T
x (k )
w j E ji f j
i i
j
x j (k )
j 1 M
T
M
wl Eli j 1, j i l 1
w j E ji T
j
T l
Elju j
q 1
(k )
ui (k ).
ui (k )
j
E ji ui (k )
80
Berdasarkan penyederhanaan diatas, masalah optimisasi FC-MPC untuk setiap subsistem i untuk setiap iterasi q dapat ditulis menjadi
2
(i) : Feasible-Cooperation MPC i
T
M
1 q T min ui ( k ) ui q ( k ) 2
q i i
u (k )
ri ( k )
ij
u
q 1 j
(k )
ui q ( k )
j 1, j i
q
kendala : ui ( k )
i
dengan i
diag Qi (0),..., Qi ( H p 1) ,
i
wi
diag Ri (0), Ri (1),..., Ri ( H p 1) ,
i
M
wi EiiT
i
i
M
w j E ji T
Eii
j
E ji ,
wl Eli T
ij
j i
l
Elj
l 1
M
wi Eii T
ri (k )
w j E ji T
i f i xi ( k )
j
f j x j (k )
j i
Ai
Eii
Bi
0
Ai Bi
Bi
Hp 1
Ai
0
0
Wij
0
0
AW i ij
Wij
, Eij
0
Hp 1
Bi
Bi
0 0
Ai
Ai 2 , fi
Ai 3
Wij
Wij
A Bentuk fungsi objektif pada optimisasi 1 2
.
2
Hp
dapat dibentuk menjadi
M T
T i
uiq (k ),
dengan
T
ri (k )
ij
u jq 1 (k ) yang merupakan
j 1, j 1
bentuk umum pemrograman kuadratik pada Subbab (2.3.4). Bentuk ini dapat mempermudah dalam penyelesaian optimisasi pada MATLAB. Berdasarkan masalah optimisasi FC-MPC, dapat dirumuskan algoritma untuk menentukan solusi optimal uiq * . Flow chart algoritma ini termuat pada Lampiran (5).
Algoritma FC-MPC : Diberikan ui0 , xi ( k ),
q
1,
i
K ,K
i
1
0,
i
0, i
M
, qmax ( k )
0,
0
81
while
for some i
i
i
do
and q
M
qmax (k )
M
uiq
arg( i )
uiq
wiuiq
(1 wi )uiq
1
Transmit u iq ke setiap subsistem terkoneksi j i uiq
i
uiq
1
end do q
q 1
end while
Pada algoritma FC-MPC, trayektori state untuk setiap subsistem i pada saat iterasi ke-q diperoleh sebagai xiq
xiq u1q , u2q ,..., uMq ; x(k ) . Iterasi dibatasi
hanya sampai qmax (k ) , tetapi iterasi dapat berhenti sebelum qmax (k ) apabila error semua subsistem sudah kurang dari toleransi error yang diberikan.
3.5 Skema Nash-Bargaining MPC Terdistribusi Sebuah permainan dapat dituliskan sebagai tupel dengan N
1, 2,..., M
N,
menyatakan himpunan para pemain,
i i N
i
,
i i N
,
menyatakan
himpunan berhingga dari strategi yang dapat dilakukan oleh pemain ke-i, dan i
:
1
2
M
menyatakan fungsi payoff pemain ke-i (Akira, 2005).
Masalah MPC terdistribusi dapat dituliskan sebagai tupel
M
,
i i
, M
i i
, M
sehingga MPC terdistribusi dapat dipandang sebuah permainan dimana para pemainnya adalah subsistem, strategi dari pemain adalah masukan kendali, dan payoff diberikan oleh nilai fungsi cost subsistem. Selanjutnya, kendali optimal MPC terdistribusi akan ditentukan berdasarkan solusi Nash-Bargaining seperti pada permainan dinamis kooperatif diskrit yang diberikan pada Subbab (2.6).
82
3.5.1 Solusi Nash-Bargaining MPC Terdistribusi Berdasarkan teori permainan dinamis kooperatif pada Subbab (2.6), solusi Nash-Bargaining dari masalah MPC terdistribusi didefinisikan sebagai masalah optimisasi
3
(i) : Nash-Bargaining MPC M
di (k )
max uˆ( k )
i
(u (k ))
i 1
kendala :
i
(u (k ))
ui ( k )
di ( k ) i
(3.5.1)
, i 1,.., M .
dengan di (k ) menyatakan titik disagreement subsistem i yang diberikan oleh di (k )
i
u d (k ) dan u d (k ) ditentukan dengan mencari solusi minimax, yaitu
4
(i) : Titik Disagreement MPC
min max ui ( k ) u i ( k )
u (k )
i
kendala : ui ( k ) u i (k )
Berdasarkan
(3.5.2)
i i
.
optimisasi
(2.6.17)-(2.6.19),
menyelesaikan
masalah
optimisasi (3.5.1) sama dengan menyelesaikan masalah optimisasi (3.5.3) berikut.
M
max u (k )
wi log di (k )
i
(u ( k ))
i 1
kendala :
i
(u ( k ))
ui ( k )
di ( k ) i
(3.5.3)
, i 1,.., M .
Berdasarkan hasil diatas, maka masalah Nash-Bargaining MPC dapat diselesaikan secara terdistribusi sebagai masalah optimisasi (3.5.4) berikut.
83
5
(i) : Nash-Bargaining MPC M
max ui ( k ) r 1
wr log d r (k )
kendala :
r
r
(ui (k ), u i (k ))
(ui ( k ), u i ( k )
ui ( k )
i
(3.5.4)
d r (k )
.
3.5.2 Model Negosiasi Nash-Bargaining MPC Model negosiasi Nash-Bargaining terdiri dari langkah-langkah yang akan memberikan solusi Nash-Bargaining secara numerik untuk subsistem i setiap waktu k, pada setiap iterasi
M
. Pada
q, q 1,2,..., qmax , algoritma untuk
menyelesaikan Nash-Bargaining MPC adalah sebagai berikut. Flow chart algoritma ini termuat pada Lampiran (6).
x(k ) , semua subsistem menghitung titik
1. Diberikan kondisi awal
d i (k )
disagreement
i
u d (k )
masing-masing
secara
terpisah
berdasarkan masalah optimisasi
ui d (k ) arg min max ui ( k ) u i ( k )
kendala : ui ( k )
i
u (k ) (3.5.5)
i
u i (k )
i
2. Setiap subsistem mengirimkan titik disagreement ke subsistem lain 3. Setiap subsistem menyelesaikan masalah optimisasi M i
max q
ui ( k ) r 1
kendala :
wr log dr (k ) r
(uiq ( k ), u qi 1 )
r
(uiq (k ), u qi 1 (k ))
d r ( k ), r
(3.5.6)
1,.., M
q
ui ( k )
i
4. Jika (3.5.6) fisibel, maka uiq* (k ) arg memenuhi
r
(uiq* (k ), u qi 1 ) dr (k ) ),
i
menjadi solusi optimal ( uiq* (k ) kemudian
subsistem
memperbarui aksi kendali masing-masing sebagai kombinasi konveks
ke-i
84
uiq ( k )
wi uiq ( k ) (1 wi )uiq 1 ( k ) . Jika tidak fisibel, subsistem ke-i
memperbarui aksi kendali uiq ( k )
wi uid ( k ) (1 wi )uiq 1 ( k ) . Pada langkah
ini, untuk q = 1, maka uid ( k ) dijadikan kondisi awal oleh subsistem i untuk menyelesaikan optimisasi (3.5.6). Untuk q 1 , uiq 1 (k ) digunakan sebagai kondisi awal (3.5.6). 5. Setiap
pemain
mengirim
uiq ( k ) uiq 1 ( k )
,(
aksi
kendali
ke
pemain
0) untuk semua i, atau jika q
lain.
Jika
qmax , atau jika
waktu maksimum untuk menghitung kendali optimal sudah tercapai, elemen pertama dari barisan kendali uiq ( k ) diaplikasikan ke sistem dan setiap subsistem kembali ke langkah 1. Selain itu, setiap subsistem kembali ke langkah 3.
Pada saat k 1, kondisi awal subsistem i untuk menyelesaikan (3.5.6) ditentukan ui0 ( k 1)
dengan
pergeseran
uiqend *T ( k 1| k ),..., u Mqend *T ( k
H p | k ), 0
barisan T
,
dengan
kendali
uiqend * (k
j | k)
menyatakan nilai kendali optimal untuk subsistem i saat iterasi q end saat waktu k
j yang diberikan kondisi awal saat k. Pada model negosiasi ini, subsistem
hanya mengkomunikasikan titik disagreement pada setiap k, dan subsistem hanya mengirim barisan kendali pada setiap iterasi q.
3.6 Sifat-sifat MPC Terdistribusi Berdasarkan algoritma FC-MPC dan model negosiasi Nash-Bargaining MPC dapat diturunkan beberapa sifat berikut.
Proposisi 3.6.1 Masalah optimisasi Nash-Bargaining MPC (3.5.4) merupakan masalah optimisasi konkaf.
Bukti : Masalah optimisasi (3.5.4) dapat ditulis kembali menjadi
85
M
min
wr log d r (k )
ui ( k )
r
kendala :
r
(ui ( k ), u i ( k )
ui ( k )
i
0 dan
di (k )
, i 1,.., M
Menurut Teorema (2.3.5), fungsi f ( x)
(ui (k ), u i ( k )
r 1
log f ( x) merupakan fungsi konveks jika
f ( x) adalah fungsi konveks. Fungsi log d r (k )
memenuhi kendala dr (k )
r
(ui (k ), u i (k )) 0 dan
merupakan fungsi konveks karena dr tertentu dan
log d r (k )
Sehingga
r
d r (k ) r
r
(ui (k ), u i (k )
r
(ui (k ), u i (k )
(ui (k ), u i (k ) konveks.
(ui (k ), u i (k ) konveks. Berdasarkan Teorema (2.3.4),
M
fungsi
wr log d r ( k )
r
(ui ( k ), u i ( k ) konveks karena merupakan kombinasi
r 1
konveks dari fungsi-fungsi konveks. Sehingga masalah optimisasi (3.5.4) adalah masalah optimisasi konkaf.
Karena masalah optimisasi (3.5.4) merupakan masalah optimsasi konkaf, maka berdasarkan Teorema (2.3.13), apabila terdapat pemaksimal lokal, maka pemaksimal lokal tersebut adalah pemaksimal global.
3.6.1 Fisibilitas Algoritma MPC Terdistribusi Untuk menganalisis fisibilitas iterasi dari algoritma MPC terdistribusi, diberikan preposisi (3.6.2) berikut.
Proposisi 3.6.2 Barisan kendali untuk subsistem i yang diperoleh dari algoritma MPC terdistribusi pada waktu k dengan kondisi fisibel awal adalah fisibel untuk waktu k , k 1, k
2,... , dengan i 1,.., M .
Bukti : (1) Algoritma FC-MPC. Pada setiap iterasi q saat k, diketahui uiq
arg( i ) dengan
i
adalah optimisasi fungsi konveks untuk setiap i. Pada
setiap iterasi, update kendali subsistem i adalah uiq
wi uiq
(1 wi )uiq 1 . Karena
86
adalah himpunan konveks, maka sebarang kombinasi konveks dari uiq dan
i
uiq
1
termuat di
i
. Sehingga barisan kendali untuk subsistem i yang di hasilkan
dari algoritma FC-MPC adalah fisibel untuk semua iterasi q 1,2,..., qend , saat k. (2) Model negosiasi Nash-Bargaining. Karena konveks, dan karena
i
u qi ( k )
i
i
dan u id( k )
dan
r
(ui (k ), u i (k ) adalah fungsi
keduanya himpunan konveks, maka
i
ada untuk semua iterasi q 1,..., qend , untuk semua
i 1,2,..., M . Misalkan ui0 ( k )
menyatakan solusi inisial saat k. pada model
i
negosiasi Nash-Bargaining, saat iterasi q = 1 subsistem i memperbarui aksi kendalinya sebagai ui1 (k )
wiuid ( q ) (k ) (1 wi )ui0 (k ) . Karena ui0 (k )
dan
i
konveks, maka sebarang kombinasi konveks dari ui0 (k ) dan uid ( k ) termuat di
i
. Sehingga ui1 (k )
i
untuk i 1,2,..., M . Jika saat iterasi q = 1 saat k subsistem
i memperbarui aksi kendali sebagai uiq (k )
ui0 (k )
i
, ui ( k )
i
, dan
i
Saat iterasi q = 2, ui1 (k )
i
wiuiq (k ) (1 wi )ui0 (k ) . Karena
konveks, maka sebarang kombinasi konveks
dari ui0 (k ) dan ui ( k ) termuat di
q
i
i
. Sehingga ui1 (k )
, maka ui2 (k )
i
i
untuk i 1,2,..., M .
. Begitu pula untuk barisan iterasi
2,3,..., qend . Sehingga barisan kendali untuk subsistem i yang di hasilkan dari
model
negosiasi
Nash-Bargaining
adalah
fisibel
untuk
semua
iterasi
q 1,2,..., qend , saat k. Untuk waktu k 1 , kondisi awal subsistem i untuk optimisasi diberikan oleh pergeseran barisan kendali ui0 ( k 1)
u iqend *T ( k 1| k ),..., u iqend *T ( k
H p | k ), 0
T
.
Hp
Tinjau kembali
i
i
dan
0
i
. Karena uiqend * (k )
fisibel, maka
j 1
uiqend (k
j)
i
untuk
j 1,2,..., H p . Akibatnya
ui0 (k 1)
i
. Sehingga
barisan masukan kendali yang dihasilkan oleh algoritma FC-MPC dan model negosiasi Nash-Bargaining MPC untuk subsistem i termuat di
i
untuk waktu
87
k 1. Analog untuk waktu k
3,... . Jadi, barisan masukan kendali yang
2, k
dihasilkan oleh algoritma MPC terdistribusi untuk subsistem i adalah fisibel untuk waktu k , k 1, k
2,... .
3.6.2 Kestabilan Sistem Lingkar Tertutup MPC Terdistribusi Misalkan pada saat k, iterasi algoritma MPC terdistribusi dinotasikan oleh q(k) dan berakhir pada iterasi ke p. Misalkan untuk setiap i ui p
(k )
T
uip
T
T
( k ),0 , uip
( k ),1 ,..., uip
M
,
( k ), H p 1
T
(3.6.1)
menyatakan solusi dari algoritma MPC terdistribusi setelah p iterasi dengan (k )
x1 ( k ), x2 ( k ),..., xM ( k )
T
menyatakan state model composite saat k.
Berdasarkan konsep MPC, kendali yang diaplikasikan ke subsistem i adalah
ui (k ) uip
(k ),0 . Nilai kendali ini diaplikasikan ke sistem sebagai state
feedback. Proposisi (3.6.2) menjamin eksistensi trayektori masukan fisibel untuk semua subsistem saat k
k
0 dan q (0)
0 , kemudian eksis juga untuk semua
0 . Salah satu pilihan trivial untuk trayektori masukan fisibel awal saat k
adalah
ui (k
j | k ) 0, j 0, i
M
0
. Hal ini mungkin dilakukan karena
berdasarkan asumsi awal, himpunan kendali admissible
i
,i
M
masing-
masing tidak kosong dan memuat titik asal didalamnya. Untuk dapat menjamin kestabilan, teknik kendali MPC terdistribusi bekerja dibawah kedua asumsi berikut.
Asumsi 3.6.1. Semua model interaksi adalah stabil, yaitu Aij stabil untuk i Asumsi 3.6.2. Untuk setiap i
M
, pasangan Aii , Bii bersifat stabilizable.
Asumsi 3.6.3. Untuk setiap i
M
, Qi (0)
Ri
0
Ri (0)
Ri (1) ...
Ri ( H p 1)
Qi (1) ... Qi ( H p 1)
Qi
j.
0 dan
88
Asumsi (3.6.1) diperlukan untuk menjamin kestabilan sistem lingkar tertutup MPC terdistribusi dibawah state feedback. Asumsi (3.6.2) merupakan syarat perlu suatu sistem dapat distabilkan. Asumsi (3.6.3) merupakan asumsi yang
diperlukan pada pendefinisian fungsi cost yang akan dijadikan fungsi
Lyapunov.
u1q , u 2q ,..., u Mq ; ( k )
Lemma 3.6.4 Barisan nilai fungsi objektif
yang
dihasilkan dari algoritma MPC terdistribusi adalah fungsi tidak naik terhadap iterasi q.
Bukti : Karena uiq * peminimal
maka
i
u1q 1 ,..., uiq 11 , ui*( q ) , uiq 11,..., uMq 1; ( k )
u1q 1, u 2q 1,..., u Mq 1; ( k )
Tinjau kembali bahwa M
u1q , u2q ,..., u Mq ; ( k )
wr
u1q , u2q ,..., u Mq ; ( k )
r
r 1
w1
u1q , u2q ,..., uMq ; (k )
1
wM
...
u1q , u2q ,..., u Mq ; ( k )
M
Berdasarkan pengambilan kendali baru uiq
wiuiq*
(1 wi )uiq 1 , maka
M
u1q , u2q ,..., uMq ; ( k )
wr
w1u1q *
r
(1 w1 )u1q 1 , w2u2q
(1 w2 )u2q 1,...,
r 1
wM u Mq *
w1
1
(1 wM )uMq 1; ( k )
w1u1q*
(1 w1 )u1q 1, w2u2q*
wM u Mq* wM
Karena
M
(1 w2 )u2q 1 ,...,
(1 wM )uMq 1; ( k )
w1u1q*
...
(1 w1 )u1q 1 , w2u2q*
wM u Mq*
(1 wM )uMq 1; (k )
w1
u1q* , u2q 1 ,..., uMq 1; ( k )
(1 w2 )u2q 1,...,
konveks, maka
u1q , u2q ,..., uMq ; (k )
wM
1
M
u1q 1 , u2q 1 ,..., uMq* ; (k )
...
89
M
wr
u1q 1 , u2q 1 ,..., urq 11 , ur*( q ) , urq 1 ,..., uMq 1; ( k )
r
r 1
u1q 1 , u2q 1 ,..., uMq 1 ; (k )
Karena
u1q , u2q ,..., u Mq ; ( k )
barisan
u1q , u2q ,..., uMq ; ( k )
u1q 1 , u2q 1,..., u Mq 1; ( k )
tidak naik.
T
Pada saat k, misalkan ui0 (0) sehingga JHp
(0)
barisan
untuk setiap q, maka
0,0,...,0 , i ui0 (0)
kendali
u10 (0),..., u M0 (0); (0)
M
. Karena 0 int( ) ,
fisibel.
Didefinisikan
menyatakan nilai fungsi cost dengan (0) . Misalkan saat k, untuk
kendali awal nol dan state awal
i
M
,
didefinisikan
ui0 (k )T dan
uiq ( k
1)
T
(k 1),1 ,..., uiq ( k
u10 ( k ), u20 ( k ),..., u M0 ( k )
T
1)
(k 1), H p 1 ,0,0,...
(3.6.2)
menyatakan himpunan kendali fisibel yang
bersesuaian dengan fungsi cost J H0 p
u10 ( k ), u20 ( k ),..., uM0 ( k ); ( k ) .
(k )
Sehingga nilai fungsi cost setelah iterasi q(k) dinotasikan oleh J Hq (pk )
u1q ( k ) (k ),..., u Mq ( k ) ( k ); ( k ) .
(k )
Lemma 3.6.5. Saat k nilai kendali ui (k
(3.6.3)
0 , misalkan algoritma MPC terdistribusi diawali dengan
j | k ) 0, j 0,
i
M
. Jika untuk semua k
0 , setiap MPC
terdistribusi diawali dengan strategi ui0 ( k )T
uiq ( k
1)
T
( k 1),1 ,..., uiq ( k
1)
( k 1), H p 1
T
,
maka k 1 M
J Hq (pk )
(k )
J H0 p
(k )
JH p
wi Li xi ( j ),0
(0) j 0 i 1
untuk
q(k )
0 dan semua k
0.
JH p
(0)
(3.6.4)
90
Bukti : Akan dibuktikan menggunakan induksi. Saat k = 0, algoritma MPC terdistribusi diawali dengan akibatnya J H0 p
(0)
(0) dan J Hq (0) p
JHp
ui (k
kendali awal (0)
j | k ) 0, j 0,
J H0 p
(0)
i
M
,
(0) .
JHp
sehingga (3.6.4) dipenuhi untuk k = 0. Untuk k = 1, didapat M
J Hq (1) p
(1)
J H0 p
J H0 p
(1)
wi Li xi (0), uiq (0) (0)
(0) i 1
M
JHp
wi Li xi (0),0
(0) i 1
JHp
(0) .
Sehingga (3.6.4) dipenuhi untuk k = 1. Diasumsikan benar untuk k > 1. Maka untuk k + 1, M
J Hq (pk
1)
(k 1)
J H0 p
J Hq (pk )
(k 1)
wi Li xi (k ), uiq ( k ) ( k )
(k ) i 1
M
J Hq (pk )
(k )
wi Li xi (k ),0 i 1 k
JHp
M
(0)
wi Li xi ( j ),0 j 0 i 1
JHp
(0) .
Terbukti untuk k + 1. Jadi, terbukti benar untuk semua k
0.
Lemma diatas akan digunakan untuk menunjukkan kestabilan dalam arti Lyapunov. Untuk membuktikan kestabilan, digunakan teorema kestabilan Lyapunov (2.3.6) pada Subbab (2.3).
Teorema 3.6.6. Diberikan algoritma MPC terdistribusi. Misalkan asumsi-asumsi (3.6.1)-(3.6.3) dipenuhi. Maka titik ekuilibrium 0 adalah stabil asimtotik untuk sistem lingkar tertutup M
xi (k 1)
Ai xi (k ) Biuiq ( k )
Wiju qj ( k )
(k ),0 j i
untuk semua
(0)
n
dan q(k ) 1,2,..., qmax (k )
(k ),0 , i
M
(3.6.5)
91
Bukti xiq
:
Kandidat
fungsi
xiq ( ,1)T , xiq ( , 2)T ,...
T
Lyapunov
J Hq (pk )
adalah
(k ) .
Misalkan
menyatakan trayektori state subsistem i
M
yang
dihasilkan oleh trayektori masukan u1q ( ),..., uMq ( ) setelah iterasi ke-q dengan state awal
. Fungsi cost diberikan oleh
J Hq (pk )
u1q ( k ) (k ),..., u Mq ( k ) ( k ); ( k )
(k )
Hp 1
M
Li xiq ( k )
wi i 1
l 0
M
Hp 1
wi i 1
l 0
(k ), l , uiq ( k ) T
q(k ) (k ), l Qi xiq ( k ) (k ), l 1 xi 2 u q ( k ) ( k ), l T R u q ( k ) ( k ), l i i i
dimana Qi , Ri kesemuanya definit positif dan xiq ( k )
(2.2.6), J Hq (pk )
(k )
(k ) kontinu di
J Hq (pk
J Hq p
1)
(k )
(k )
xi (k ), i
M
.
merupakan fungsi
berbentuk kuadratik, maka menurut Teorema n
. Berdasarkan lemma (3.6.5), diperoleh
J Hq (pk )
(k 1)
J Hq (pk )
(k ),0
J Hq (pk )
Pertama, akan ditunjukkan bahwa Lyapunov. Karena J Hq (pk )
(k ), l
J Hq (pk )
(k )
(k )
0 untuk semua
n
(k )
. Sehingga, menurut definisi (2.4.4), J Hq (pk )
merupakan fungsi Lyapunov di fungsi cost J Hq (pk )
n
. Karena Qi , Ri kesemuanya definit positif, dan
(k ) merupakan bentuk kuadratik dari matriks-matriks Qi dan
Ri , maka menurut Definisi (2.3.7) dan Teorema (2.3.8), J Hq (pk ) positif, yaitu J Hq (pk ) semua k
(k )
(k )
0 . Karena J Hq (pk )
0 untuk semua
n
(k )
, (k ) 0 ,
(k )
definit
q(k )
0 dan
(k ) merupakan fungsi Lyapunov definit positif,
maka menurut Teorema (2.4.5), titik ekuilibrium 0 adalah stabil. Selanjutnya akan dibuktikan bahwa J Hq p
J Hq (pk
1)
(k 1)
J Hq (pk )
(k )
0
92
Karena Karena Qi , Ri kesemuanya definit positif, maka menurut Teorema (2.3.9),
J Hq (pk )
u1q , u2q ,..., uMq ; ( k )
(k ) konveks tegas. Sehingga
merupakan barisan
turun, akibatnya
J Hq (pk
1)
J Hq (pk )
(k 1)
(k )
yang artinya J Hq p
J Hq (pk J Hq (pk )
1)
J Hq (pk )
(k 1) (k )
J Hq (pk )
(k )
(k )
0 untuk
p(k )
0 dan semua k
0 Sehingga, menurut Teorema (2.4.5), titik
ekuilibrium 0 adalah stabil asimtotik.
BAB IV APLIKASI KENDALI MODEL PREDIKTIF TERDISTRIBUSI PADA KANAL IRIGASI Kanal irigasi merupakan systemwide plant yang memuat beberapa reach yang saling berinteraksi dan mencakup area yang luas. Untuk mengoperasikan kanal irigasi sehingga paling aman dan efisien, hal yang paling krusial adalah mengatur ketinggian air sehingga sesuai dengan yang dikehendaki, terutama pada kondisi ekstrim. Perangkat utama kanal irigasi adalah pompa dan pintu air (gate). Teknik kendali MPC terdistribusi akan diaplikasikan pada kanal irigasi yang memuat 5 reach.
4.1 Model Matematika Kanal Irigasi Kanal irigasi yang akan dibahas terdiri dari 5 reach yang diilustrasikan pada Gambar (4.1.1). Diasumsikan tidak ada gangguan pada kanal irigasi seperti gangguan penambahan volume air karena hujan, kebocoran pada aliran air, atau gangguan lainnya.
Gambar 4.1.1 Kanal irigasi 5 reach.
Sistem kanal irigasi di bagi menjadi 5 subsistem, masing-masing subsistem terdiri dari reach dan kontroler lokal yang terletak di upstream gate masing-masing reach. Subsistem ke-4 dan ke-5 memiliki dua kontroler lokal,
93
94
yaitu upstream gate dan pompa pada downstream. Model dari masing-masing subsistem dibentuk berdasarkan volume air di reservoir yang dipengaruhi oleh upstream in-flow (aliran masuk dari hulu) dan downstream out-flow (aliran keluar menuju hilir).
Gambar 4.1.2 Ilustrasi Masing-masing Reach
Misalkan hi (k ) menyatakan tinggi air (m) pada reach i saat k, Ai s menyatakan luas permukaan air (m2) untuk reach i, Qi in dan Qi out secara berturutturut menyatakan in-flow (m3/s) dan out-flow (m3/s) untuk reach i yang diukur pada upstream dan downstream. Didefinisikan qi menyatakan debit aliran air yang melewati gate ke-i, p4 dan p5 berturut-turut menyatakan debit aliran air yang melewati pompa di reach 4 dan 5. Berdasarkan hukum konservasi massa, Qiout
Qiin1 untuk i 1, 2,3 dan Q4out
p4 , Q5out
p5 . Misalkan ketinggian air
(meter) pada reach i saat k adalah hi(k), maka volume air reach i saat k adalah hi ( k ). Ais (m3). Misalkan Ts menyatakan waktu sampling. Perubahan volume air
pada reach i pada saat k ke k +1 dinyatakan oleh persamaan Ais hi (k 1) hi ( k )
Ts Qiin (k ) Qiout (k ) .
(4.1.1)
Sehingga model yang berlaku pada reach i dinyatakan oleh hi (k 1) hi ( k )
Ts Qiin ( k ) Qiout ( k ) , s Ai
(4.1.2)
hi ( k 1)
Ts in Qi (k ) Ais
(4.1.3)
hi ( k )
Ts out Qi (k ). Ais
95
4.1.1 Model Terdesentralisasi Kanal Irigasi Model terdesentralisasi menyatakan model pada masing-masing subsistem secara independent. Satu buah subsistem pada kanal irigasi terdiri dari gate dan
xii (k ) hi (k ) , untuk masing-masing subsistem i, model
reach. Misalkan
terdesentralisasi dinyatakan oleh 1) Subsistem 1
x11 (k 1) y1 (k )
x11 (k )
Ts q1 (k ) A1s
A11 x11 (k ) B11q1 (k ),
x11 (k ) C11x11 (k ),
dengan A11 1, B11
Ts , dan C11 1 . A1s
2) Subsistem 2
x22 (k 1) y2 (k )
x22 (k )
Ts q2 ( k ) A2s
A22 x22 (k ) B22q2 (k ),
x22 (k ) C22 x22 (k ).
dengan A22
1, B22
Ts , dan C22 A2s
1.
3) Subsistem 3
x33 (k 1) y3 (k )
x33 (k )
Ts q3 (k ) A3s
A33 x33 (k ) B33q3 (k ,
x33 (k ) C33 x33 (k ).
dengan A33
1, B33
Ts , dan C33 1 . A3s
4) Subsistem 4 x44 ( k 1)
x44 ( k ) x44 ( k )
y4 ( k )
dengan A44
Ts q4 ( k ) A4s Ts A4s
Ts p4 ( k ) A4s
Ts A4s
q4 ( k ) p4 ( k )
A44 x44 ( k ) B44
x44 ( k ) C44 x44 (k ).
1, B44
Ts A4s
Ts , dan C44 A4s
1.
q4 (k ) p4 ( k )
,
96
5) Subsistem 5 x55 (k 1)
x55 (k ) x55 (k )
y5 ( k )
dengan A55
x55 ( k )
1, B55
Ts q5 (k ) A5s Ts A5s
Ts p5 ( k ) A5s Ts A5s
q5 ( k ) p5 ( k )
A55 x55 ( k )
B55
q5 ( k ) , p5 ( k )
C55 x55 ( k ).
Ts A5s
Ts , dan C55 1 . A5s
4.1.2 Model Interaksi Kanal Irigasi Model interaksi menyatakan efek yang ditimbulkan dari suatu subsistem. Misalkan xij ( k ) menyatakan efek pada subsistem i yang ditimbulkan oleh subsistem j, dalam kasus ini menyatakan efek perubahan ketinggian air pada subsistem i. Berdasarkan persamaan (4.1.3), model interaksi antar subsistem pada kanal irigasi dinyatakan oleh persamaan-persamaan berikut. 1) Interaksi subsistem 1 dengan subsistem 2 adalah pengurangan volume air di reach 1 yang berpindah ke reach 2 melalui gate 2, dinyatakan oleh x12 ( k 1)
dengan B12
Ts q2 ( k ) A1s
B12 q2 ( k ),
Ts . A1s
2) Interaksi subsistem 2 dengan subsistem 3 adalah pengurangan volume air di reach 2 yang berpindah ke reach 3 melalui gate 3, dinyatakan oleh x23 (k 1)
dengan B23
Ts q3 ( k ) A2s
B23q3 ( k ),
Ts . A2s
3) Interaksi subsistem 2 dengan subsistem 5 adalah pengurangan volume air di reach 2 yang berpindah ke reach 3 melalui gate 5, dinyatakan oleh x25 (k 1)
Ts q5 ( k ) A5s
B25 q5 (k ),
97
dengan B25
Ts . A5s
4) Interaksi subsistem 3 dengan subsistem 4 adalah pengurangan volume air di reach 3 yang berpindah ke reach 4 melalui gate 4, dinyatakan oleh x34 ( k 1)
dengan B34
Ts q4 ( k ) A3s
B34 q4 ( k ),
Ts . A3s
Untuk interaksi antara subsistem i dan j yang tidak dituliskan diatas, artinya tidak terjadi interaksi xij
0 .
4.1.3 Model Composite Kanal Irigasi Model composite menggabungkan model desentralisasi dan model interaksi. Model composite untuk sistem kanal irigasi untuk masing-masing subsistem adalah sebagai berikut.
1) Model composite untuk subsistem 1 diberikan oleh
x11 (k 1) x12 (k 1)
1 0 0 0 0
x11 (k )
0 0 0 0 0
x12 (k )
x13 (k 1) x14 (k 1)
0 0 0 0 0
x13 (k )
0 0 0 0 0
x14 (k )
x15 (k 1)
0 0 0 0 0
x15 (k )
A1x1 (k ) B1q1 (k ) W12 q2 (k ), x11 (k ) x12 (k ) y1 (k )
1 1 0 0 0 x13 (k ) x14 (k ) x15 (k )
C1 x1 ( k ),
Ts A1s 0 0
q1 (k )
0 Ts A1s 0
0
0
0
0
q2 (k )
98
dengan x11 (k ) x12 (k ) x1 (k )
x13 (k ) , A1 x14 (k ) x15 (k )
dan C1
Ts A1s
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 , B1 0 0 0 0 0 0 0 0 0 0
0 , W12 0
0
0 0
0 0
1 1 0 0 0.
2) Model composite untuk subsistem 2 diberikan oleh
x21 (k 1) x22 (k 1) x23 (k 1)
0 0 0 0 0 0 1 0 0 0
x24 (k 1)
0 0 0 0 0 0 0 0 0 0
x25 (k 1)
0 0 0 0 0
x21 (k ) x22 (k ) x23 (k ) x24 (k ) x25 (k )
0 0
0 0
Ts q (k ) A2s 3
0 0
0
Ts A2s
0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 Ts A2s 0 0
q2 (k )
0
q5 (k )
x21 (k ) x22 (k ) x23 (k ) x24 (k ) x25 (k )
0 Ts A2s 0 0 0
0 Ts A1s
q2 (k )
,
99
0
0
0
0 Ts q (k ) A2s 3
0 0
0 0
q5 (k )
0
0
p5 (k )
0
Ts A2s
0
0
A2 x2 (k ) B2 q2 (k ) W23q3 (k ) W25
q5 (k ) p5 (k )
x21 (k ) x22 (k ) y2 ( k )
0 1 1 0 1 x23 (k ) x24 (k )
C2 x2 ( k ),
x25 (k )
dengan
x21 (k ) x2 ( k )
W23
x22 (k ) x23 (k ) , A2 x24 (k )
0 0 0 0 0 , B2 0 0 0 0 0
x25 (k )
0 0 0 0 0
0 0
0 0
0 0
Ts , W25 A2s
0 0
0 , 0
0
Ts A2s
0
0 dan C2
0 0 0 0 0 0 1 0 0 0
0 1 1 0 1.
0 Ts A2s 0 0 0
,
,
100
3) Model composite untuk subsistem 3 diberikan oleh x31 (k 1) x32 (k 1)
0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
x33 (k 1) x34 (k 1)
0 0 0 0 0 0 0 0 0 0
x35 (k 1)
0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
x31 (k ) x32 (k ) x33 (k ) x34 (k )
x33 (k ) x34 (k )
Ts q (k ) A3s 3
0 q4 ( k ) Ts A3s 0
0
0
0
0 Ts q (k ) A3s 3
0 0 Ts A3s
0 0
0
0
x35 ( k )
A3 x3 (k ) B3q3 (k ) W34
0 0
0 0
x35 (k ) x31 (k ) x32 (k )
0 0
q4 (k ) p4 (k )
0 0
0
q4 ( k ) p4 ( k )
,
x31 (k ) x32 (k ) y3 (k )
0 0 1 1 0 x33 (k ) x34 (k )
C3 x3 ( k ),
x35 (k ) dengan x31 (k ) x32 (k ) x3 (k )
x33 (k ) , A3 x34 (k ) x35 (k )
dan C3
0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 , B3 0 0 0 0 0 0 0 0 0 0
0 0 1 1 0.
0 0
0 0
0 0
Ts , W34 A3s
0 Ts A3s
0
0
0
0 0
0
,
101
4) Model composite untuk subsistem 4 diberikan oleh x41 (k 1) x42 (k 1) x43 (k 1) x44 (k 1) x45 (k 1)
x41 (k ) x42 (k )
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
x43 (k ) x44 (k )
0 0 0 1 0 0 0 0 0 0
x45 ( k )
x42 (k ) x43 (k ) x44 (k )
0 0 0 1 0 0 0 0 0 0
x45 (k )
q4 (k )
A4 x4 (k ) B4
p4 (k )
0 0
0 q4 ( k ) Ts A4s
0 p4 ( k ) Ts A4s
0
0
0
0
0 0
0 0
Ts A4s
Ts A4s
0
0
x41 (k )
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0
q4 (k ) p4 (k )
,
x41 (k ) x42 (k ) y4 ( k )
0 0 0 1 0 x43 (k ) x44 (k )
C4 x4 ( k ),
x45 (k ) x41 (k ) x42 (k ) dengan x4 (k )
x43 (k ) , A4 x44 (k ) x45 (k )
dan C4
0 0 0 1 0.
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 , B4 0 0 0 1 0 0 0 0 0 0
0 0
0 0
0 Ts A4s
0 , Ts A4s
0
0
102
5) Model composite untuk subsistem 5 diberikan oleh x51 (k 1) x52 (k 1) x53 (k 1) x54 (k 1) x55 (k 1)
x51 (k ) x52 (k )
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
x53 (k ) x54 (k )
0 0 0 0 0 0 0 0 0 1
x55 ( k )
0 0 0 0 0
x51 (k )
0 0 0 0 0 0 0 0 0 0
x52 (k )
0 0 0 0 0
x54 (k )
0 0 0 0 1
x55 (k )
A5 x5 (k )
x53 (k )
0 0
0 0
0 q5 (k ) 0
0 0
Ts A5s
Ts A5s
0
0
0
0
0
0
q5 (k )
0 Ts A5s
0 Ts A5s
p5 (k )
p5 (k )
q5 (k ) , p5 (k )
B5
x51 (k ) x52 (k ) y5 ( k )
0 0 0 0 1 x53 (k ) x54 (k )
C5 x5 ( k ),
x55 (k ) x51 (k ) x52 (k ) dengan x5 (k )
x53 (k ) , A5 x54 (k ) x55 (k )
dan C5
0 0 0 0 1.
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 , B5 0 0 0 0 0 0 0 0 0 1
0 0
0 0
0 0
0 0
Ts A5s
Ts A5s
,
103
Selanjutnya, didefinisikan state dan keluaran untuk subsistem ke-i, yaitu
xi1 (k 1) qi ( k ),
xi 2 (k 1) xi (k 1)
xi 3 (k 1) , ui ( k ) xi 4 (k 1)
i 1,2,3
qi ( k ) pi ( k )
, i
4,5,
xi 5 (k 1) sehingga, state model composite dapat ditulis kembali sebagai x1 (k 1)
A1x1 (k ) B1u1 (k ) W12u2 (k )
x2 (k 1)
A2 x2 (k )
x3 (k 1)
A3 x3 ( k ) B3u3 (k ) W34u 4 (k )
x4 (k 1)
A4 x4 ( k )
B4u4 (k )
x5 ( k 1)
A5 x5 ( k )
B5u5 ( k )
B2u2 ( k ) W23u3 ( k ) W25u5 ( k )
atau dapat dituliskan sebagai x1 ( k 1)
x1 ( k )
x2 (k 1) x3 (k 1) x4 (k 1)
A1 0 0 0
0 A2 0 0
0 0 A3 0
0 0 0 A4
0 0 0 0
x2 ( k ) x3 ( k ) x4 ( k )
x5 (k 1)
0
0
0
0
A5
x5 (k )
0
0 W23 0 0 y1 (k ) y2 (k ) y3 ( k ) y4 (k ) y5 ( k )
0 0 u1 ( k ) 0 0
0
0
B33 u3 (k )
B11
W34 u4 (k ) B44
W25 0 u5 (k ), 0 B55
0
C1
0
0
0
0
x1 (k ) x2 ( k )
0 0
C2 0
0 C3
0 0
0 0
0
0
0
C4
0
x3 (k ) . x4 ( k )
0
0
0
0
C5
x5 (k )
W12 B22 0 u2 ( k ) 0 0
104
Berdasarkan uraian diatas, diperoleh matriks-matriks sistem
A1
A4
1 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 , A2 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 , A3 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 , 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 , A5
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 .
0 0 0 1 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 1
Ts A1s B1
0 , B2 0 0 0
W12
0 Ts A2s 0
, B3
0 0 0 Ts A1s 0 0 0
, W23
0 0
0 0
0 0
0 0
0 0
Ts , B4 A3s
0 Ts A4s
0 , B5 Ts A4s
0 0
0 0
0
0
Ts A5s
Ts A5s
0 0 0 0
0 0
0 0
0 0
0 0
Ts , W25 A2s
0 0
0 , W34 0
0
0 0
Ts A2s
0
0 Ts A3s 0
0
.
.
0
4.2 Deskripsi Systemwide Kanal Irigasi MPC terdistribusi akan diaplikasikan pada kanal irigasi yang diberikan oleh P.J. Overloop (2006). Luas permukaan masing-masing reach diberikan oleh Tabel (4.2.1).
Tabel 4.2.1 Luas permukaan reach kanal irigasi Reach 1 2 3 4 5 Ais (m2) 397 653 503 1530 1614
105
Masing-masing gate memiliki kapasitas maksimum 2,8 m3/s dan masingmasing pump memiliki kapasitas maksimum 2 m3/s. Waktu sampling yang digunakan adalah Ts = 240 detik dan horison prediksi sepanjang 5 step (20 menit). Berdasarkan parameter-parameter tersebut, matriks-matriks model composite untuk sistem ini diberikan oleh
B1
240 397 0 0 0 0
0 240 653 , B3 0 0 0
, B2
0 240 397 , W23 0 0 0
W12
0 0
0 0
240 , B4 503 0 0
0 0
0 240 1530 0
0 0
0 , B5 240 1530 0
0 0
0 0
0 0
240 , W25 653 0 0
0 0 , W34 0 0 240 0 653
0 0
0 0 240 1614
0 0
0 0
0 240 503 0
0
0 , 0 240 1614
.
0 0
4.3 Fungsi Objektif Kanal Irigasi Matriks pembobotan cost masing-masing subsistem adalah 100 untuk pembobotan state, dan 40 untuk pembobotan input. Stage cost masing-masing subsistem didefinisikan oleh Li xi , ui
dengan
Qi
1 T xi (k )Qi xi ( k ) uiT (k ) Ri ui (k ) 2
100, Ri
(4.2.1)
40, i 1,2,3,4,5 . Fungsi objektif subsistem ke-i untuk
horison berhingga Hp = 5 didefinisikan oleh Hp 1 i
Li xˆiT (k
xˆi (k ), uˆi (k ); xi (k )
j | k ), uˆiT (k
j | k)
j 0 4
j
1 T xˆi (k 02
j | k )Qi xˆi (k
j | k ) uˆiT (k
j | k ) Riuˆi (k
j | k) (4.2.2)
106
Setiap gate memiliki kapasitas maksimum 2,8 m3/s dan setiap pump memiliki kapasitas maksimum 2 m3/s, sehingga 0 ui 0 0
2,8 ,i 2
ui
2,8, i 1,2,3 dan
4,5 untuk setiap detik. Karena waktu sampling Ts = 240 s,
maka untuk setiap gate dan setiap pump selama waktu sampling k memiliki kapasitas
0 ui (k j
maksimum
672
j) 672, i 1,2,3 dan
m3/s 0 0
ui (k
dan j)
480 672 ,i 480
m3/s,
sehingga
4,5 untuk setiap
0 . Parameter-parameter tersebut akan menjadi kendala pada optimisasi FC-
MPC dan juga optimisasi pada model negosiasi Nash-Bargaining MPC.
4.4 Aplikasi FC-MPC Pada Kanal Irigasi Algoritma FC-MPC pada Bab III akan diaplikasikan pada persamaan model composite dan fungsi objektif dari kanal irigasi. Untuk melihat hasil aplikasi ini, dilakukan simulasi dengan MATLAB yang akan menghasilkan grafik input dan output dari sistem kanal irigasi dan juga total cost yang dihasilkan. Dari hasil simulasi tersebut, dapat dianalisis performansi dari teknik kendali ini.
4.4.1 Masalah Optimisasi FC-MPC Pada Kanal Irigasi Pada aplikasi ini, state masing-masing subsistem yaitu ketinggian air akan dibawa ke titik 5 meter. Karena fungsi objektif yang digunakan pada aplikasi ini adalah fungsi objektif sebagai regulator (yaitu membawa state ke titik nol), maka titik 5 meter dianggap sebagai titik nol. Sebagai konsekuensinya, nilai state awal akan menjadi negatif. Nilai state awal sesungguhnya adalah 2 untuk subsistem 1, 2, dan 3, sedangkan state awal subsistem 4 dan 5 adalah 0. Sehingga nilai awal yang
akan
digunakan
x1 (k )
x2 (k )
x3 (k )
pada
masalah
2 , dan x4 (k )
x5 (k )
optimisasi
FC-MPC
adalah
5 . Berdasarkan FC-MPC pada
Bab III, masalah optimisasi FC-MPC kanal irigasi untuk setiap subsistem i saat iterasi q dapat ditulis sebagai
107
(i) : FC-MPC Kanal Irigasi i
kendala : 0 ui (k 0
q i i
u (k )
ri ( k )
ij
u
q 1 j
uiq ( k )
(k )
j 1, j i
j | k ) 672, i 1, 2,3
ui (k
0
T
M
1 q T min ui ( k ) q ui ( k ) 2
j | k) 0
672
, i
480 j
4,5
Hp 1
dengan i
diag Qi (0), Qi (1), Qi (2), Qi (3), Qi (4) ,
i
diag Ri (0), Ri (1), Ri (2), Ri (3), Ri (4) ,
i
wi
5 i
wi EiiT
i
5
w j E ji T
Eii
j
E ji ,
wl EliT
ij
l
Elj ,
l 1
j i 5
ri ( k )
wi EiiT
w j E ji T
i f i xi ( k )
j
f j x j (k ),
j i
Eii
Bi Ai Bi Ai4 Bi
0 Bi
0
0 0 Bi
, Eij
Wij AW i ij
0 Wij
0
Ai4Wij
0 0
Ai Ai 2 ,
Ai 3
fi
.
Wij A5
4.4.2 Hasil Simulasi FC-MPC Pada Kanal Irigasi Simulasi dilakukan menggunakan program MATLAB (Lampiran 7). Simulasi ini menggunakan nilai state awal x1 (k )
x4 (k )
x5 (k )
x2 (k ) x3 (k )
2 , dan
5 untuk dengan periode kendali 30 step atau 120 menit. Simulasi
dilakukan dengan dua pendekatan. Pendekatan pertama menggunakan nilai bobot kombinasi konveks dari kelima fungsi objektif adalah sama yaitu wi
1
5
untuk
setiap i. Pendekatan kedua menggunakan bobot kombinasi linier konveks yang proporsional terhadap beban kendali, dalam hal ini adalah volume air yang harus ditambahkan subsistem ke reach sehingga mencapai ketinggian yang diinginkan,
108
yaitu
w1
0,3871, w2
0,2354, w3
0,3056, w4
0,0402, dan w5
0,0317 .
Hasil simulasi berupa grafik input-output masing-masing subsistem diberikan oleh Gambar (4.4.1).
(a)
(b)
(c)
(d)
(e) Gambar 4.4.1 Hasil simulasi FC-MPC dengan nilai wi sama
109
Grafik pada Gambar (4.4.1) adalah grafik input-output masing-masing subsistem untuk kasus 1, yaitu nilai wi sama untuk semua subsistem. Grafik-grafik tersebut adalah hasil simulasi pada program MATLAB yang termuat pada lampiran (7). Grafik dengan tipe stairs menunjukkan nilai input kendali masingmasing subsistem, sedangkan grafik dengan tipe dot menunjukkan nilai output, dalam hal ini adalah ketinggian air yang sudah ditransformasikan dimana titik 5 meter dianggap sebagai titik nol. Hal ini juga berlaku untuk grafik-grafik pada Gambar (4.4.2), Gambar (4.5.1), dan Gambar (4.5.2).
(a)
(b)
(c)
(d)
110
(e) Gambar 4.4.2 Hasil simulasi FC-MPC dengan wi proporsional terhadap beban kendali (volume reach)
Menurut Teorema (3.6.4), sistem lingkar tertutup hasil kendali MPC terdistribusi memiliki titik equilibrium 0 yang stabil asimtotoik. Berdasarkan hasil simulasi, yaitu Gambar (4.4.1) dan (4.4.2), dapat dilihat bahwa state masingmasing subsistem untuk nilai wi sama (kasus 1) dan nilai wi proporsional dengan volume reach (kasus 2) semuanya menuju ke state equilibrium 0 (stabil asimtotik), yang artinya hasil aplikasi ini sesuai dengan hasil analitik Teorema (3.6.4). Selanjutnya akan dilihat waktu yang dibutuhkan oleh subsistem untuk mencapai state equilibrium untuk kasus 1 dan kasus 2. Pada Gambar (4.4.1.d), state mencapai titik nol sekitar menit ke-120, sedangkan untuk Gambar (4.4.2.d), state belum mencapai titik nol pada menit yang sama. Sebaliknya, pada Gambar (4.4.2.c), state mencapai titik nol sekitar menit ke-20, sedangkan pada Gambar (4.4.1.c), state mencapai titik nol sekitar menit ke-120. Hal ini menunjukkan bahwa untuk beberapa subsistem, kasus 1 mencapai state equilibrium lebih cepat dibanding kasus 2. Sebaliknya, untuk beberapa subsistem kasus 2 mencapai state equilibrium lebih cepat dibanding kasus 2. Berdasarkan hasil ini, tidak dapat ditentukan performansi mana yang lebih baik diantara kedua kasus ini jika dilihat dari waktu yang dibutuhkan untuk mencapai state ekuilibrium untuk masingmasing subsistem. Pada akhir bab ini, akan diselidiki nilai total cost yang dihasilkan.
111
4.5 Aplikasi Nash-Bargaining MPC Pada Kanal Irigasi Model negosiasi Nash-Bargaining MPC akan diaplikasikan pada kanal irigasi. Hasil aplikasinya disimulasikan dalam bentuk grafik input-output. Nilainilai awal yang digunakan adalah sama seperti pada aplikasi FC-MPC. Pertama, pembobotan kombinasi linier konveks pada fungsi objektif subsistem adalah sama, yaitu 1
5
untuk semua subsistem (kasus 3). Kedua, pembobotan untuk
setiap subsistem proporsional terhadap beban kendali (volume reach) yaitu
w1 0,3871, w2 0,2354, w3 0,3056, w4 dan w5 0,09 (kasus 4).
4.5.1
0,0402, dan w5
0,0317
Model Negosiasi Nash-Bargaining MPC Pada Kanal Irigasi Model negosiasi Nash-Bargaining MPC yang diberikan pada pembahasan
subbab (3.5.2) akan diaplikasikan pada kanal irigasi. Fungsi objektif subsistem yang digunakan adalah fungsi objektif (4.2.2).
4.5.2 Hasil Simulasi NB- MPC Pada Kanal Irigasi Simulasi dilakukan menggunakan program MATLAB (Lampiran 8). Simulasi ini menggunakan nilai state awal x1 (k )
x4 (k )
x5 (k )
x2 (k )
x3 (k )
2 , dan
5 untuk dengan periode kendali 30 step atau 120 menit. Hasil
simulasi berupa grafik input-output masing-masing subsistem diberikan oleh Gambar (4.5.1).
(a)
(b)
112
(c)
(d)
(e) Gambar 4.5.1 Hasil simulasi model negosiasi Nash-Bargaining MPC dengan proporsi wi sama Grafik-grafik pada Gambar (4.5.1) adalah grafik input-output masingmasing subsistem untuk kasus 3, yaitu nilai wi sama untuk semua subsistem. Grafik-grafik tersebut adalah hasil simulasi pada program MATLAB yang termuat pada Lampiran (8). Grafik dengan tipe stairs menunjukkan nilai input kendali masing-masing subsistem, sedangkan grafik dengan tipe dot menunjukkan nilai output, dalam hal ini adalah ketinggian air yang sudah ditransformasikan dimana titik 5 meter dianggap sebagai titik nol. Hal ini juga berlaku untuk grafik-grafik pada Gambar (4.5.2) tetapi menggunakan nilai wi yang proporsional dengan volume reach.
113
(a)
(b)
(c)
(d)
(e) Gambar 4.5.2 Hasil simulasi Nash-Bargaining MPC dengan wi proporsional terhadap volume reach
114
Berdasarkan Gambar (4.5.1) dan (4.5.2), dapat dilihat bahwa state semua subsistem adalah stabil asimsotik, yaitu konvergen ke titik equilibrium 0, hasil ini menunjukkan hasil yang sesuai dengan hasil analitik Teorema (3.6.4). Sama halnya pada hasil simulasi aplikasi FC-MPC pada kanal irigasi, untuk kedua kasus diatas, terdapat beberapa subsistem untuk kasus 3 yang mencapai titik equilibrium 0 yang lebih cepat dibanding kasus 4, sebaliknya, ada pula yang lebih lambat. Sehingga performansi kedua kasus ini tidak dapat dibandingkan dengan cara menyelidiki waktu yang dibutuhkan untuk mencapai titik ekuilibrium dari masing-masing subsistem. Pada subbab selanjutnya, diselidiki nilai total cost yang dihasilkan dari kasus 1 sampai kasus 4. Berdasarkan grafik input-output yang dihasilkan dari simulasi diatas, dapat diamati bahwa teknik kendali MPC terdistribusi bekerja menghitung nilai kendali optimal yang diaplikasikan ke sistem sedemikian sehingga state sistem menuju titik nol. Sebagai ilustrasi bagaimana teknik kendali ini bekerja, khusus untuk subsistem ke-4 (Gambar 4.5.2.d) dan subsistem ke-5 (Gambar 4.5.2.e) yang masing-masing memiliki dua input, yaitu aliran masuk melalui gate (input 1) dan aliran keluar melalui pump (input 2), dapat di lihat bahwa apabila state masih kurang dari nol, maka input 1 bekerja menambahkan air, sedangkan input 2 sama dengan nol. Begitu pula sebaliknya, apabila state melebihi titik nol, maka input 1 sama dengan nol dan input 2 bekerja mengeluarkan air, walaupun pada gambar-gambar diatas tidak terlihat karena nilainya sangat kecil dibandingkan skala grafik yang digunakan. Hal inilah yang menunjukkan peran kendali pada kasus ini, bagaimana teknik kendali MPC terdistribusi bekerja mengendalikan state ketinggian air dari kanal irigasi 5 subsistem yang saling berinteraksi.
4.6 Perbandingan Performansi FC-MPC dan NB-MPC Pada Kanal Irigasi Pada subbab ini, akan dibandingkan nilai cost yang dihasilkan oleh FCMPC dan NB-MPC untuk kasus wi sama dan kasus wi proporsional. Perbandingan cost untuk masing-masing subsistem dan total cost untuk masing-masing kasus diberikan oleh Tabel (4.6.1) berikut.
115
Tabel 4.6.1 Perbandingan cost FC-MPC dan NB-MPC kanal irigasi Sub Sistem (i)
FC-MPC (wi sama)
1 2 3 4 5 Total
26,3812 55,7845 37,1724 2.212,1000 1.036,0000 3.367,4000
FC-MPC (wi proporsional terhadap volume reach) 51,0607 71,0134 56,5085 1077,4000 297,0759 1553,0 000
NB-MPC (wi sama)
46,1974 146,7652 73,4052 918,9722 974,7592 2.160,100
NB-MPC (wi proporsional terhadap volume reach) 63,0056 202,3558 638,4223 200,4227 172,8764 1.277,1000
Jika ditinjau dari cost yang dihasilkan pada Tabel (4.6.1), khusus untuk aplikasi ini, NB-MPC memberikan total cost yang lebih kecil dibandingkan FCMPC. Berdasarkan tabel diatas, dapat dilihat bahwa ada beberapa cost subsistem dimana FC-MPC lebih kecil dibanding NB-MPC, tetapi ada yang sebaliknya. Hal ini menunjukkan bahwa kita tidak dapat membandingkan cost masing-masing subsistem, tetapi hanya dapat dibandingkan total cost dari semua subsistem.
Hasil simulasi diatas merupakan hasil simulasi aplikasi MPC terdistribusi pada kanal irigasi 5 subsistem dengan konfigurasi seperti pada Gambar (4.1.1). Untuk konfigurasi yang berbeda ataupun banyaknya subsistem yang berbeda, perbandingan cost yang dihasilkan belum tentu sama dengan hasil diatas, sehingga perlu dilakukan simulasi khusus untuk sistem yang dikendalikan.
BAB V PENUTUP
5.1 Kesimpulan Systemwide control merupakan masalah kendali yang memuat beberapa subsistem yang saling berinteraksi. Persamaan state space dari systemwide dapat dinyatakan sebagai model composite yang merupakan gabungan dari model terdesentralisasi dan model interaksi antar subsistem. Masalah kendali ini memuat M subsistem dan M fungsi objektif. Dengan memandang subsistem beserta fungsi objektifnya sebagai pemain, maka masalah ini dapat dipandang sebagai masalah permainan dinamis kooperatif, yang dapat diselesaikan dengan menentukan solusi Pareto, kemudian dilanjutkan dengan menentukan solusi Nash-Bargaining. Teknik kendali MPC diaplikasikan pada masing-masing subsistem, kemudian diselesaikan secara terdistribusi, yaitu semua subsistem bekerjasama untuk menentukan kendali optimal, dalam hal ini menentukan solusi Pareto dan Nash-Bargaining, dengan melakukan pertukaran informasi nilai kendali pada setiap iterasi yang dilakukan. Pada pembahasan ini, algoritma FC-MPC yang akan menghasilkan solusi Pareto, kemudian dilanjutkan dengan model negosiasi NashBargaining MPC yang akan menghasilkan solusi Nash-Bargaining. Berdasarkan Teorema (3.6.4), algoritma MPC terdistribusi akan menghasilkan state equilibrium 0 yang stabil asimsotik. Teknik kendali MPC terdistribusi tersebut diaplikasikan pada kanal irigasi 5 subsistem dengan konfigurasi pada Gambar (4.1.1). Berdasarkan simulasi menggunakan program MATLAB, output untuk semua subsistem konvergen ke state equilibrium 0, hal ini sesuai dengan hasil analitik yang diberikan pada Teorema (3.6.4). Jika dilihat dari total cost yang dihasilkan, model negosiasi Nash-Bargaining MPC menghasilkan total cost yang lebih kecil dari pada FCMPC baik menggunakan wi sama ataupun wi proporsional terhadap beban kendali (volume reach). Sedangkan total cost terkecil diberikan oleh Nash-Bargaining MPC dengan wi proporsional terhadap volume reach. Akan tetapi, hasil ini adalah 116
117
hasil simulasi numerik untuk sistem kanal irigasi khusus untuk 5 subsistem dengan konfigurasi Gambar (4.1.1), sehingga perlu simulasi lanjut untuk menyelidiki hasil untuk banyak subsistem berbeda ataupun konfigurasi berbeda, ataupun untuk systemwide yang berbeda.
5.2 Saran Pada tesis ini, model systemwide yang digunakan pada teknik kendali MPC terdistribusi diasumsikan tidak ada gangguan. Model gangguan dapat ditambahkan pada model ini untuk membentuk teknik kendali MPC terdistribusi yang memuat gangguan. Teknik kendali MPC terdistribusi hanya diaplikasikan pada kanal irigasi yang terdiri dari 5 reach dan konfigurasi Gambar (4.1.1). Pada penelitian lebih lanjut, dapat dikembangkan aplikasi pada kanal irigasi untuk reach yang lebih banyak dan dengan konfigurasi yang berbeda. Total cost yang dihasilkan merupakan hasil simulasi khusus untuk konfigurasi Gambar (4.1.1), sehingga masalah ini masih terbuka untuk menyelidiki secara teoritis bagaimana pengaruh pembobotan kombinasi linier konveks dari fungsi objektif. Selain itu, juga dapat diselidiki secara teoritis apakah total cost yang dihasilkan oleh solusi Nash-Bargaining selalu lebih kecil dari solusi Pareto atau tidak. Pada tesis ini digunakan kanal irigasi yang memuat cabang, tetapi efek dari percabangan ini diabaikan, sehingga dapat diteliti lebih lanjut pengaruh dari percabangan terhadap performansi sistem kanal irigasi. Teknik kendali MPC terdistribusi juga dapat diaplikasikan pada systemwide lain, seperti pada sistem reaktor kimia yang memuat beberapa tangki, sistem lalu-lintas yang memuat beberapa persimpangan, atau sistem-sistem lain yang memuat beberapa subsistem yang saling berinteraksi.
118
DAFTAR PUSTAKA
Anton, H., 1997, Aljabar Linier Elementer, Alih bahasa oleh Pantur Silaban dan I Nyoman Susila, Erlangga, Jakarta. Bazaraa, M.S., Sherali, H.D., Shetty, C.M., 2006, Nonlinear Programming, Theory and Algorithms, 3rd edition, John Wiley and Sons, New Jersey. Boyd, S., Vandenberghe, L., 2004, Convex Optimization, Cambridge University Press. Camacho, E.F., Bordons, C., 1999, Model Predictive Control 2nd edition, Springer-verlag, London. Doan, M.D., T. Keviczky, and B. De Schutter, “An improved distributed version of Han’s method for distributed MPC of canal systems,” Proceedings of the 12th IFAC Symposium on Large Scale Systems: Theory and Applications, Villeneuve d’Ascq, France, 6 pp., July 2010. Engwerda, Jacob, 2005, LQ Dynamic Optimization and Differential Games, John wiley & sons, Tilburg University, Netherlands. Elaydi, S., 2005, An Introduction to Difference Equations, 3rd ed., Springer Science and Business Media, inc., USA. Fletcher, R., 2000, Practical Methods of Optimization, A Wiley-Interscience Publication, UK. Larson, R., Falvo, D.C., 2009, Elementary Linear Algebra, 6th ed., Houghton Mifflin Harcourt Publishing Company, USA. Maciejowski, J.M., 2001, Predictive Control with Constraints, Prentice Hall, USA. Myerson, R.B., 1997, Game Theory Analysis of Conflict, Harvard University Press, England. Narahari, Y., 2009, Game Theory (Lecture Note Chapter 2), Department of Computer Science and Automation, Bangalore, India. Ogata, Katsuhiko, 1995, Discrete Time Control System, Prentice-Hall, USA. Overloop, V., 2006, Model Predictive Control on Open Water Systems, Ph.D Thesis, Delft University Press, The Netherlands. Trench, W. F., 2003, Introduction to Real Analysis, Pearson Education, USA.
119
Valencia, F., Espinosa, J.J., Schutter, B.D., Stankova, K., Feasible-Cooperation Distributed Model Predictive Control Scheme Based on Game Theory, 18th IFAC World Congress, Milano (Italy) August 28 - September 2, 2011. Venkat, Aswin N., 2006, Distributed Model Predictive Control: Theory and Applications, UNIVERSITY OF WISCONSIN–MADISON, UK. Wang, L., 2009, Model Predictive Control System Design and Implementation Using MATLAB, Springer, Australia. Nocedal, J., Wright, S.J., 1999, Numerical Optimization, Springer-Verlag, New York, inc.
LAMPIRAN
121
Lampiran 1 : Program MATLAB Contoh (2.5.1) A=0.5; B=[1.5 -1.5]; C=1; Qi=10; Ri=100; Hp=5; Hc=12; r=10; xpast=0; upast=[0;0]; QQ=blkdiag(Qi,Qi,Qi,Qi,Qi); RR=blkdiag(Ri,Ri,Ri,Ri,Ri,Ri,Ri,Ri,Ri,Ri); z12=zeros(1,2); Tk=[r;r;r;r;r]; Psi=blkdiag(C,C,C,C,C)*[A;A^2;A^3;A^4;A^5]; Upsilon=blkdiag(C,C,C,C,C)*[B;B+A*B;B+A*B+A^2*B;B+A*B+A^2*B+A^3*B; B+A*B+A^2*B+A^3*B+A^4*B]; Theta=blkdiag(C,C,C,C,C)*[B,z12,z12,z12,z12;B+A*B,B,z12,z12,z12; B+A*B+A^2*B,B+A*B,B,z12,z12; B+A*B+A^2*B+A^3*B,B+A*B+A^2*B,B+A*B,B,z12; B+A*B+A^2*B+A^3*B+A^4*B,B+A*B+A^2*B+A^3*B, B+A*B+A^2*B,B+A*B,B]; for p=1:Hc Ek = Tk - Psi*xpast - Upsilon*upast; G = 2.*Theta'*QQ*Ek; H = Theta'*QQ*Theta + RR; DeltaUkopt = (1/2).*inv(H)*G; Deltau1kopt=DeltaUkopt(1); Deltau2kopt=DeltaUkopt(2); u1(p)=upast(1)+ Deltau1kopt; u2(p)=upast(2)+ Deltau2kopt; xk(p)=A*xpast + B*upast; yk(p)=xk(p); xpast=xk(p); upast=[u1(p);u2(p)]; end k=1:Hc; figure subplot(3,1,1) stairs(k,u1) xlabel('Waktu Instant') ylabel('Input 1') subplot(3,1,2) stairs(k,u2) xlabel('Waktu Instant') ylabel('Input 2') subplot(3,1,3) plot(k,yk,'ob') hold on plot(k,r,'.g') xlabel('Waktu Instant') ylabel('Output')
122
Lampiran 2. Program MATLAB Contoh (2.6.1) A=1; x0=4; J1=@(u)(x0^2+u(1)^2+(A*x0+u(1)+u(2))+u(3)^2); J2=@(u)(x0^2+u(2)^2+(A*x0+u(1)+u(2))+u(4)^2); p=1; for a=0.1:0.1:0.9 u0=[10;10;10;10]; J=@(u)(a.*(x0^2+u(1)^2+(A*x0+u(1)+u(2))+u(3)^2)+(1a).*(x0^2+u(2)^2+(A*x0+u(1)+u(2))+u(4)^2)); [uopt,Jopt] = fminsearch(J,u0); J1opt(p)=J1(uopt); J2opt(p)=J2(uopt); p=p+1; end figure plot(J1opt,J2opt) xlabel('Nilai Fungsi Objektif J1') ylabel('Nilai Fungsi Objektif J2')
123
Lampiran 3. Program MATLAB Contoh (2.6.2) A=1; x0=4; J1=@(u)(x0^2+u(1)^2+(A*x0+u(1)+u(2))+u(3)^2); J2=@(u)(x0^2+u(2)^2+(A*x0+u(1)+u(2))+u(4)^2); J3=@(u)(u(5)); u0=[2,2,2,2,2]; option=optimset('Algorithm','active-set'); dis1=fmincon(J3,u0,[],[],[],[],[],[],@mynonlincon,option) dis2=fmincon(J3,u0,[],[],[],[],[],[],@mynonlincon2,option) d1=J1(dis1) d2=J2(dis2) function [c,ceq]=mynonlincon(u) A=1; x0=4; J1=@(u)(x0^2+u(1)^2+(A*x0+u(1)+u(2))+u(3)^2); J2=@(u)(x0^2+u(2)^2+(A*x0+u(1)+u(2))+u(4)^2); c=u(5)-(x0^2+u(1)^2+(A*x0+u(1)+u(2))+u(3)^2); ceq=[]; end function [c,ceq]=mynonlincon2(u) A=1; x0=4; J1=@(u)(x0^2+u(1)^2+(A*x0+u(1)+u(2))+u(3)^2); J2=@(u)(x0^2+u(2)^2+(A*x0+u(1)+u(2))+u(4)^2); c=u(5)-(x0^2+u(2)^2+(A*x0+u(1)+u(2))+u(4)^2); ceq=[]; end A=1; x0=4; J1=@(u)(x0^2+u(1)^2+(A*x0+u(1)+u(2))+u(3)^2); J2=@(u)(x0^2+u(2)^2+(A*x0+u(1)+u(2))+u(4)^2); u0=[10,10,10,10]; d=[32;32]; Jnash=@(u)(-1.*((3/4).*log(d(1)(x0^2+u(1)^2+(A*x0+u(1)+u(2))+u(3)^2))+(1/4).*log(d(2)(x0^2+u(2)^2+(A*x0+u(1)+u(2))+u(4)^2)))); option=optimset('Algorithm','active-set'); solNash=fmincon(Jnash,u0,[],[],[],[],[],[],@nonlincon,option) J1nash=J1(solNash) J2nash=J2(solNash)
124
Lampiran 4. Flow chart algoritma active-set
MULAI
Diberikan
dan
0
active-set
0
,r=0
( r 1)
r
r
untuk min r 1
q r
(r )
(r)
r r 1
Selesaikan T
(r )
r 1
r
T
H
H 0
0
0
a
0
0
0
a
Tidak Ya
r
0 ,
r
Tidak r
{q}
0
Ya Tidak q
0 Ya
* opt
r r
SELESAI
r
1
r 1
a r
min 1, a
r
,
r
{a}
a
r
0
125
Lampiran 5. Flow chart algoritma FC-MPC MULAI
ui0 , xi ( k ), 0, i
i
qmax ( k ) q
u1q
w1u1q* (1 w1 )u1q
1
u1q
u1q
1,
u2q*
u1q* arg( 1 )
1
u2q
0,
i
M
,
0,
0, ,
i
arg(
2
1
uMq*
)
w2u2q* (1 w2 )u2q
1
u2q u2q
2
q
uMq
1
uMq
M
untuk
beberapa i
M
dan
qmax (k )
Tidak
uiopt : uiq ,
i
SELESAI
M
M
)
wM uMq* (1 wM )uMq 1
q 1
i
q
1
arg(
Ya
uMq
1
126
Lampiran 6. Flow chart model negosiasi Nash-Bargaining MULAI
ui0 , xi ( k ), i
q
u1d (k )
arg min max u1 ( k ) u 1 ( k )
d1 (k )
1
1
1, qmax ( k )
, 0,
0
uM d ( k )
arg min max
d M (k )
u d (k )
(ui (k ), u i )
d1 ( x), d2 (k ),..., d3 (k )
uM ( k ) u
M
M
(k )
M
u (k )
u d (k )
T
r
d r (k )
1
M
u (k )
d (k )
r
0, i
0,
i
(ui (k ), u i )
dr (k )
2
3
4
5
127
1
u1q*
u1q
arg( 1 )
3
2
uiq (k )
w1u1q* (1 w1 )u1q
uMq*
wiuid (k ) (1 wi )uiq 1 (k )
uMq
1
1
u1q u1q
arg(
M
i
M
uiq (k )
)
uMq
uMq 1
q 1
untuk beberapa
i
M
wM uMq* (1 wM )uMq 1
1
q
4
dan q
qmax (k ) Tidak
uiopt : uiq ,
i
SELESAI
M
Ya
wiuid (k ) (1 wi )uiq 1 (k )
5
128
Lampiran 7. Program MATLAB FC-MPC kanal irigasi
129
130
131
132
133
134
135
136
137
Lampiran 8. Program MATLAB Nash-Bargaining MPC kanal irigasi
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
