Document Type : Research Article
Authors
1 Department of Electrical Engineering, Yazd Branch, Islamic Azad University, Yazd, Iran
2 Department of industrial Engineering, Meybod Univercity, Meybod, Iran
Abstract
Keywords
Main Subjects
ریزشبکه نوع جدیدی از شبکه برق مبتنی بر DG (تولید توزیع‌شده) است که می‌تواند به‌صورت متصل به شبکه یا خارج از شبکه راه‌اندازی شود. همان‌طور که منابع جدید انرژی ازجمله فتوولتائیک، توربین‌های بادی و فناوری‌های ذخیره انرژی با گذشت زمان به بلوغ می‌رسند و هزینه‌های آنها به‌سرعت درحال کاهش است، ریزشبکه‌ها نیز به‌تدریج شروع به جلب‌توجه بیشتر محققان می‌کنند. برای یک شبکه بزرگ، ریزشبکه را می‌توان به‌عنوان یک واحد کامل قابل کنترل در نظر گرفت. هر منبع میکروگرید می‌تواند هنگامی که به شبکه متصل است، برق را به شبکه منتقل کند. هنگامی که شبکه از کار می‌افتد یا نیاز به تعمیر و نگهداری دارد، ریزشبکه می‌تواند مستقل عمل کند و برای تأمین برق بار خود خودکفا باشد.
ریزشبکه‌ها معمولاً ظرفیت تولید محدودی دارند؛ به‌طوری‌که نمی‌توانند کل بار شبکه را تأمین کنند و پایداری مبدل‌های داخلی آنها برای پایداری کل ریزشبکه مهم است. هنگامی که به شبکه متصل است، ولتاژ و فرکانس سیستم به‌راحتی تحت‌تأثیر حضور شبکه بزرگ‌تر قرار نمی‌گیرد؛ اما زمانی که ریزشبکه به‌طور مستقل کار می‌کند، ولتاژ و فرکانس آن به‌راحتی تحت‌تأثیر بار سایر عوامل قرار می‌گیرد؛ بنابراین، مطالعه عمیق پایداری و مقاومت ریزشبکه در مقابله با اختلالات و بازگشت به حالت پایدار ضروری است [2،1 [.
در مرجع [3]، یک استراتژی کنترلی جدید مبتنی بر مفهوم کنترل فیدبک حالت برای کنترل قوی مبدل‌های توان دو طرفه در ریزشبکه‌های هیبریدی ارائه شده است. در این روش نسبت به روش کنترل افتی نتایج پذیرفتنی و مؤثری به دست آمده است و تنظیم ولتاژ و اشتراک توان به نحو مناسب‌تری انجام می‌شود.
در مرجع [4]، ابزارها و تکنیک‌های مختلف تجزیه و تحلیل پایداری برای ریزشبکه‌ها بحث می‌شود. در این مرجع تعاریفی از مسائل پایداری ریزشبکه ارائه شده و ناپایداری در ریزشبکه‌ها براساس ساختار آنها طبقه‌بندی شده است. در مرجع [5]، یک روش کنترل پیشنهاد شده است که DGها را قادر می‌کند تغییرات آتی در تقاضای بار را سریع‌تر و پیشگیرانه‌تر جبران کنند. در مرجع [6]، روشی مطرح شده است تا عوامل گذرا را مشخص و تأثیرگذاری آنها را بر سیستم محدود کند. در مرجع [7]، برخی از الگوریتم‌های پیشرفته برای انتقال میکروگرید به حالت جزیره‌ای و سنکرون‌سازی هنگام برگشت میکروگرید از حالت جزیره به شبکه معرفی می‌شود. در مرجع [8]، طرح یک اینورتر متصل به شبکه یکپارچه (SGI)[i] برای میکروگرید پیشنهاد می‌شود. SGI پیشنهادی قادر است در حالت متصل به شبکه یا حالت مستقل با کم‌کردن اثرات اتصال و جداشدن از شبکه کار کند.
در این مقاله یک روش کنترلی با استفاده از کنترلکنندة فیدبک حالت[ii] برای کنترل کلیدزنی مبدل‌های واسط و جبران بارهای نامتعادل و غیرخطی ارائه می‌شود. این کنترل‌کننده از یک تنظیم‌کنندة درجه 2 خطی (LQR)[iii] برای تولید یک حاشیه بهره مناسب استفاده کرده است و به تخصیص بهینه بار بین منابع تولید پراکنده خواهد پرداخت. درنهایت، برای بهبود عملکرد گذرای میکروگرید‌های شامل یک و دو منبع تولید پراکنده، کنترلکنندة PID را با استفاده از الگوریتم تکاملی بهینه‌سازی اجتماع ذرات[iv] طراحی خواهیم کرد. در مرحله شبیه‌سازی شرایط مختلفی برای بارگذاری میکروگرید در نظر گرفته می‌شود. بارهای نامتعادل و غیرخطی که تعادل شبکه را مختل و کیفیت توان را دچار مشکل می‌کنند، توسط کنترلکنندة پیشنهادی جبران می‌شوند و منابع تولید پراکنده توسط یک منبع ولتاژ DC ایدئال مدل‌سازی می‌شوند که به‌واسطة مبدل‌های منبع ولتاژ به میکروگرید متصل می‌شوند و در هر دو حالت عملکرد میکروگرید، یعنی حالت‌های متصل به شبکه و جزیره‌ای بررسی خواهد شد. مراحل کار در این مقاله شامل موارد زیر است:1- طراحی کنترلکنندة فیدبک حالت با استفاده از رگولاتور درجه دوم خطی و استفاده از قائده بریسون در انتخاب اولیه ماتریس‌ها؛ 2- مطالعه میکروگرید شامل یک DG در دو حالت متصل به شبکه و جزیره‌ای، بررسی تأثیر انتخاب پارامترهای شاخص LQR بر پاسخ گذرای شبکه و به دست آوردن پارامترهای شاخص LQR با استفاده از الگوریتم اجتماع ذرات؛3- طراحی کنترلکنندةPID در حلقه فیدبک با استفاده از الگوریتم اجتماع ذرات؛ 1- طراحی کنترلکنندة فیدبک حالت با استفاده از رگولاتور درجه دوم خطی[v].در طراحی کنترلکنندة فیدبک حالت به کمک رگولاتور درجه دوم خطی شکل (1)، فرض می‌کنیم حالت x کاملاً قابل اندازه‌گیری و برای کنترل در دسترس باشد.
|
|
|
+ |
|
xref |
|
کنترل کننده |
|
|
شکل (1): کنترلکنندة فیدبک حالت با استفاده از رگولاتور درجه دوم خطی در پیکربندی فوق مدل فضای حالت فرایند به‌صورت زیر است]6[:
| (1) |
در بلوک دیاگرام شکل (1)، z(t) خروجی کنترل‌شده[vi] نامیده می‌شود و مربوط به سیگنالی است که می‌خواهیم در کمترین زمان ممکن به حداقل مقدار ممکن خود برسانیم. همچنین، در بلوک فوق m تعداد ورودی‌ها، n تعداد حالت‌های سیستم و l تعداد خروجی‌های کنترلی است.مسئله LQR بهینه شامل یافتن ماتریس انتقال کنترلکنندة C(s) به نحوی است که معیار زیر را حداقل کند [7]:
| (2) |
در رابطه فوق، ماتریس وزن‌دهی[vii]Q یک ماتریس معین یا نیمه‌معین مثبت، حقیقی و متقارن و ماتریس جریمه کنترل[viii] R نیز ماتریس مثبت معین، حقیقی و متقارن است و ρ یک عدد ثابت مثبت است. 1-1- قاعده بریسون[ix]انتخاب اولیه برای ماتریس‌های R و Q توسط قانون بریسون به شکل انتخاب ماتریس‌های قطری به‌صورت زیر مطرح می‌شود [8]: i= {1, 2,…, l}
| (3) | ={1, 2,…, m} j |
که متناظر با معیار زیر است [8]:
| (4) |
در رابطه فوق، l تعداد خروجی‌های کنترلی و در اینجا تعداد حالت‌ها است و m تعداد ورودی‌ها را نشان می‌دهد.قاعده بریسون اساساً متغیرهایی را تنظیم می‌کند که در JLQR ظاهر می‌شوند؛ به‌طوری‌که بیشترین مقدار پذیرفتنی برای هر عبارت یکی باشد. اگرچه این قاعده برخی اوقات نتیجه خوبی به دست می‌دهد،بیشتر به‌صورت یک نقطه شروع برای روش طراحی آزمون و خطا با هدف به دست آوردن مشخصه مطلوب برای سیستم حلقه بسته استفاده می‌شود.اساس عملکرد کنترلکنندة فیدبک حالت LQR برای معیار رابطه (4)، به دست آوردن یک ماتریس بهره ساده به‌صورت زیر است [12،11،10،9،8]:
| (5) | u= -Kx |
K یک ماتریس m×n به شکل زیر است:
| (6) |
P یک ماتریس مثبت معین[x] و یک جواب یکتا برای معادله زیر است [9]:
| (7) |
معادله فوق به معادله جبری ریکاتی[xi] معروف است. شکل کلی شاخص LQR بهینه به‌صورت زیر است [6،5،4،3]:
| (8) |
| (9) |
در روابط فوق، ابتدا مقادیر اولیه R و Q توسط قاعده بریسون تعیین می‌شوند و سپس با استفاده از روابط فوق ماتریس‌های و محاسبه می‌شوند [5]:
| (10) |
در این حالت K به‌صورت زیر محاسبه می‌شود [9]:
|
(11) |
پاسخ معادله جبری ریکاتی به‌صورت زیر است [9]:
| (12) |
در محیط نرم‌افزار MATLAB، دستور [K, P, e]=LQR (A, B, Q, R, N) مقدار بهینه کنترلکنندة فیدبک حالت LQR را برای فرایند زیر محاسبه می‌کند [9]:
| (13) |
در رابطه بالا، K ماتریس بهره، P پاسخ معادله ریکاتی و e مقادیر ویژه حلقه بسته از رابطه زیر به دست می‌آید [9]:
| (14) | e= eig (A-B*K) |
برای انتخاب ماتریس‌های Q و R گزینه‌های دیگری نیز وجود دارد؛ اما باید به این نکته توجه کرد که همواره این ماتریس‌ها باید مقادیر مثبت معین داشته باشند. 2-1- کنترلکنندة فیدبک حالتدر کنترلکنندة فیدبک حالت، حالت‌های انتخابی سیستم با مقادیر مرجع خود مقایسه می‌شوند تا نحوه کلیدزنی مبدل را مشخص کنند. تولید مقادیر مرجع برای ولتاژ خروجی vcf و جریان ic از روی شرایط پخش بار آسان است؛ اما تولید مقدار مرجع برای جریان i1 امری دشوار است. برای سهولت این امر، بردارهای حالت جدید را به‌صورت زیر تعریف می‌کنیم [16]:
| (15) |
بنابراین، ماتریس تبدیل حالت زیر را خواهیم داشت]11 [:
| (16) | = x= |
معادلات حالت تبدیل‌یافته با ترکیب معادلات 11 و 12 به‌صورت زیر به دست می‌آیند ]11[:
| (17) |
قانون کنترلی به‌صورت زیر تعریف می‌شود[9]:
| (18) |
در رابطه بالا K یک ماتریس بهره[xii] و xref بردار مرجع است. ماتریس بهره را می‌توان با استفاده از رگولاتور درجه دوم خطی گسسته[xiii] (LQR) با ماتریس وزن‌دهی[xiv] حالت Q و جریمه کنترلی R به دست آورد.یک رگولاتور درجه دو خطی برای تولید یک حاشیه بهره نامحدود و حاشیه فاز حداقل 60 درجه به کار گرفته می‌شود. یک ویژگی مهم دیگر LQR این است که در مقابل ورودی غیرخطی بسیار انعطاف‌پذیر است. این موضوع را از شکل (2) می‌توان فهمید. وقتی خطا بزرگ باشد و کنترل بین 1و -1 محدود شده باشد، درایه‌های ماتریس بهره K باید کوچک باشند. برای مجموعه‌هایی از مقادیر R که درحال کاهش‌اند، مجموعه‌ای متناظر از مقادیر درحال افزایش K به دست خواهیم آورد؛ بنابراین، همیشه یک مقدار برای R وجود دارد که Kx به‌طور مناسبی محدود شود.
|
u |
|
x |
|
z |
|
+ |
|
_
|
|
xref |
|
K
|
|
|
شکل (2): استفاده از رگولاتور درجه دوم خطی در فیدبک حالت شکل (2) بلوک دیاگرام یک کنترلکنندة فیدبک حالت با حضور رگولاتور درجه دوم خطی را نشان می‌دهد. در این بلوک دیاگرام پارامتر x نمایانگر حالت‌های سیستم است و پارامتر z خروجی کنترلی سیستم نامیده می‌شود. مقدار z برابر است با [10،9،8]:
| (19) |
بلوک دیاگرام شکل (2) می‌تواند دارای دو نوع خروجی مجزا باشد:خروجی اندازه‌گیری‌شده[xv] y(t): مربوط به سیگنالی است که قابل اندازه‌گیری است و بنابراین برای اهداف کنترلی در دسترس است. اگر C(s) ماتریس انتقال کنترل‌کننده باشد، خواهیم داشت [7]:
| (20) |
که U(s) و Y(s) به‌ترتیب تبدیل لاپلاس ورودی u(t) و خروجی اندازه‌گیری y(t) را نشان می‌دهند.خروجی کنترل‌شدة z(t): مربوط به سیگنالی است که می‌خواهیم در کمترین زمان ممکن به حداقل مقدار ممکن خود برسانیم.در بعضی موارد داریم: z(t) = y(t)؛ به این معنی که هدف از کنترل ما این است که کل خروجی اندازه‌گیری‌شده بسیار کوچک باشد. برای زمانی که خروجی y(t) یک بردار باشد، اغلب یکی از نیازهای مسئله کوچک‌کردن یکی از خروجی‌های اندازه‌گیری‌شده، مثلاً y1(t) است. در این حالت خواهیم داشت [11]:
| (21) | z(t) = y1(t) |
ماتریس بهره K از رابطه زیر محاسبه می‌شود [11]:
| (22) | , |
P یک ماتریس حقیقی متقارن مثبت معین است که معادله جبری ریکاتی زیر را برآورده می‌کند [11]:
| (23) |
همچنین، مقدار شاخص عملکرد برابر است با [11]:
| (24) |
2- مطالعه میکروگرید شامل یک DGهمان‌طور که در شکل (3) نشان داده شده است، سیستم در نظر گرفته شده در این بخش شامل یک منبع DG است. میکروگرید در نقطه PCC به شبکه اصلی متصل شده است. منبع DG دارای بار محلی است که شامل بارهای نامتعادل و غیرخطی است.در حالت متصل به شبکه بار مشترک به‌طور کامل توسط شبکه تغذیه می‌شود و توان مورد نیاز بار محلی به‌طور مشترک توسط منبع DG-1 و شبکه تأمین می‌شود. در حالت جزیره‌ای، یعنی حالتی که CB-1، باز و میکروگرید از شبکه اصلی تأمین برق جدا شود، DG-1 باید تمامی توان مورد نیاز را برای بار محلی و نیز بار مشترک فراهم آورد. در بخش‌های بعدی هر دو حالت را در نظر خواهیم گرفت. شکل(3): دیاگرام تک‌خطی سیستم شامل میکروگرید، شبکه برق اصلی و بارهای محلی و مشترک2-1- تغذیه بار محلی به‌صورت مشترک با شبکه برق اصلی (حالت متصل به شبکه)در این بخش تغذیه بار محلی DG به‌وسیلة اشتراک‌گذاری این بار توسط شبکه اصلی برق و منبع DG نشان داده خواهد شد. در این حالت به‌دلیل اینکه منبع DG بار مشترک را تأمین نمی‌کند، خواهیم داشت ]12[:
| (25) |
فرض کنید در این حالت بخواهیم با توجه به برنامه‌ریزی‌های انجام‌شده، سیستم به‌گونه‌ای رفتار کند که DG-1، 10 درصد از توان حقیقی و 20 درصد از توان راکتیو بار محلی خود را تأمین کند؛ بنابراین، در این حالت خواهیم داشت [12]:
| (26) |
داده‌های سیستم شامل مشخصات بارها، مبدل منبع ولتاژ و منابع ولتاژ در جدول (1) آمده‌اند. جدول (1): پارامترهای سیستم شکل (3)
|
کمیت‌های سیستم |
مقادیر |
|
فرکانس سیستم |
50 HZ |
|
ولتاژ منبع (Vs) |
|
|
امپدانس فیدر |
Rs=1.025Ω, Ls=57.75 mH |
|
بار محلی نامتعادل |
|
|
|
|
|
بار محلی غیرخطی |
یک یکسوساز سه‌فاز شامل بار مقاومتی - سلفی با مقادیر R=200Ω و L=200mH |
|
امپدانس بار مشترک |
|
|
امپدانس خط میکروگرید |
|
|
DGها و جبران‌سازها |
|
|
ولتاژ DC (Vdc1) |
3.5 KV |
|
نسبت ترانسفورماتور |
3KV/11KV, 0.5 MVA, 2.5% reactance (Lf) |
|
تلفات VSC |
1.5 Ω |
|
خازن فیلتر (Cf) |
50 μF |
فرض کنید در زمان t=0.5 (s) امپدانس بار مشترک به نصف مقدار اولیه کاهش یابد. شکل (4) توان‌های حقیقی و غیرفعال را برای منبع DG-1 نشان می‌دهد. ولتاژ نقطه (PCC1) و جریان ig1 در شکل (5) نشان داده شده‌اند. تخصیص بار در همان نسبت‌های دلخواه و نیز تعادل ولتاژ حتی پس از تغییر مقدار بار مشترک نشان‌دهندة یک عملکرد پایدار برای شبکه است. گفتنی است برای برقراری تعادل در شبکه و تخصیص بهینه بار از قانون کنترلی فیدبک حالت بیان‌شده در بخش‌های قبل استفاده شده است. قانون کنترلی سیستم مطابق بلوک دیاگرام نشان داده شده در شکل (3) است. (الف) توان حقیقی (MW) (ب) توان غیرفعال (MVAR)شکل (4): تخصیص توان‌های حقیقی و غیرفعال در حالت متصل به شبکه (برای سیستم شامل یک DG) شکل (5): (الف) ولتاژ سه‌فاز در نقطه PCC1 شکل (5): (ب) جریان تزریقی شبکه به نقطه PCC1 در این مثال به‌منظور عملکرد مناسب شبکه پارامترهای شاخص LQR به‌صورت زیر انتخاب می‌شوند:
| (27) | K= [8.1109 3.2170 0.1170]. |
که در آن، diag یک ماتریس قطری است. دو متغیر که برای کنترل مهم‌تر از همه هستند، جریان ic و ولتاژ vcf هستند. انتخاب ماتریس Q اهمیت این حالت‌ها را نشان می‌دهد. در شکل‌های صفحه قبل افت ولتاژ[xvi] و درنتیجه افت توان‌های حقیقی و راکتیو در لحظه تغییر بار یعنی زمان 5/0 ثانیه به وضوح دیده می‌شود. پس از مدت زمان کوتاهی (حدود 025/0 ثانیه)، ولتاژ، جریان و توان‌ها مجدداً به حالت تعادل بازمی‌گردند. در جدول (2) برخی از نتایج عددی حاصل از این نمودارها آمده‌اند. جدول (2): نتایج عددی حالت متصل به شبکه
|
شکل (4) |
توان حقیقی |
مقدار اولیه |
مقدار نهایی (MW) |
مقدار مینیمم (t= 0.5 s) |
ماکزیمم فروجهش (MW) |
|
PL1 |
0.7 |
0.7 |
0.56 |
0.14 |
|
|
PG1 |
0.628 |
0.63 |
0.44 |
0.17 |
|
|
P1 |
0.072 |
0.07 |
0.003 |
0.067 |
|
|
شکل (4) |
توان غیرفعال |
مقدار اولیه (MVAr) |
مقدار نهایی (MVAr) |
مقدار مینیمم (t= 0.5 s) |
ماکزیمم فروجهش (MVAr) |
|
QL1 |
0.76 |
0.76 |
0.62 |
0.14 |
|
|
QG1 |
0.61 |
0.61 |
0.5 |
0.11 |
|
|
QL1 |
0.15 |
0.15 |
0.12 |
0.03 |
|
|
شکل (5) |
افت ولتاژ در نقطه PCC1(%) |
دامنه ولتاژ ( KV) |
دامنه جریان (Amp) |
||
|
4.5% |
5 |
175 |
|||
2-2- بررسی تغذیه بار محلی در حالت جزیره‌ای
سیستم مطرح‌شده در بخش قبل را در نظر بگیرید. فرض می‌کنیم در زمان t=0.4 s، CB-1 باز شود و درحالی‌که بار مشترک همچنان به سیستم متصل است، حالت جزیره‌ای اتفاق بیفتد. در این حالت باید هم بار مشترک و هم بار محلی به‌طور کامل توسط واحد DG تأمین شود. در این حالت، خواهیم داشت [13،12،11،10]:
|
(28) |
|
در این حالت نیز پارامترهای شاخص LQR مشابه قسمت قبل انتخاب می‌شوند. در شکل (6)، توزیع توان‌های حقیقی و در شکل (7) ولتاژ نقطه PCC1 نشان داده شده است. به محض اینکه حالت جزیره‌ای‌شدن اتفاق می‌افتد، افت ولتاژ نقطه PCC1 موجب ایجاد افت اندکی در مقدار توان مورد نیاز بار محلی PL1 می‌شود. در همین حال، توان تولیدی توسط واحد DG-1 یعنی P1 به‌منظور تأمین توان مورد نیاز بار مشترک افزایش می‌یابد که در شکل با مقدار منفی PG1 نشان داده شده است. مقدار منفی برای توان PG1 در شکل (6) نشان‌دهندة شارش توان در جهتی مخالف با جهت نشان داده شده در شکل (3) است.
شکل (6): تخصیص توان‌های حقیقی در DG-1 در حالت جزیره‌ای (MW)
شکل (7): ولتاژ سه‌فاز در نقطه PCC1در زمان t=1 s، CB-1 دوباره بسته می‌شود و سیستم به حالت متصل به شبکه بازمی‌گردد. همان‌طور که در شکل (5) مشاهده می‌شود، پس از وصل دوبارة مدارشکن در لحظه 1 ثانیه، سیستم بلافاصله متعادل نمی‌شود و تا رسیدن توان‌ها به مقادیر اولیه مدت زمانی (حدود 17/0 ثانیه) وقت لازم است. علاوه بر آن، مقادیر توان‌های حقیقی در این لحظات، اندکی فراجهش و فروجهش نسبت به حالت‌های تعادل خود دارند. در جدول (3) برخی از ویژگی‌های این نمودارها آمده‌اند. جدول (3): نتایج عددی حالت جزیره‌ای
|
|
توان حقیقی |
مقدار اولیه (MW) |
مقدار میانی در حالت جزیره‌ای (MW) |
مقدار نهایی (MW) |
زمان بازیابی مجدد (ثانیه) (پس از اتصال دوباره CB-1) |
||
|
شکل |
PL1 |
0.7 |
0.63 |
0.695 |
0.17 |
||
|
(6) |
PG1 |
0.629 |
-0.13 |
0.625 |
0.16 |
||
|
|
P1 |
0.071 |
0.76 |
0.07 |
0.17 |
||
|
(7) |
افت ولتاژ (Voltage droop) در(%)PCC1 |
دامنه ولتاژ KV |
|||||
|
مقدار اولیه |
مقدار میانی (در حالت جزیره‌ای) |
||||||
|
3.9% |
5 |
4.8 |
|||||
3-2- بررسی تأثیر انتخاب پارامترهای شاخص LQR بر پاسخ گذرای شبکهبا انتخاب مناسب پارامترهای شاخص LQR یعنی ماتریس‌های وزن‌دهی Q و جریمه کنترلی R می‌توان تا حدی این افت ولتاژها و افت توان‌ها در لحظه تغییر مقدار بار مشترک را بهبود بخشید. به کمک قاعده بریسون می‌توانیم مقادیر زیر را برای Q و R انتخاب کنیم [13]:
| (29) |
با انتخاب این مقادیر برای پارامترهای شاخص LQR مقدار ماتریس بهره به این صورت به دست می‌آید:
K= [9.536 3.217 0.155]
(الف) توزیع توان حقیقی (MW)
(ب) توزیع توان راکتیو (MVAR)شکل (8): توزیع توان‌های حقیقی و راکتیو در DG-1 و DG-2 (حالت دوم)جدول (4)، برخی از نتایج عددی حاصل از شبیه‌سازی‌ها را برای حالت‌های اول و دوم انتخاب R و Q نشان می‌دهد. با توجه به این جدول و با استفاده از شکل موج‌های ولتاژ و توان، می‌توان مشاهده کرد که در حالت دوم ولتاژها و توان‌ها در لحظه کاهش بار مشترک نسبت به حالت اول افت کمتری خواهند داشت. جدول (4): (الف) نتایج عددی توان حقیقی
|
شکل (8) توزیع توان منابع DG-1 و DG-2 را برای این مقادیر ماتریس بهره نشان می‌دهد. در شکل (9)، شکل موج‌های ولتاژ در نقاط PCC1 و PCC2 و جریان ig1 نشان داده شده است.وان حقیقی |
مقدار اولیه (MW) |
مقدار نهایی (MW) |
مقدار کمینه (t= 0.5 s) |
بیشترین افت (MW) |
|||||
|
حالت |
1 |
2 |
1 |
2 |
1 |
2 |
1 |
2 |
|
|
شکل‌های (6) و (8) |
PL1
|
1.07 |
1.06 |
1.06 |
1.06 |
0.89 |
0.96 |
0.17 |
0.1 |
|
PG1
|
0.86 |
0.84 |
0.85 |
0.85 |
0.66 |
0.75 |
0.19 |
0.1 |
|
|
P1
|
0.21 |
0.22 |
0.21 |
0.21 |
0.13 |
0.19 |
0.08 |
0.02 |
|
جدول(4): (ب) نتایج عددی توان راکتیو
|
توان غیرفعال |
مقدار اولیه (MV Ar) |
مقدار نهایی (MV Ar) |
مقدار کمینه (t= 0.5 s) |
بیشترین افت (MV Ar) |
|||||
|
حالت |
1 |
2 |
1 |
2 |
1 |
2 |
1 |
2 |
|
|
شکل‌های (6) و (8) |
QL1 |
1.15 |
1.16 |
1.16 |
1.16 |
0.96 |
1.05 |
0.2 |
0.2 |
|
QG1 |
0.92 |
0.93 |
0.93 |
0.93 |
0.77 |
0.84 |
0.16 |
0.09 |
|
|
Q1 |
0.23 |
0.23 |
0.23 |
0.23 |
0.19 |
0.21 |
0.04 |
0.02 |
|
جدول (4): (پ) نتایج عددی ولتاژ و جریان
|
شکل‌های (5) و (7) |
|
افت ولتاژ PCC1 (%) |
افت ولتاژ PCC1 (%) |
دامنه ولتاژ (KV) |
دامنه جریان ig1 (Amps) |
||||
|
حالت |
1 |
2 |
1 |
2 |
1 |
2 |
1 |
2 |
|
|
|
13% |
6% |
12.8% |
6.2% |
6 |
6 |
210 |
210 |
|
شکل (9): (الف) ولتاژ سه‌فاز در نقاط PCC1 و PCC2 شکل (9): (ب) جریان تزریقی به نقطه PCC1 (حالت دوم)
4-2- راهکارهای پیشنهادی برای بهبود پاسخ گذرای میکروگرید
در بخش قبل پس از معرفی کنترلکنندة فیدبک حالت LQR، چگونگی تأثیر این کنترل‌کننده بر عملکرد میکروگرید‌ها بررسی شد. در آنجا دیدیم کنترلکنندة پیشنهادی در حفظ تعادل و پایداری میکروگرید در حالت متصل به شبکه عملکرد مطلوبی داشت؛ اما در هنگام وقوع جزیره‌ای‌شدن و اتصال مجدد به شبکه قدرت، تا رسیدن میکروگرید به حالت تعادل و رسیدن توزیع توان‌ها به حالت پایدار، مدت زمانی به طول می‌انجامد و مقادیر توان در این حالت دارای اندکی فراجهش و فروجهش خواهند بود. در اینجا برای بهبود پاسخ گذرای میکروگرید، دو راهکار را معرفی می‌کنیم. این راهکارها عبارت‌اند از:
روش اول: به دست آوردن پارامترهای شاخص LQR با استفاده از الگوریتم اجتماع ذرات؛
روش دوم: طراحی کنترلکنندة PID در حلقه فیدبک با استفاده از الگوریتم اجتماع ذرات.
1-2-4- الگوریتم بهینه‌سازی اجتماع ذرات[xvii]
در الگوریتم PSO تعدادی از موجودات وجود دارند که به آنها ذره گفته می‌شود و در فضای جست‌وجوی تابعی پخش شده‌اند که قصد کمینه کردن (یا بهینه‌کردن) مقدار آن را داریم. هر ذره مقدار تابع هدف را در موقعیتی از فضا محاسبه می‌کند که در آن قرار گرفته است. سپس با استفاده از ترکیب اطلاعات محل فعلی‌اش و بهترین محلی که در گذشته در آن بوده است و نیز اطلاعات یک یا چند ذره از بهترین ذرات موجود در جمع، جهتی را برای حرکت انتخاب می‌کند. همه ذرات جهتی برای حرکت انتخاب می‌کنند و پس از انجام حرکت، یک مرحله از الگوریتم به پایان می‌رسد. این مراحل چندین‌بار تکرار می‌شوند تا جواب مدنظر به دست بیاید. درواقع، انبوه ذرات که مقدار کمینة یک تابع را جست‌وجو می‌کنند، همانند دسته‌ای از پرندگان عمل می‌کنند که به‌دنبال غذا می‌گردند.
برای فرموله‌کردن الگوریتم PSO دو متغیر x و v به‌ترتیب با عنوان موقعیت و سرعت ذره تعریف می‌شوند. بهترین موقعیت ذره (ازلحاظ شایستگی در تابع هدف) با P_best و بهترین موقعیت بهترین ذره در کل گروه با عنوان g_best شناخته می‌شود [14]. برای اطمینان از همگرایی PSO استفاده از ضرایبی موسوم به ضرایب انقباض برای تنظیم مناسب پارامترهای PSO لازم است. روابط اصلاح‌شدة سرعت و موقعیت ذره را با توجه به ضرایب انقباض می‌توان به‌صورت زیر نوشت [14]:
|
(30) |
|
|
|
(31) |
+ |
|
d شمارنده تکرار، موقعیت ذره در تکرار d، موقعیت ذره در تکرار d+1، سرعت ذره در تکرار d، w وزن اینرسی و و ضرایب شتاب نام دارند.
در رابطه بالا، rand یک تابع تولید عدد تصادفی با توزیع یکنواخت در بازه [1 0] است. α تابعی از و است که توسط Clerk برای تقویت همگرایی و پایداری PSO پیشنهاد شد [14]
|
(32) |
|
|
|
(33) |
|
|
انتخاب مناسب w سبب تعادل در جست‌وجوی محلی و جست‌وجوی کل فضای مسئله خواهد شد. عموماً برای عملکرد بهینه الگوریتم، w به‌صورت دینامیک تغییر داده می‌شود [14]:
|
(34) |
|
نشان‌دهندة حد تعداد تکرار الگوریتم و iter شمارندة تکرار است. در روابط (1) و (2) برای جلوگیری از واگرایی الگوریتم و کوچک‌کردن گام‌های جست‌وجو، سرعت ذره با مقدار محدود می‌شود [14،13،12،11،10]:
|
(35) |
|
جست‌وجوی محلی را بهبود می‌دهد و درواقع فرایند یادگیری تدریجی ذره را مدل می‌کند. اگر بزرگ تنظیم شود، ممکن است ذره از نقطه بهینه جهش کند و اگر کوچک تنظیم شود، همگرایی الگوریتم کند می‌شود و ذره ممکن است در دام کمینه محلی بیفتد. معمولاً برای هر متغیر تصمیم بین 10 درصد تا 20 درصد محدوده متغیر مربوطه تنظیم می‌شود.
2-4-2- به دست آوردن پارامترهای شاخص LQR با استفاده از الگوریتم اجتماع ذرات
با انتخاب مناسب پارامترهای شاخص LQR می‌توان عملکرد گذرای سیستم را بهبود بخشید. برای این منظور، پارامترهای شاخص LQR برای میکروگرید شامل دو منبع تولید پراکنده را با استفاده از الگوریتم تکاملی اجتماع ذرات محاسبه می‌کنیم.
برای این منظور، فرض کنید پارامترهای شاخص LQR یعنی R و Q به‌صورت زیر در نظر گرفته شوند [15،14]:
|
(36) |
|
در این حالت، تابع هدف که باید مینیمم شود، به‌صورت زیر تعریف می‌شود [14]:
|
(37) |
|
در رابطه (35)، P1ss، PL1ss، PG1ss و P2ss، PL2ss و PG2ss مقادیر نهایی (ماندگار) به‌دست‌آمده برای توان‌های حقیقی در جدول (4) هستند. Q1ss، QL1ss، QG1ss و Q2ss، QL2ss و QG2ss نیز مقادیر نهایی (ماندگار) توان‌های راکتیو هستند. درواقع، این تابع هدف به‌گونه‌ای تعریف شده است که مقادیر توان‌ها کمترین میزان انحراف را از مقدار نهایی یا مقدار حالت ماندگار خود داشته باشند و در کمترین زمان ممکن به این حالت ماندگار برسند. در این حالت چون میکروگرید شامل دو منبع تولید پراکنده است، نیاز به حضور دو جبران‌ساز خواهیم داشت و درنتیجه باید برای هر جبران‌ساز یک کنترلکنندة LQR طراحی کنیم.
مراحل اجرای الگوریتم PSO برای حل مسئله به‌صورت زیر انجام می‌شود:
۱- موقعیت و سرعت اولیه ذرات در محدودة مجاز به‌صورت تصادفی تولید می‌شود.
۲- تابع هدف رابطه (35) محاسبه می‌شود.
3- P_best برای هر ذره برابر موقعیت اولیه‌اش در نظر گرفته می‌شود (در تکرار اول). همچنین، g_best یعنی برترین ذره در میان ذرات ازنظر میزان شایستگی در تابع هدف رابطه (7) تعیین می‌شود.
۴- سرعت و موقعیت جدید هر ذره با استفاده از روابط (1) و (2) محاسبه می‌شود.۵- تابع هدف کلی از (36) برای هر ذره محاسبه می‌شود.6- اگر مقدار تابع هدف کلی جدید برای هر ذره بهتر از مقدارش در P_best باشد، P_best جایگذین می‌شود. به همین نحو، چنانچه در کل جمعیت ذره‌ای بهتر ازنظر میزان شایستگی در تابع هدف پیدا شد، جایگزین g_best می‌شود.7- مراحل 4 تا 7 تا رسیدن به بیشینه تعداد تکرار مشخص‌شده (itermax) انجام می‌گیرد. جواب نهایی مسئله از تعیین g_best در تکرار آخر به دست می‌آید.مقادیر پارامترهای الگوریتم PSO برای تابع هدف رابطه (36) به‌صورت زیر انتخاب شده است:
| (38) | , Itermax=100 Popsize=40, Npar=50 |
پس از اجرای الگوریتم PSO روی تابع هدف رابطه (36) در طی 100 تکرار نتایج زیر برای پارامترهای شاخص LQR حاصل شد: جدول (5): نتایج به‌دست‌آمده برای پارامترهای شاخص LQR با استفاده از الگوریتم اجتماع ذرات
| R | Q3 | Q2 | Q1 | متغیر |
| 0.11 | 9.8950 | 2.1242 | 0.0012 |
در شکل (10) تخصیص توان‌های حقیقی و در شکل (11) ولتاژ نقاط PCC1 و PCC2 در هنگام استفاده از مقادیر به‌دست‌آمده در جدول (5) برای پارامترهای شاخص LQR نشان داده شده‌اند. جدول (6) نتایج عددی حاصل از این نمودارها را نشان می‌دهد. شکل (10): (الف) تخصیص توان‌های حقیقی در DG-1 شکل (10): (ب) تخصیص توان‌های حقیقی در DG-2 با مقایسة شکل موج‌های توان‌های حقیقی و ولتاژ به‌دست‌آمده در شکل‌های (20) و (21) با شکل موج‌های مشابه به‌دست‌آمده در بخش‌های قبل و بدون استفاده از الگوریتم اجتماع ذرات مشاهده می‌شود که در این حالت عملکرد گذرای شبکه بهبود چشمگیری خواهد داشت. همان‌طور که از این شکل‌ها و جدول (9) پیداست، در این حالت پس از وصل مجدد میکروگرید به شبکه اصلی مقادیر توان‌ها به‌سرعت به حالت پایدار (ماندگار) خود می‌رسند و دیگر از نوسانات توان خبری نخواهد بود. شکل (11): (الف) ولتاژ نقاط PCC1(با انتخاب مقادیر بهینه برای پارامترهای R و Q) شکل (11): (ب) ولتاژ نقاط PCC2(با انتخاب مقادیر بهینه برای پارامترهای R و Q) جدول (6): نتایج عددی توزیع بار مشترک توسط منبع تولید پراکنده با انتخاب پارامترهای بهینه
3- طراحی کنترلکنندةPID در حلقه فیدبک با استفاده از الگوریتم اجتماع ذرات برای اصلاح عملکرد گذرای سیستم، می‌توانیم یک کنترلکنندة PID در حلقه فیدبک طراحی کنیم. کارایی بالا، طراحی مناسب و سادگی پیاده‌سازی بعد از طراحی ازجمله عواملی‌اند که باعث کاربرد وسیع کنترلکنندة PID در صنعت شده است. استفاده از روش‌های هوشمند در طراحی پارامترهای این کنترل‌کننده در سال‌های اخیر افزایش چشمگیری داشته است. در گذشته این طراحی توسط روش‌های کلاسیک مانند مکان هندسی ریشه‌ها، روش زیگلر - نیکولز و ... انجام می‌گرفت. در این پایان‌نامه برای به دست آوردن ضرایب کنترلکنندة PID از روش بهینه‌سازی اجتماع ذرات کمک خواهیم گرفت. در شکل (12)، بلوک دیاگرام کنترلکنندة فیدبک حالت با حضور PID Controller نشان داده شده است [15،14].
|
x= [i1 icf]] Vcf]
|
|
x |
|
uc (k)= - [x(k)-xref(k)] |
|
xref= [i1ref icfref Vcfref]
|
|
PID controller |
|
xref |
|
_
|
|
|
|
|
شکل (12): بلوک دیاگرام کنترل فیدبک حالت با حضور کنترلکنندة PID
| (39) |
در این حالت قانون کنترلی به‌صورت زیر تعریف می‌شود:
| (40) |
uc (k)= - [x(k)-xref(k)] If uc (k)>h u= +1 If uc (k)>h u= - 1 If uc (k) <-h u=-1 |
متغیر u به‌صورت زیر تعریف می‌شود:
| (41) |
مقدار h برای حل این مسئله بهینه‌سازی برابر با 005/0 در نظر گرفته شده است. 1-3- طراحی کنترلکنندة PID برای میکروگرید شامل یک منبع تولید پراکندهدر این قسمت درصدد خواهیم بود که برای بهبود حالت گذرای سیستم، کنترلکنندة PID در حلقه فیدبک را طراحی کنیم. ضرایب کنترل‌کننده با استفاده از الگوریتم PSO محاسبه می‌شود و تابع هدف که باید مینیمم شود، به‌صورت زیر تعریف می‌شود:
| (42) |
در رابطه (9)، P1ss، PL1ss و PG1ss مقادیر نهایی (ماندگار) به‌دست‌آمده برای توان‌ها در جدول (3) هستند. Q1ss، QL1ss و QG1ss نیز مقادیر نهایی (ماندگار) توان‌های راکتیو هستند. درواقع، این تابع هدف به‌گونه‌ای تعریف شده است که مقادیر توان‌ها کمترین میزان انحراف را از مقدار نهایی یا مقدار حالت ماندگار خود داشته باشند و در کمترین زمان ممکن به این حالت ماندگار برسند. مقادیر پارامترهای الگوریتم PSO برای تابع هدف رابطه به‌صورت زیر انتخاب شده‌اند:
| (43) | Popsize=40, Npar=50, Itermax=100 |
پس از اجرای الگوریتم PSO روی تابع هدف نشان داده شده در رابطه (41) در طی 100 تکرار نتایج زیر برای ضرایب کنترلکنندة PID حاصل شد. جدول (7): پارامترهای بهینه به‌دست‌آمده برای کنترلکنندة PID با استفاده از الگوریتم PSO
| Kp | Ki | Kd |
| 1.8722 | 1.5880 | 0.0220 |
شکل (13) نمودار همگرایی تابع هدف را نشان می‌دهد. در شکل (14) تخصیص توان حقیقی و در شکل (15) ولتاژ نقطه PCC1 با وجود کنترلکنندة PID نشان داده شده است. جدول (8) نتایج عددی حاصل از این نمودارها را نشان می‌دهد. شکل (13): نمودار همگرایی الگوریتم PSO برای تابع هدف رابطه (9) شکل (14): توزیع توان‌های حقیقی DG-1 باوجود کنترلکنندة PID (به کمک الگوریتم PSO) شکل (15): ولتاژ سه‌فاز در نقطه PCC1 باوجود کنترلکنندة PID (به کمک الگوریتم PSO)با توجه به شکل (15-13) و اطلاعات نشان داده شده در جدول (11)، مشاهده می‌شود با حضور کنترلکنندة PID پاسخ گذرای سیستم تا حد مناسبی بهبود یافته است و پس از وصل شبکه، توزیع توان در مدت زمان کمی به حالت اولیه خود بازمی‌گردد. جدول (8): نتایج عددی حالت جزیره‌ای با حضور کنترلکنندة PID (با استفاده از روش PSO)
|
شکل (14) |
توان حقیقی |
مقدار اولیه |
مقدار میانی در حالت جزیره‌ای MW |
مقدار نهایی (MW) |
زمان بازیابی مجدد (ثانیه) پس از وصل مجدد شبکه |
|
|
PL1 |
0.7 |
0.64 |
0.7 |
0.03 |
||
|
PG1 |
0.63 |
-0.14 |
0.63 |
0.03 |
||
|
P1 |
0.07 |
0.76 |
0.07 |
0.04 |
||
| شکل(15) | افت ولتاژ در(%)PCC1
)Voltage droop( |
دامنه ولتاژ KV |
||||
|
مقدار اولیه |
مقدارمیانی (در حالت جزیره‌ای) |
|||||
|
2.5% |
5 |
4.89 |
||||
با مقایسة نتایج جدول (8) با نتایج حاصل از تخصیص توان در حالت جزیره‌ای به کمک کنترلکنندة LQR که در جدول (3) آمده است، می‌توان دید باوجود کنترلکنندة PID، رفتار گذرای سیستم تا حد زیادی بهبود می‌یابد و پس از اتصال دوباره میکروگرید به شبکه اصلی برق (در زمان 1 ثانیه)، توزیع توان به‌سرعت به حالت پایدار خود می‌رسد و دیگر از نوسانات توان خبری نخواهد بود.
1-3- عملکرد روش پیشنهادی در شرایط مختلف بهره‌برداری و مقایسه عملکرد این روش با روش‌های دیگر
در جدول 9، نتایج حاصل از مقایسة روش پیشنهادی با روش کنترل LQR با پارامترهای ثابت برای سه حالت بارگذاری مختلف نشان داده شده‌اند. همان‌طور که مشاهده می‌شود، روش پیشنهادی در هر سه حالت بارگذاری دارای مقادیر کمتری برای فراجهش توان، نوسان ولتاژ و زمان بازیابی توان است. این امر نشان می‌دهد روش پیشنهادی می‌تواند پاسخ گذرای میکروگرید را در شرایط مختلف بهره‌برداری به‌طور چشمگیری بهبود بخشد.
جدول (9): مقایسه عملکرد روش پیشنهادی با روش کنترل LQR با پارامترهای ثابت برای سه حالت بارگذاری مختلف
|
پارامتر |
حالت بارگذاری کم |
حالت بارگذاری متوسط |
حالت بارگذاری زیاد |
|
فراجهش توان پیشنهادی (MW) |
0.0001 |
0.0001 |
0.0001 |
|
فراجهش توان LQR (MW) |
.0002 |
.0002 |
0.0002 |
|
نوسان ولتاژ پیشنهادی (V) |
0.002 |
0.002 |
0.002 |
|
نوسان ولتاژ LQR (V) |
0.003 |
0.003 |
0.003 |
|
زمان بازیابی توان پیشنهادی (s) |
0.01 |
0.01 |
0.01 |
|
زمان بازیابی توان LQR (s) |
0.02 |
0.02 |
0.02 |
در جدول 10، نتایج حاصل از این تحقیق با نتایج حاصل از چند مقالة دیگر درزمینة کنترل میکروگرید مقایسه شده‌اند. همان‌طور که مشاهده می‌شود، روش پیشنهادی در مقایسه با روش‌های دیگر، مقادیر کمتری برای فراجهش توان، نوسان ولتاژ و زمان بازیابی توان دارد. این امر نشان می‌دهد روش پیشنهادی می‌تواند پاسخ گذرای میکروگرید را به‌طور چشمگیری بهبود بخشد.
جدول (10): مقایسة عملکرد روش پیشنهادی با روش‌های دیگر درزمینة کنترل میکروگرید
|
مقاله |
روش کنترل |
فراجهش توان (MW) |
نوسان ولتاژ (V) |
زمان بازیابی توان (s) |
|
این تحقیق |
پیشنهادی |
0.0001 |
0.002 |
0.01 |
|
[1] |
LQR با پارامترهای ثابت |
0.0002 |
0.003 |
0.02 |
|
[2] |
LQR با پارامترهای متغیر |
0.0002 |
0.003 |
0.02 |
|
[3] |
PID |
0.0003 |
0.004 |
0.03 |
4- نتیجه‌گیری
در این مقاله پس از معرفی دو ساختار مختلف برای میکروگرید، چگونگی تأثیر کنترلکنندة پیشنهادی بر عملکرد این میکروگرید‌ها در حالت‌های جزیره‌ای و متصل به شبکه بررسی شده است. پس از شبیه‌سازی مشاهده کردیم که کنترلکنندة LQR در هر دو حالت متصل به شبکه و جزیره‌ای، تعادل سیستم را تا حد مناسبی تأمین می‌کند و تخصیص توان دلخواه را در میکروگرید برقرار می‌کند؛ اما در هنگام جدایی شبکه و وقوع حالت جزیره‌ای و رفع این حالت (اتصال مجدد به شبکه قدرت)، تا رسیدن سیستم به شرایط پایدار، مدت زمانی طول کشیده است و توان‌های حقیقی و راکتیو در این حالت‌ها دارای مقادیری فراجهش و فروجهش نسبت به حالت پایدار خود بودند. در ادامه، به‌منظور اصلاح عملکرد گذرای میکروگرید‌ها در هنگام وقوع حالت جزیره‌ای، دو راهکار کنترلی معرفی شدند که عبارت بودند از: 1- به دست آوردن پارامترهای شاخص LQR با استفاده از الگوریتم اجتماع ذرات و 2- طراحی کنترلکنندةPID در حلقه فیدبک با استفاده از الگوریتم اجتماع ذرات. در راهکار اول، پس از محاسبه پارامترهای شاخص LQR با استفاده از الگوریتم اجتماع ذرات برای میکروگرید شامل دو منبع تولید پراکنده، مشاهده کردیم که با انتخاب مناسب پارامترهای شاخص LQR می‌توان عملکرد گذرای میکروگرید را در هنگام وقوع حالت جزیره‌ای بهبود بخشید. در راهکار دوم، کارایی بالا، طراحی مناسب و سادگی پیاده‌سازی بعد از طراحی ازجمله عواملی بودند که باعث استفاده از کنترلکنندة PID برای رسیدن به هدف مدنظر شدند. به‌منظور بررسی صحت عملکرد کنترل‌کننده، کارایی آن در بهبود کیفیت توان میکروگرید‌های شامل 1 منبع تولید پراکنده و 2 منبع تولید پراکنده سنجیده شد. ضرایب کنترلکنندة PID با استفاده از الگوریتم تکاملی اجتماع ذرات محاسبه شدند و درنهایت شکل موج‌های ولتاژ و تخصیص توان حقیقی برای هر دو حالت ترسیم شد. با مقایسة این شکل موج‌ها با شکل‌های به‌دست‌آمده با حضور کنترلکنندة LQR، مشاهده کردیم که کنترلکنندة PID پیشنهادی پاسخ گذرای سیستم و کیفیت توان شبکه را تا حد زیادی بهبود می‌بخشد. ذکر این نکته شایان توجه است که کنترلکنندة LQR پیشنهادی در برقراری تعادل و پایداری میکروگرید در حالت متصل به شبکه عملکرد موفقی داشت.
[1] تاریخ ارسال مقاله: 19/10/1401
تاریخ پذیرش مقاله: 21/12/1402
نام نویسندۀ مسئول: حمیدرضا اکبری
نشانی نویسندۀ مسئول: ایران، یزد، دانشگاه آزاد اسلامی، واحد یزد، گروه مهندسی برق
[i] Seamless Grid-Connected Inverter
[ii] State feedback controller
[iii] Linear quadratic regulator
[iv] Particle Swarm Optimization
[v] State Feedback LQR
[vi] The controlled output
[vii] weighting matrix
[viii] control penalty
[ix] Bryson's rule
[x] positive-definite
[xi] Algebraic Riccati Equation
[xii] Gain matrix
[xiii] Discrete time linear quadratic regulator
[xiv] Weighting matrix
[xv] The measured output
[xvi] Voltage Droop
[xvii] Particle Swarm Optimization (PSO)
)