Authors
1 Department of Electrical Engineering, Kermanshah Branch, Islamic Azad University, Kermanshah, Iran
2 Department of Electrical Engineering, Engineering Faculty, Razi University, Kermanshah, Iran
Abstract
Keywords
نگرانی‌های زیست محیطی و افزایش مصرف انرژی موجب شده است توجه به منابع انرژی تجدیدپذیر افزایش یابد. در میان منابع تجدیدپذیر، انرژی بادی و خورشیدی نفوذ زیادی در سیستم‌های قدرت سراسر جهان یافته‌اند. ورود انرژی بادی و خورشیدی در سیستم قدرت، به دلیل ماهیت غیرقطعی آنها، چالش‌های عمده‌ای را در برنامه‌ریزی و بهره‌برداری سیستم‌های قدرت سبب شده است [1].
با توجه به بحران جهانی در افزایش آلودگی‌های زیست‌محیطی برای تولید انرژی الکتریکی، استفاده از منابع انرژی‌های تجدیدپذیر به‌منزلۀ راه‌حلی مؤثر، شایان توجه فراوانی قرار گرفته است؛ چون این منابع نسبت به منابع فسیلی، رایگان بوده و دارای عدم آلودگی‌اند. با وجود این، مشکل اصلی به‌کارگیری این منابع وجود عدم‌قطعیت آنها در مقدار توان تولیدی است که چالش‌های بزرگی را در برنامه‌ریزی و بهره‌برداری سیستم قدرت، فراروی برنامه‌ریزان و بهره‌برداران سیستم قدرت قرار داده است. بنابراین ارائۀ مدل‌های پیش‌بینی دقیق رفتار این منابع اهمیت فراوانی دارد [2]. دربارۀ مسئلۀ پخش توان بهینه[1]، مطالعات زیادی دربارۀ مدل پخش توان بهینۀ سیستم قدرت و روش حل با توجه عدم‌قطعیت‌ها انجام شدهاند [4,3]. روش‌های احتمالاتی عوامل تصادفی مانند ساختار توپولوژی شبکه، بارها و خروجی ژنراتورها و درنتیجه، مقدار احتمالی ولتاژ و توان شاخه‌ها در شبکه‌های قدرت را نشان می‌دهد [5].
روش‌های احتمالاتی عموماً به دو دستۀ روش‌های عددی و روش‌های تحلیلی تقسیم شدهاند. ازجمله روش‌های عددی روش شبیه‌سازی مونت کارلو[2] است که اجرای آن ساده و دقت محاسبات آن زیاد است؛ اما زمان‌بر است. همچنین، روش‌های تحلیلی مبتنی بر خطی‌سازی‌اند و دقت محاسباتی آنها کمتر است. یکی از روش‌های استفاده‌شده در روش تحلیلی برای حل پخش توان احتمالاتی با توجه به تقاضای بار، استفاده از سری گرام چالیر[3] است که در مرجع [6] به آن اشاره شده است. ازجمله روش‌های تحلیلی دیگر روش تخمین نقطه‌ای[4] برای حل پخش توان بهینۀ احتمالاتی[5] است که در مراجع [9-7] پیشنهاد شده است. این روش، محاسبات کمتری دارد و برای محاسبات مسائل مقیاس بزرگ مشکلاتی ندارد [1].
روش تبدیل بی بو[6] رویۀ جدید مدسازی عدم‌قطعیت‌ها است که عملکرد مناسبی در حالت غیرخطی و نسبت به روش تخمین دارد که در مرجع [10] به آن اشاره شده است. این روش دارای دقت زیاد محاسباتی بوده و به‌طور چشمگیری سریع‌تر از روش شبیه‌سازی مونت کارلو است. یکی از ویژگی‌های خوب روش تبدیل بی بو، امکان مدلسازی متغیرهای عدم‌قطعیت است.
روش نقطۀ داخلی[7] برای حل مسئلۀ بهینه‌سازی براساس ساخت یک مسیر از یک نقطۀ اولیه داده شده است تا راه‌حل نهایی را در یک منطقه امکان‌پذیر جستجو کند. این روش برای نخستین‌بار در سال 1984 پیشنهاد شده [11] و از آن زمان تا کنون در پژوهش‌های مختلفی بررسی شده است [12]. این روش در مسائل مختلف حوزۀ بهره‌برداری مانند پخش بار اقتصادی [13]، پخش توان بهینۀ محدود به قیود امنیتی [14]، ارزیابی حاشیه‌های بارگذاری [15] و کاهش تلفات [16] به‌کار گرفته شده است. روش نقطۀ داخلی، بیشتر برای حل الگوریتم پخش توان بهینۀ خطی‌شده برای پشتیبانی برنامه‌ریزی توسعۀ خطوط انتقال استفاده می‌شود [17]. ضعف عمدۀ روش نقطۀ داخلی هنگامی است که محدودیت‌های فعال، مانع از رسیدن نقطۀ مدنظر به صفر یا مرز پذیرفتنی شوند؛ بنابراین چالش بزرگ در به‌کارگیری این روش در حل مسئلۀ پخش توان بهینه، طراحی روش همگرایی قوی برای نشان‌دادن محدودیت‌های فعال برای جلوگیری از ناکارآمدی این روش است. برای مقابله با این مشکل در مرجع [18] راهکارهایی برای اصلاح بیان شدهاند. در مرجع [19] روش تابع جریمه برای حل مسئلۀ بهینه‌سازی غیرخطی به کار گرفته شده که در آن از تقریب درجه دوم استفاده شده است که موجب حل مسئلۀ روش نقطۀ داخلی می‌شود.
روش شبیه‌سازی مونت کارلو به‌طور گسترده‌ای برای حل مسئلۀ پخش توان احتمالی در حضور عدم‌قطعیت‌ها به‌عنوان متغیرهای تصادفی ورودی به‌منظور ایجاد فرآیند احتمالاتی مطرح شده است [5]. این روش برای دستیابی به جواب هدف نهایی باید ده‌ها هزار نتیجۀ نمونه‌برداری در مدت زمان زیادی را محاسبه کند؛ بنابراین، برای محاسبۀ پخش توان بهینۀ منابع تجدیدپذیر با مقیاس بزرگ مناسب نیست؛ زیرا مدت زمان انجام محاسبات بسیار زیاد است. این روش با توجه به تابع توزیع احتمال ورودی، شبیه‌سازی نمونه‌ها را برای سناریوهای مختلف در سیستم قدرت تولید می‌کند و معادلات پخش توان قطعی را برای ولتاژ شین، زوایای فاز، پخش توان اکتیو و راکتیو حل می‌کند. به این ترتیب، توابع توزیع خروجی‌های مختلف محاسبه می‌شوند.
پیشنهاد دیگر، استفاده از روش تخمین نقطه‌ای با تقریب گشتاورهای متغیر است. روش تخمین نقطه‌ای، مقدار میانگین و انحراف استاندارد هر متغیر سیستم قدرت را به‌طور مؤثر به دست می‌آورد؛ در حالی که روش شبیه‌سازی مونت کارلو نیاز به آگاهی از توابع توزیع احتمال دارد. روش تخمین نقطه‌ای به دو دلیل به‌طور گسترده‌ای برای حل پخش توان احتمالی استفاده می‌شود: خطی‌سازی معادلات پخش بار ضروری نیست و ازنظر محاسبات کارآمد است. رایج‌ترین روش تخمین نقطه‌ای براساس عملکرد هانگ[8] است؛ با این حال، هار[9] و مورالس[10] برخی از اصلاحات را در روش تخمین نقطه‌ای انجام داده‌اند [5].
با نفوذ منابع تجدیدپذیر به سیستم‌های قدرت مدرن، لازم است روش‌های پخش توان بهینۀ احتمالاتی مناسب برای مدلسازی عدم‌قطعیت‌ها در معادلات پخش توان به کار گرفته شوند. در پخش توان بهینۀ احتمالاتی، متغیرهای تصادفی ورودی در نظر گرفته می‌شوند که هدف اصلی این روش، دستیابی به اطلاعات آماری خروجی است. تا کنون الگوریتم‌های مختلفی برای محاسبات پخش توان بهینۀ احتمالاتی ارائه شدهاند که به روش‌های تحلیلی و روش‌های عددی طبقه‌بندی می‌شوند.
یکی از مهم‌ترین موضوعات در سیستم‌های مدرن قدرت امروزی، تولید انرژی الکتریکی برای سیستم‌های قدرت به‌منظور بهینه‌سازی هزینۀ تولیدی برای واحدهای فعال موجود در شبکۀ قدرت است. پخش بار، ابزار مناسب برای تنظیم میزان توان واحدهای تولیدی است که محدودیت‌های غیرخطی آنها و شبکۀ قدرت را در نظر می‌گیرد. پخش بار اقتصادی، مسئله‌ای غیرخطی، غیرمحدب و چالش‌برانگیز است که برای حل آن، با توجه به مشخصات پیچیدۀ موجود در مسئله، از الگوریتم‌های ابتکاری استفاده می‌شود. در این مقاله، مسئلۀ پخش بار اقتصادی با محدودیت‌های غیرخطی، به مسئلۀ بهینه‌سازی تبدیل شده و با استفاده از الگوریتم بهینه‌سازی یادگیری ردیابی بازگشتی، به حل آن پرداخته شده است [20].
برای بررسی عدم‌قطعیت در سیستم تولید باد از توزیع ویبول[11] [21] و توزیع ریلی[12] [22] استفاده می‌شود. سپس توان اکتیو تزریق‌شده از سیستم تولید باد با استفاده از رابطۀ سرعت و قدرت ژنراتور محاسبه می‌شود. در مراجع [24,23] همبستگی بین توربین بادی در نظر گرفته می‌شود. همچنین در برخی از پژوهش‌ها، عواملی ازقبیل عدم‌قطعیت ژنراتور معمولی نیز شایان توجه قرار گرفته است. از بحث فوق دیده می‌شود در ادبیات پژوهش، تقریباً همۀ انواع عدم‌قطعیت‌ها در پخش توان بهینۀ احتمالاتی با سیستم تولید باد در نظر گرفته شدهاند. در سیستم تولید باد، توان اکتیو تزریقی، عدم‌قطعیت است و هیچ مدل دقیقی از سیستم تولید باد در هیچ کدام از پژوهش‌ها شایان توجه قرار نگرفته است [25]. برای بررسی عدم‌قطعیت ناشی از انرژی خورشیدی از توزیع بتا[13] استفاده می‌شود که با استفاده از این توزیع، توان تزریقی سلول خورشیدی محاسبه می‌شود [26].
در حال حاضر برای مدلسازی عدم‌قطعیت بار مبتنی بر رویکرد مشاهدات از داده‌های سیستم قدرت برای انواع دوره‌های زمانی مختلف، جمع‌آوری و تجزیه و تحلیل می‌شوند تا اطلاعات مربوطه در شبیه‌سازی‌ها اعمال شوند. موضوع بار یا پیش‌بینی منابع انرژی موضوع مهمی است که روش‌های متعددی برای پیش‌بینی تغییرات بار وجود دارند که بیشتر آنها مبتنی بر تحلیل رفتار کاربران یا بررسی داده‌های تاریخی‌اند که برای بهبود تحلیل احتمالی استفاده می‌شوند. پس برای بررسی عدم‌قطعیت ناشی از تغییرات بار از توزیع نرمال[14] استفاده می‌شود [27].
نویسندگان در مرجع [28] توان بهینۀ احتمالاتی براساس مقایسۀ بین روش‌های تخمین دو و سه نقطه‌ای با روش شبیه‌سازی مونت کارلو بررسی شده است. این مقایسه براساس بررسی توان اکتیو تولیدی واحدهای تولید با در نظر گرفتن عدم‌قطعیت‌های توان تزریقی باد و خورشید بوده و تحلیل نتایج روی سیستم استاندارد 98 شینهPEGASE اروپا نشان داده شد که در یک شبکۀ مقیاس بزرگ روش تخمین سه نقطه‌ای در رسیدن به جواب بهینه، چه ازنظر زمان اجرا و چه ازنظر هزینۀ بهینۀ بهره‌برداری، کارایی بهتری دارد.
نوآوری این مقاله آن است که مسئلۀ پخش توان بهینۀ احتمالاتی در یک شبکۀ استاندارد 118 شینه IEEE با نرم‌افزارهای MATLAB و DigSilent با روش‌های عددی و تحلیلی متفاوت انجام شده است و نتایج آنها به تفصیل با همدیگر مقایسه و بررسی شدهاند. این روش‌ها ازنظر دقت و زمان همگرایی حصول جواب بهینه، تحلیل شدهاند. همچنین نتایج پروفیل ولتاژ تولیدی ژنراتورها و توان اکتیو تولیدی ژنراتورها مقایسه شدهاند. در انتها، تمامی نتایج به‌دست‌آمده ازنظر مقدار تابع هدف مسئله مقایسه و بررسی شدهاند که هزینۀ بهره‌برداری بر حسب دلار بر ساعت و زمان مورد نیاز برای رسیدن به جواب بهینه است.
ساختار مقاله به شرح زیر است: در بخش دوم، فرمول‌بندی مسئله ازنظر پخش توان بهینه، محدودیت‌ها، تابع هدف، انواع روش‌های پخش توان احتمالاتی و در آخر همین بخش، مدل‌های تابع توزیع بیان شدهاند. در بخش سوم، شبکۀ نمونه استاندارد 118 شینه IEEE با استفاده از روش‌های یادشده در نرم‌افزارهای MATLAB و DigSilent، شبیه‌سازی و نتایج به‌دست‌آمده بررسی شدهاند. سرانجام در بخش چهارم، نتایج حاصل‌شده بیان شدهاند.
پخش توان بهینۀ مسئله، بهینه‌سازی غیرخطی است که هدف اصلی آن حداقل‌کردن توان در سیستم قدرت با در نظر گرفتن تعدادی از قیدهای برابری و نابرابری است. ازنظر ریاضی، پخش توان بهینه به‌صورت معادله (1) نشان داده می‌شود [29]:
|
(1) |
|
در اینجا، بردار متغیرهای کنترل یا مستقل، بردار متغیرهای حالت یا وابسته‌اند. تابع هدف پخش توان بهینه، مجموعه‌ای از محدودیت‌های نابرابری و مجموعه‌ای از محدودیت‌های برابری است.
مجموعه‌ای از متغیرهایی که پخش توان بهینه را در شبکۀ قدرت کنترل می‌کنند، در فرم برداری این متغیرها به‌صورت معادله (2) نشان داده می‌شوند:
|
(2) |
|
در اینجا، توان اکتیو ژنراتور ام است. انتخاب شین مبنا اختیاری است و این شین می‌تواند هرکدام از شین‌های ژنراتورها باشد. دامنه ولتاژ شین ام، تپ ترانسفورماتور شاخه ام، جبران‌ساز موازی شین ام، ، و به‌ترتیب نشان‌دهندۀ ژنراتور، جبران‌سازی توان راکتیو موازی و ترانسفورماتور است.
وضعیت سیستم قدرت با متغیرهای حالت تعریف می‌شود که این متغیرها با بردار نشان داده می‌شوند.
|
(3) |
|
در اینجا، توان اکتیو ژنراتور شین اسلک، توان راکتیو ژنراتور متصل‌شده به شین ، ولتاژ شین ام بار، بارگیری خط ام با نشان داده می‌شوند. و به‌ترتیب شمارۀ شین‌های بار و خطوط انتقال‌اند.
همان‌طور که پیش از این ذکر شد مسئلۀ پخش توان بهینه دارای هر دو نوع محدودیت‌های برابری و نابرابری است. این محدودیت‌ها در این بخش از مقاله به‌صورت جداگانه توضیح داده می‌شوند.
محدودیت‌های برابری
در پخش توان بهینه، معادلات تعادل توان جز محدودیت‌های برابری‌اند که به‌صورت معادلات (4) و (5) نشان داده می‌شوند:
|
(4) |
|
|
(5) |
|
که نشان‌دهندۀ زاویۀ اختلاف فاز ولتاژ بین شین‌های و است. شمارۀ شین‌ها، و به‌ترتیب تقاضای بار اکتیو و راکتیو، کندوکتانس منتقل‌شده و سوسپتانس بین شین و است.
محدودیت‌های نا برابری
محدودیت‌های نابرابری در پخش توان بهینه منعکس‌کنندۀ محدودیت‌های بهره‌برداری تجهیزات موجود در سیستم قدرت و همچنین محدودیت‌های اعمال‌شده بر خطوط و شین‌های بار برای تضمین امنیت سیستم قدرت‌اند.
محدودیت‌های نابرابری به‌صورت زیر معرفی می‌شوند [30]:
|
(6) |
|
در اینجا و به‌ترتیب محدودیت‌های حداقل و حداکثر ولتاژ شین ژنراتوری ام است. این محدودیت‌ها بیشتر با استانداردهای مرتبط ارائه می‌شوند. افزایش یا کاهش بیش از حد ولتا، به واردشدن آسیب‌های جدی ازجمله بی‌ثباتی ولتاژ در سیستم قدرت منجر می‌شود که ازنظر فنی و اقتصادی پذیرفتنی نیست.
نقاط ممنوعه در بهره‌برداری از ژنراتورها
نقاط ممنوعۀ ژنراتور نقاطی‌اند که در صورت بهره‌برداری از واحد در این نقاط، ممکن است به آنها آسیب‌های مکانیکی وارد شود. این نقاط برای واحد ام به‌صورت معادله (7) مدل می‌شوند:
|
(7) |
|
در اینجا، و به‌ترتیب مقادیر حداکثر و حداقل محدوده امین نقاط ممنوعه در واحد ام است. و بیان‌کنندۀ حداقل و حداکثر محدودۀ توان اکتیو در ژنراتور ام است.
|
(8) |
|
در اینجا، و به‌ترتیب حداقل و حداکثر محدودۀ توان راکتیو در ژنراتور ام است.
|
(9) |
|
|
(10) |
|
در اینجا، و به‌ترتیب محدودیت‌های حداقل و حداکثر ام تپ چنجر ترانسفوماتور است. و به‌ترتیب حداقل و حداکثر محدودیت امین شیفت‌دهندۀ فاز ترانسفوماتور است. و شماره تغییرات تپ ترانسفورماتور و شیفت‌دهندۀ فاز نصب‌شده در سیستم قدرت است.
|
(11) |
|
در اینجا، و به‌ترتیب حداقل و حداکثر محدودیت امین جبران‌ساز موازی است. شمارۀ جبران‌سازهای خازنی متصل‌شده به سیستم قدرت است.
|
(12) |
|
در اینجا، و محدودیت حداکثر امین خط انتقال است. شمارۀ خط انتقال در سیستم قدرت است.
قابلیت انتقال توان یک خط انتقال به حد پایداری آن محدود می‌شود. حد حرارتی یا پایداری با ظرفیت حمل جریان هادی، براساس اطلاعات ارائه‌شدۀ سازندگان، مشخص می‌شود. اگر ظرفیت حمل جریان با جریان حرارتی نشان داده شود، حد بارگذاری حرارتی خط برابر است با:
|
(13) |
|
معادله (13) توان حقیقی را در خط بدون تلفات نشان می‌دهد.
حداکثر توان انتقالی در 90 به دست می‌آید. در عمل زاویۀ بار برای خط به‌تنهایی محدود به 30 تا 45 درجه است. در خط بدون تلفات است و رابطه (14) نوشته می‌شود:
|
(14) |
|
در اینجا، دو جمله نخست داخل پرانتزها، ولتاژهای مبنای واحدند که با و نشان داده شدهاند و جمله سوم، حد حرارتی است. معادله (14) به‌صورت (15) نوشته می‌شود.
|
(15) |
|
|
(16) |
|
مجموع راکتانس القایی استاتور و روتور امین توربین بادی است. و به‌ترتیب توان راکتیو و ولتاژ شین محل اتصال امین ژنراتور بادی است.
|
(17) |
|
یکی از مهم‌ترین هدف‌های پخش توان بهینه، حداقل‌کردن هزینۀ سوخت ژنراتورها و توان تولیدی ژنراتورها است. برای دستیابی به این هدف از معادلۀ درجه دوم (18) استفاده می‌شود [29]:
|
(18) |
|
در اینجا، ، و ضرایب هزینه‌ای ژنراتور ام تولیدکنندۀ توان خروجی هستند.
به دلیل گستردگی سیستم قدرت در حوزه‌های مکانی و زمانی مختلف، پارامترهای متعددی بر آن تأثیر گذارده‌اند و ممکن است عدم‌قطعیت‌هایی را به سیستم اعمال کنند. برای ارزیابی عدم‌قطعیت یک سیستم از تئوری احتمال یعنی از روش‌های احتمالاتی استفاده می‌شود که هدف اصلی آنها تعیین وضعیت سیستم به‌منزلۀ یک تابع از متغیرهای ورودی نامعین به‌صورت معادله (19) است:
|
(19) |
|
شامل تقاضای ورودی نامعین است. توان تولیدی با منابع انرژی توزیع‌شده مانند باد، انرژی خورشیدی و ... به‌صورت معادله (20) بیان شده است:
|
(20) |
|
همان‌طور که پیش از این ذکر شد پخش توان بهینه دارای برخی متغیرهای حالت و کنترل است؛ ازاین‌رو، در پخش توان بهینۀ احتمالاتی بردار عدم‌قطعیت خروجی است که به‌صورت معادله (21) محاسبه می‌شود:
|
(21) |
|
عدم‌قطعیت متغیرهای ورودی باعث می‌شود متغیرهای خروجی نیز نامعین باشند؛ به این معنی که حتی اگر یک متغیر ورودی نامعین باشد، تمام متغیرهای خروجی نامعین خواهند شد. پارامترهای استفاده‌شده در معادلات (19) تا (21) در فهرست علایم تشریح شدهاند.
شبیه‌سازی مونت کارلو، روش شبیه‌سازی مبتنی بر ترکیب اعداد تصادفی برای حل مسائل با وجود عدم‌قطعیت‌ها و حالت‌های احتمالاتی است. این روش محاسبات مدل قطعی را با استفاده از مجموعه‌ای عدد تصادفی به شکل ورودی اجرا می‌کند. از این روش، بیشتر در مواردی استفاده می‌شود که مدل دارای پیچیدگی‌های زیادی ازجمله غیرخطی‌بودن یا بیش از دو پارامتر عدم‌قطعیت باشد. روش شبیه‌سازی مونت کارلو از تکرار چندین باره حل مسئله برای رسیدن به جواب نهایی استفاده می‌کند؛ در حالی که این روش نتایج دقیقی از مدل ارائه می‌دهد؛ اما چون زمان محاسبه در آن زیاد و وقت‌گیر است، امروزه خیلی شایان توجه قرار نمی‌گیرد [31].
ازطریق روش تخمین دو نقطه‌ای به محاسبه دو مقدار قطعی در هر متغیر نامعیناستفاده می‌شود که در هر دو سمت مقدار میانگین آن قرار گرفته است. بنابراین، پخش توان بهینه در هر متغیر نامعین دو مرتبه محاسبه می‌شود؛ یکبار در مقدار زیر میانگین و بار دیگر در مقدار بالاتر از میانگین؛ در حالی که سایر متغیرها در مقدار میانگین خود حفظ می‌شوند.
ممکن است این دو نقطه در حوزۀ میانگین متغیر ارائه‌شده به‌صورت متقارن یا نامتقارن باشند. سپس هر مجموعه از نقاط انتخاب‌شده آزمایش تابع غیرخطی می‌شوند. با توجه به اینکه مقدار متغیرهای ورودی نامعین در نمونه‌گیری‌های مختلف متفاوت است، راه‌حل بهینه در هریک از این مسائل تغییر می‌یابد [32]. درخور ذکر است در این روش، تعداد نمونه‌های انتخابی با افزایش تعداد متغیرهای نامعین افزایش می‌یابد. در سیستم‌های مقیاس بزرگ که دارای تعداد متغیر نامعین بیشتری‌اند، روش تخمین دو نقطه‌ای به‌خوبی نمونه‌های مقیاس کوچک اجرا نمی‌شود.
فرمول‌بندی اصلی روش تخمین دو نقطه‌ای با وضعیت نامتقارن دو نقطه نمونه‌گیری به شرح زیر است [33]:
مرحله اول) تعیین تعداد متغیرهای نامعین
مرحله دوم) تنظیم
مرحله سوم) تنظیم
مرحله چهارم) تعیین
، و با معادلات (22) تا (25):
|
(22) |
|
|
(23) |
|
|
(24) |
|
|
(25) |
|
مرحله پنجم) تعیین دو تجمیع و
|
(26) |
|
|
(27) |
|
مرحله ششم) اجرای پخش توان بهینۀ قطعی برای هر دو تجمیع با استفاده از که است.
مرحله هفتم) محاسبۀ مقادیر
|
(28) |
|
|
(29) |
|
مرحله هشتم) محاسبۀ میانگین و انحراف استاندارد متغیرهای خروجی
|
(30) |
|
|
(31) |
|
روش تخمین سه نقطه‌ای برای تجزیه و تحلیل متغیرهای تصادفی و به دست آوردن داده‌های منحنی ولتاژ و پخش توان استفاده می‌شوند. هنگامی که تابعی از متغیرهای تصادفی باشد، چنانچه متغیر تصادفی ام انتخاب شده باشد، مقدار نمونه‌گیری باید دارای سه نقطه نمونه‌گیری باشد تا جایگزین متغیر تصادفی اصلی شود. روش تخمین سه نقطه‌ای از نقطه نمونه‌گیری استفاده می‌کند؛ اما نقاط نمونه‌گیری دارای داده‌های یکسانی‌اند. درواقع، تنها نقاط نمونه‌گیری نیاز به محاسبه‌شدن دارند؛ بنابراین، الگوریتم روش تخمین سه نقطه‌ای به‌صورت زیر توصیف می‌شود [34]:
مرحلۀ اول) نقاط نمونه‌گیری در توزیع نرمال استاندارد انتخاب می‌شود. هر جز در توزیع نرمال استاندارد بیان شده است. میانگین برابر است با صفر، انحراف استاندارد برابر است با 1، ضریب چولگی برابر است با صفر و کشیدگی برابر 3 است. اگر نقاط نمونه‌گیری مؤلفه باشد، معادله آن به شکل معادله (32) فرمول‌بندی می‌شود:
|
(32) |
|
در اینجا، برابر با 1 و 1، و ضرایب چولگی و کشیدگی و تعداد نقاط نمونه‌گیری را نشان می‌دهند. زمانی که برابر 3 است، میانگین جایگزین تمام مؤلفه های اصلی می‌شود.
مرحله دوم) در این مرحله، گشتاور و نتایج پخش توان بهینۀ احتمالاتی به دست می‌آید. هنگامی که برابر با 1 و 2 است، وزن‌های مقدار تابع در تمام نقاط نمونه‌گیری در معادله (33) به دست می‌آید. زمانی که برابر 3 است، وزن مقدار تابع با معادله (34) به دست می‌آید؛ ازاین‌رو، گشتاور خروجی متغیر با معادله (35) به دست می‌آید.
|
(33) |
|
|
(34) |
|
|
(35) |
|
مرحله سوم) مقدار پیش‌بینی و انحراف استاندارد متغیر خروجی را با معادلات (36) و (37) به دست می‌آید:
|
(36) |
|
|
(37) |
|
روش تبدیل بی بو برای غلبه بر نقایص مرتبط با روش‌های احتمالاتی سنتی معرفی شده است، به‌ویژه آنهایی که از فرآیند خطی‌سازی استفاده می‌کنند. روش‌های تحلیلی براساس برخی مفروض‌های ریاضی و الگوریتم‌های پیچیده ایجاد می‌شوند. این روش، رویکرد قدرتمندی است که در محاسبۀ مسائل تصادفی با و بدون متغیرهای عدم‌قطعیت استفاده می‌شود [10].
قلب روش تبدیل بی بو این است که چگونه می‌توان نمونه‌های مناسب از متغیرهای ورودی را تولید کرد تا بتواند اطلاعات مناسب دربارۀ تابع چگالی احتمال متغیرهای ورودی را حفظ کند. فرض برداری از متغیرهای تصادفی بعدی با میانگین و کوواریانس است. متغیر تصادفی دیگر وابسته به است که با معادله (19) محاسبه می‌شود، می‌تواند مجموعه‌ای از توابع غیرخطی باشد.
با استفاده از روش تبدیل بی بو میانگین و کوواریانس متغیرهای خروجی و با مراحل زیر به دست می‌آید:
مرحله 1) با استفاده از
|
(38) |
|
|
(39) |
|
|
(40) |
|
(1+n2) نمونه یا نقاط نمونه به دست می‌آید.
مرحله 2) با استفاده از
|
(41) |
|
|
(42) |
|
|
(43) |
|
وزن‌های مرتبط با هر محاسبه می‌شوند. باید تأکید کرد وزن‌های مرتبط باید شرایط معادله (44) را داشته باشند:
|
(44) |
|
در معادلات (40) و (41)، مقدار ردیف یا ستون امین ریشۀ دوم ماتریس است. ریشۀ دوم ماتریس به این معنا است که ماتریس و باشد. در اینجا وزن اختصاص‌یافته به نقطه است که به نام نقطۀ صفر شناخته می‌شود. این نقطه موقعیت سایر نقاط اطراف مقدار میانگین را کنترل می‌کند.
مرحله 3) در این مرحله هرکدام از نقاط نمونه به یک تابع غیرخطی می‌پردازند تا مجموعه‌ای از نقاط نمونه تبدیل‌شده با معادله (45) به دست آیند.
|
(45) |
|
باید تأکید کرد در روش تبدیل بی بو تابع غیرخطی، جعبه سیاه در نظر گرفته می‌شود؛ ازاین‌رو، ساده‌سازی یا خطی‌سازی ضروری نیست. شکل (1) این روش را نشان می‌دهد.
شکل (1): نتیجۀ حاصل از روش تبدیل بی بو مربوط به مرحله 3
مرحله 4) در این مرحله با معادلات (46) و (47) به‌ترتیب میانگین و کوواریانس متغیر خروجی محاسبه می‌شوند.
|
(46) |
|
|
(47) |
|
روش تبدیل بی بو، دو ویژگی استثنایی دارد؛ نخست، نقاط نمونه به‌طور تصادفی انتخاب نمی‌شوند؛ بلکه به طریقی انتخاب می‌شوند که یک معیار پیش‌بینی‌شده میانگین و کوواریانس داشته باشند. دوم، وزن‌های مربوط به نقاط انتخاب‌شده در محدوده [1و0] نباید باشند؛ بلکه می‌توانند مقادیر مثبت یا منفی داشته باشند. این نقاط باید شرط معادله (44) را برآورده کنند که وزن‌ها باید برابر با واحد باشد.
در این بخش، مراحل اصلی روش متداول نقطۀ داخلی ارائه شدهاند که این روش به‌طور دقیق در مرجع [35] توضیح داده شده است. فرمول‌بندی بهینه‌سازی غیرخطی محدودشده در معادله (48) ارائه شده است.
|
(48) |
|
در این فرمول‌بندی، محدودیت‌های عملکردی (محدودیت‌های نابرابری) به محدودیت‌های برابری تبدیل می‌شوند و با استفاده از متغیرهای باقیمانده در وارد می‌شوند؛ همان‌طور که در [36] دیده می‌شود. سپس بردار شامل هر دو متغیرهای اصلی و باقی مانده است. متغیرهای اسلک و به محدودیت‌های برابری تبدیل می‌شوند.
محدودیت‌های معادله (48) با اضافه‌کردن یک تابع مانع لگاریتمی به تابع هدف حذف می‌شوند که این تغییرات در معادله (49) نشان داده شدهاند. برای نتیجه، متغیرهای و باید از صفر بیشتر باشند و متغیرهای هرگز نمی‌توانند روی مقادیر مرزها فرض شوند.
|
(49) |
|
در ابتدا پارامتر مانع باید بیشتر از صفر فرض شود، مقدار در پایان روند تکرار باید نزدیک به صفر باشد. پس تابع لاگرانژ تعریف‌شده در معادله (50) نشان داده شده است.
|
(50) |
|
شرایط بهینگی روش کروش کان تاکر[15] مرتبۀ اول برای مسئلۀ بهینه‌سازی با نتایج روش نیوتون رافسون در معادله (51) نشان داده شده است.
|
(51) |
|
که در این معادله
|
(52) |
|
بسط ریاضی برای دستیابی به معادله (51) در ضمیمۀ مقاله ارائه شده است.
معادلات (52) نشان‌دهندۀ پارامترهای مانع و متغیرهای اسلک به عناصر قطری ماترس و همچنین بردار گرادیان است. مشاهده می‌شود مقادیر صفر را نمی‌توان برای متغیرهای اسلک در نظر گرفت.
از حل معادله (51) و به دست می‌آید. همچنین ، ، و به‌ترتیب از معادلات (70) و (71) به دست می‌آیند که در ضمیمۀ مقاله نشان داده شدهاند. طول گام و با معادلات (53) محاسبه می‌شوند. از این معادلات برای حفظ مثبت‌بودن و و همچنین سیگنال مناسب و استفاده می‌شود.
|
(53) |
|
پارامتر با مقدار 9995/0 معمولاً به معادلات (50) برای حفظ مقادیر و اضافه می‌شود. سپس متغیرهای مسئلۀ بهینه‌سازی با معادلات (50) بروزرسانی می‌شوند.
|
(54) |
|
پارامتر مانع با توجه به شکاف دوگانگی با معادلۀ (55) بروزرسانی می‌شود [37]. پارامتر در معادله (55) برای کنترل و بهبود فرآیند همگرایی معرفی شده است [15].
|
(55) |
|
از معادله (55) چنین استنباط می‌شود که حاصل ضرب‌های یا ، هنگامی که محدودیت مرتبط - با منتهی‌شدن یا به یک مقدار بزرگ وقتی یا نمی‌توانند صفر شوند - فعال می‌شود، بسیار بزرگ خواهند شد؛ درنتیجه، حاصل ضرب‌های یا بزرگ خواهند شد، در حالی که باید صفر باشند. این ویژگی باعث شده است مسئلۀ پخش توان بهینه از لحاظ عددی، ناپایدار یا ناهمگرا شود.
از معادله (49) دیده می‌شود متغیرهای نمی‌توانند حد مجاز خود را در نظر بگیرند؛ زیرا متغیرهای اسلک نمی‌توانند دقیقاً برابر صفر باشند که این قضیه بر کیفیت راه‌حل تأثیر می‌گذارد. با توجه به مقدمۀ مقاله، بسیاری از آثار پیشنهادی در ادبیات برای غلبه بر مشکلات یادشده بیان شدهاند.
برای حل مشکل مسئلۀ ابتدا، پارامترهای جدید در معادله (49) اضافه می‌شود. مقدار همیشه مثبت است. در این فرمولاسیون جدید، فرآیند همگرایی پارامتر به جای و کنترل می‌شود. واردکردن به تابع مانع لگاریتمی موجب می‌شود متغیرهای اسلک، اکنون به مقدار صفر برسند. در طی فرآیند تکرار نیوتون رافسون مقدار به‌صورت مستقل از مقدار پارامتر است که در معادلۀ (56) تعریف می‌شود.
|
(56) |
|
در معادله (56)، تکرار شمارش، ضریب کاهش است. پارامتر تعریف‌شده برای سرعت‌بخشیدن به سرعت همگرایی است.
درنهایت پارامتر اصلاح می‌شود. معادله (69) در ضمیمۀ مقاله برای کاهش پارامتر مانع استفاده می‌شود و پارامتر برای ها به‌صورت معادله (57) بازنویسی می‌شود.
|
(57) |
|
جملۀ اول جمع در سمت راست معادله (52) شبیه GAP شناخته‌شده در منابع [15] و [36] است. علاوه بر این، با توجه به مثبت‌بودن و منفی‌بودن اصلاح جدیدی ظاهرشده در انتهای معادله (58) همیشه مثبت است و تمایل دارد یک مقدار بسیار کوچک باشد؛ زیرا در همگرایی فرض می‌شود یک مقدار نزدیک به صفر است که در معیارهای اولیه و همگرایی دیده می‌شود. مقدار جدید با معادله (58) محاسبه می‌شود.
|
(58) |
|
پارامتر معمولاً در محاسبات برای تسریع روند همگرایی معرفی می‌شود. زمانی که متغیر اسلک صفر باشد، ضریب لاگرانژ یک مقدار بزرگ فرض می‌شود؛ زیرا متغیر بهینه‌سازی در حد پایین قرار دارد؛ درنتیجه ضرایب و تبدیل به صفر می‌شوند که این عمل موجب می‌شود مسئلۀ پخش توان بهینه از لحاظ عددی پایدار و ساده‌تر شود تا به توان همگرایی دست یابد.
مدل احتمالاتی توان بادی
برای مدلسازی سرعت باد از تابع توزیع ویبال استفاده می‌شود. این مدلسازی به‌صورت زیر انجام می‌شود [38]:
مرحلۀ اول) سرعت باد با تابع چگالی احتمال مناسب همچون ویبال مدلساز شده است که این تابع سرعت باد و احتمال وقوع را در هر سرعت باد تعیین می‌کند. تابع توزیع ویبال با معادله (59) مدلسازی می‌شود.
|
(59) |
|
در اینجا، سرعت باد بر حسب متر بر ثانیه، پارامتر شکل در توزیع ویبال و پارامتر مقیاس در توزیع ویبال بر حسب متر بر ثانیه است.
مرحلۀ دوم) برای بررسی عدم‌قطعیت، مسئلۀ چندین مرتبه ارزیابی‌شده تا مهم‌ترین عامل یا شرایط احتمالی پوشش داده می‌شود. به‌منظور مدلسازی عدم‌قطعیت توان توربین بادی، نمونه‌های سرعت باد تولیدشده در هر ارزیابی به شیوۀ مناسب بررسی می‌شوند.
مرحلۀ سوم) نمونه‌های سرعت باد تولیدشده با استفاده از منحنی توان ـ سرعت به توان خروجی توربین تبدیل می‌شوند. این تبدیل با معادله (60) توضیح داد شده است [39].
|
(60) |
|
در اینجا، توان اکتیو خروجی توربین بادی بر حسب مگاوات، توان اکتیو اندازه‌گیری‌شده توربین بادی برحسب مگاوات، سرعت شروع تولید توان اکیتو بر حسب متر بر ثانیه، سرعت توقف تولید توان اکتیو بر حسب متر بر ثانیه، سرعت اندازه‌گیری‌شده توربین بادی بر حسب متر بر ثانیه و سرعت باد بر حسب متر بر ثانیه است. تغییرات تابع چگالی احتمال ویبال در شکل (2) نشان داده شدهاند.
شکل (2): تابع چگالی احتمال ویبال
مدل احتمالاتی توان فتوولتائیک
تولید توان فتوولتائیک به میزان زیادی از شرایط طبیعی تأثیر می‌گیرد؛ بنابراین، توان خروجی با توجه به شدت نور تغییر می‌کند. پژوهش‌ها نشان می‌دهند شدت نور در مدت زمان کوتاه با تابع توزیع بتا توصیف می‌شود [40]. تابع چگالی احتمال با معادله (61) توصیف می‌شود:
|
(61) |
|
در معادله (56)، و پارامترهای شکل توزیع بتا، تابع گاما، شدت نور واقعی در یک دوره و حداکثر شدت نور در خلال دورۀ زمانی‌اند.
پارامترهای و با مقدار میانگین و انحراف استاندارد شدت نور با معادلات (62) و (63) بیان می‌شوند:
|
(62) |
|
|
(63) |
|
پس از به دست آوردن تابع چگالی احتمال نور، توان خروجی و حداکثر توان خروجی محاسبه می‌شوند. و با معادلات (64) و (65) به دست می‌آیند:
|
(64) |
|
|
(65) |
|
در اینجا، کل مساحت آرایۀ فتوولتائیک و بازده تبدیل باتری است.
بنابراین، تابع چگالی احتمال توان خروجی با معادله (66) به دست می‌آید:
|
(66) |
|
درنتیجه، تابع چگالی احتمال توان خروجی با توزیع بتا توصیف می‌شود.
مقادیر تابش خورشید و دمای هوا دو پارامتر مهم در تولید توان الکتریکی نیروگاه خورشیدی‌اند. این پارامترها بر حسب زمان متغیرند. همان‌طور که اشاره شد در این مقاله تابش خورشید با استقاده از تابع توزیع بتا مدل شده است. شرایط واقعی بهره‌برداری نیروگاه خورشیدی با وضعیت استاندارد متفاوت است؛ بنابراین میزان توان تولیدی سلول فتوولتائیک در شرایط واقعی به‌صورت معادله (67) محاسبه می‌شود:
|
(67) |
|
که اندازۀ توان سلول فتوولتائیک، میزان تابش خورشید، مقدار تابش در شرایط آزمون استاندارد، ضریب حداکثر دما، دمای سلول بر حسب درجۀ سانتیگراد و دمای مرجع بر حسب درجۀ سانتیگراد است. تغییرات تابع چگالی احتمال بتا در شکل (3) نشان داده شدهاند.
شکل (3): تابع چگالی احتمال بتا
مدل احتمالاتی بار
در مدل قطعی تقاضای بار در هر شین، ثابت در نظر گرفته می‌شود. زمان و آب و هوا دو عامل برای مؤلفۀ قطعی از تغییرات است، در حالی که اجزای تصادفی متغیرهای تصادفی مستقل‌اند. الگوهای رفتاری مصرف‌کنندگان انرژی به تقاضا متغیر در هر شین منجر می‌شود. این تغییرات با تجزیه و تحلیل آماری محاسبه می‌شوند؛ درنتیجه، تقاضا به‌طور مداوم با درجۀ بالایی از عدم‌قطعیت متفاوت است.
بار با تابع چگالی احتمال تعریف می‌شود. توابع چگالی مختلفی برای مدل‌کردن این متغیرهای تصادفی وجود دارد؛ برای مثال، توابع چگالی احتمال ویبال، بتا و نرمال [41]. در این مقاله تقاضای بار با توزیع نرمال با مقدار متوسط و انحراف استاندارد در نظر گرفته می‌شود. تابع توزیع نرمال با معادله (68) محاسبه می‌شود. فرم تغییرات این تابع در شکل (4) نشان داده شده است.
|
(68) |
|
شکل (4): تابع چگالی احتمال نرمال
2-12- ضریب همبستگی
در روش‌های پیشنهادی، هزینۀ احتمالی با راه‌حل‌های احتمالاتی محاسبه می‌شود. روش‌های احتمالاتی برای لحاظ عدم‌قطعیت مربوط به مقدار بار مورد نیاز و انرژی تولیدی از منابع انرژی باد و خورشید در نظر گرفته می‌شوند. در هر سیستم قدرت مقدار بارها، سرعت باد و تابش خورشید به وضوح مشخص نیستند؛ زیرا هر کدام از این عدم‌قطعیت‌ها یک متغیر تصادفی‌اند. به‌طور کلی استفاده از یک تابع چگالی، احتمال روش مناسب برای مدلسازی رفتار تصادفی است. در مطالعات انجام‌شده در این مقاله، فرمول‌بندی تابع چگالی احتمال باد، خورشید و بار به‌ترتیب با معادلات (59)، (61) و (63) انجام می‌شود؛ با این حال، زمانی که همبستگی بین عدم‌قطعیت‌ها وجود دارد، تغییرات آنها وابسته به یکدیگرند که این وابستگی نیز باید در محاسبۀ پخش توان بهینۀ احتمالاتی در نظر گرفته شود.
اگر واریانس (انحراف استاندارد) بین متغیرهای عدم‌قطعیت و باشد، ماتریس واریانس (انحراف استاندارد) متغیر همبسته به‌صورت معادله (69) نشان داده می‌شود [42]:
|
(69) |
|
که در آن ضریب همبستگی متغیرهای و با معادله (70) محاسبه می‌شود.
|
(70) |
|
ضریب همبستگی مقداری در بازه 1- تا 1+ است. هنگامی که دو متغیر مستقل باشند، ضریب همبستگی آنها صفر است. اگر دو متغیر به همدیگر وابسته باشند، دارای ضریب همبستگی مثبت‌اند. وابستگی درخلاف جهت که نشان‌دهندۀ افزایش یک متغیر هنگام کاهش دیگری است، وابستگی معکوس، نامیده و با ضریب همبستگی منفی نشان داده می‌شود.
برای تحلیل پخش توان بهینۀ احتمالاتی، شبکۀ استاندارد 118 شینه IEEE، دارای 118 شین، 54 ژنراتور و 187 خط انتقال، شبکه مورد مطالعه انتخاب شده است [42]. از تحلیل نتایج این شبکه برای سنجش دقت روش‌های شبیه‌سازی مونت کارلو و روش‌های تحلیلی استفاده می‌شود. شبیه‌سازی در یک رایانه شخصی با پردازنده 8/1 گیگا هرتز و حافظه 8 گیگا بایت اجرا شده است.
اطلاعات مربوط به این شبکه در [45-43] بیان شده است. همچنین اطلاعات مربوط به توربین بادی و نیروگاه خورشیدی در پیوست ارائه شدهاند. با توجه به اینکه در مرجع‌های مختلف شین‌های متفاوتی برای استقرار نیروگاههای خورشیدی و توربین بادی پیشنهاد شدهاند، در این مقاله محل استقرار نیروگاههای فوق به‌ترتیب در شین‌های 7 و 8 لحاظ شده است. بدیهی است با تغییر محل استقرار، امکان بروز تغییراتی در نتایج بیان‌شده، دور از انتظار نخواهد بود.
نرم‌افزارهای استفاده‌شده برای سنجش روش‌های عددی و تحلیلی، MATLAB و DigSilent است. روش‌های شبیه‌سازی مونت کارلو، تخمین نقطه‌ای و تبدیل بی بو با نرم‌افزار MATLAB و روش نقطۀ داخلی با نرم‌افزار DigSilent اجرا شدهاند. دیاگرام تک‌خطی مورد مطالعه در شکل (5) نشان داده شد.
شکل (5): نمودار تک‌خطی شبکه 118 شینه IEEE مورد مطالعه
با مقایسه همزمان پنج روش عددی و تحلیلی نشان داده شده در جدول (1)، در روش‌های شبیه‌سازی مونت کارلو، تخمین دو و سه نقطه‌ای و روش تبدیل بی بو بازه تغییرات نزدیک‌به‌هم است؛ اما با توجه به نتایج، در روش نقطۀ داخلی به دلیل وجود پارامترهای مانع و جریمه که بازۀ تولیدات را محدود می‌کند، نتایج مقادیر متفاوت‌تری نسبت به سایر روش‌ها دارد؛ اما همچنان تمامی نتایج طبق قیدهای تولید توان در بازه تعیین شده است و و هیچ‌کدام از توان‌های تولیدشده از محدودۀ حداقل یا حداکثر خود خارج نشده‌اند. این تغییرات در نمودار شکل (6) نشان داده شدهاند.
جدول (1): مقایسۀ همزمان نتایج میانگین تغییرات توان اکیتو تولیدی (بر حسب مگاوات) در چهار روش محاسباتی
|
میانگین |
شماره ژنراتور |
||||
|
IPM |
UT |
3PEM |
2PEM |
MCS |
|
|
97/49 |
0 |
0 |
0 |
0 |
1 |
|
100 |
0 |
0 |
0 |
0 |
2 |
|
5 |
0 |
0 |
0 |
0 |
3 |
|
03/5 |
100 |
100 |
100 |
100 |
4 |
|
15 |
33/286 |
17/279 |
50/284 |
20/282 |
5 |
|
99/89 |
0 |
0 |
0 |
0 |
6 |
|
99/23 |
15/81 |
76/76 |
98/79 |
58/78 |
7 |
|
01/25 |
0 |
0 |
0 |
0 |
8 |
|
75/26 |
0 |
0 |
0 |
0 |
9 |
|
99/29 |
100 |
100 |
100 |
100 |
10 |
|
98/29 |
320 |
320 |
320 |
320 |
11 |
|
30 |
0 |
0 |
0 |
0 |
12 |
|
03/8 |
0 |
0 |
0 |
0 |
13 |
|
20 |
18/73 |
90/69 |
34/72 |
28/71 |
14 |
|
100 |
0 |
0 |
0 |
0 |
15 |
|
99/79 |
75/70 |
52/67 |
92/69 |
88/68 |
16 |
|
99/99 |
0 |
0 |
0 |
0 |
17 |
|
03/8 |
0 |
0 |
0 |
0 |
18 |
|
01/8 |
06/72 |
39/70 |
63/71 |
09/71 |
19 |
|
98/99 |
119 |
119 |
119 |
119 |
20 |
|
39/242 |
304 |
304 |
304 |
304 |
21 |
|
97/99 |
41/53 |
25/52 |
11/53 |
73/52 |
22 |
|
99/99 |
54/52 |
41/51 |
24/52 |
88/51 |
23 |
|
03/25 |
100 |
100 |
100 |
100 |
24 |
|
02/5 |
255 |
255 |
255 |
255 |
25 |
|
99/199 |
36/39 |
92/37 |
98/38 |
52/38 |
26 |
|
99/99 |
100 |
100 |
100 |
100 |
27 |
|
99/41 |
490 |
40/485 |
63/490 |
19/490 |
28 |
|
99/41 |
29/260 |
52/258 |
52/259 |
25/259 |
29 |
|
99/99 |
93/37 |
51/37 |
82/37 |
68/37 |
30 |
|
33 |
0 |
0 |
0 |
0 |
31 |
|
03/5 |
0 |
0 |
0 |
0 |
32 |
|
5 |
0 |
0 |
0 |
0 |
33 |
|
03/5 |
96/66 |
59/65 |
60/66 |
16/66 |
34 |
|
99/99 |
03/68 |
95/66 |
75/67 |
40/67 |
35 |
|
25 |
100 |
100 |
100 |
100 |
36 |
|
50 |
91/52 |
67/51 |
58/52 |
18/52 |
37 |
|
100 |
0 |
0 |
0 |
0 |
38 |
|
01/10 |
104 |
104 |
104 |
104 |
39 |
|
99/24 |
02/258 |
29/257 |
83/257 |
54/257 |
40 |
|
99/41 |
0 |
0 |
0 |
0 |
41 |
|
98/19 |
0 |
0 |
0 |
0 |
42 |
|
40 |
100 |
100 |
100 |
100 |
43 |
|
98/99 |
100 |
100 |
100 |
100 |
44 |
|
40 |
48/262 |
74/261 |
28/262 |
05/262 |
45 |
|
80 |
0 |
0 |
0 |
0 |
46 |
|
72/216 |
65/45 |
25/45 |
54/45 |
41/45 |
47 |
|
98/129 |
67/50 |
28/50 |
57/50 |
44/50 |
48 |
|
100 |
0 |
0 |
0 |
0 |
49 |
|
17 |
0 |
0 |
0 |
0 |
50 |
|
40 |
87/39 |
63/39 |
81/39 |
73/39 |
51 |
|
40 |
66/55 |
44/55 |
61/55 |
53/55 |
52 |
|
94/79 |
12/64 |
12/60 |
09/63 |
81/61 |
53 |
|
20 |
0 |
0 |
0 |
0 |
54 |
بر طبق روش شبیه‌سازی مونت کارلو که براساس تکرار عمل می‌کند، طی روند زمان طولانی‌تر به جواب بهینه دست می‌یابد؛ اما جواب‌های به‌دست‌آمده به واقعیت نزدیک‌تر است؛ به همین دلیل این روش، مبنای محاسبات در نظر گرفته می‌شود.
در روش تخمین دو نقطه‌ای به دلیل اینکه این روش دارای نقاط ورودی کمتری نسبت به روش شبیه‌سازی مونت کارلو است، طبق معادلات (17) تا (22) گشتاورهای مربوط به این نقاط محاسبه می‌شود و درنهایت میانگین و انحراف استاندارد توان ژنراتورها محاسبه می‌شوند. این روش برای رسیدن به جواب بهینه سریع‌تر از روش شبیه‌سازی مونت کارلو است.
برای اصلاح روش تخمین دو نقطه‌ای، در روش تخمین سه نقطه‌ای علاوه بر محاسبۀ وزن برای هر متغیر از معادله (30) استفاده می‌شود که مربوط به ضرایب چولگی و کشیدگی است. این معادلات سبب افزایش دقت این روش می‌شود. ضرایب چولگی نشان‌دهندۀ نامتقارن‌بودن توزیع احتمالاتی، برابر با گشتاور سوم نرمال فرض شده است و در صورتی که داده‌ها نسبت به میانگین متقارن باشند، برابر صفر است و مقدار آن برابر یک توزیع نامتقارن با کشیدگی به سمت بالاتر مثبت و برای توزیع نامتقارن با کشیدگی به سمت مقادیر کوچک‌تر، منفی است.
همچنین ضریب کشیدگی، توصیف‌کنندۀ میزان قله‌ای‌بودن و یا مسطح‌بودن توزیع احتمالی و برابر با گشتاور چهارم نرمال است. مقدار این ضریب برای توابع چگالی احتمال تیز و دنبالۀ پهن، بیشتر خواهد بود و معیاری از تیزی منحنی در نقطۀ حداکثر است. بنا بر دلایل یادشده به‌کارگیری روش تخمین سه نقطه‌ای توصیه می‌شود.
با توجه به نتایج به‌دست‌آمده در روش تبدیل بی بو، براساس به‌کارگیری معادلات (33) تا (42)، این روش برخلاف روش تخمین نقطه‌ای، نقاط نمونه را تصادفی انتخاب نمی‌کند؛ بلکه طوری این نقاط انتخاب می‌شوند که دارای مقادیر میانگین و کوواریانس باشند. وزن‌های مربوط به نقاط باید طبق معادله (39) مثبت و منفی باشند. نتایج نشان می‌دهند دقت این روش در حد روش شبیه‌سازی مونت کارلو است.
شکل (6): نمودار مقایسۀ نتایج میانگین تغییرات توان اکتیو تولیدی (بر حسب مگاوات) در روش‌های مورد مطالعه
اگر شبکۀ استاندارد از دیدگاه پروفیل ولتاژ تولیدی ژنراتورها بررسی شود، مطابق شکل (7) مشاهده می‌شود که کمترین تغییرات ولتاژ مربوط به روش نقطۀ داخلی است و بازۀ تغییرات این روش بین 1 تا 01/1 پریونیت است، پس دارای ثبات بیشتری است. بیشترین بازۀ تغییرات مربوط به روش شبیه‌سازی مونت کارلو است. کمترین میزان آن 943/0 پریونیت مربوط به ژنراتور شماره 35 است و بیشترین مقدار آن 05/1 پریونیت مربوط به ژنراتورهای 11 و 29 است. مقادیر تغییرات روش‌های تخمین نقطه‌ای و روش تبدیل بی بو مشابه همدیگرند که این روش‌ها دارای بازۀ تغییرات حداقل 0239/1 پریونیت و حداکثر 0599/1 هستند.
شکل (7): نمودار تغییرات پروفیل ولتاژ تولیدی (بر حسب کیلو ولت) در روش‌های مورد مطالعه
با توجه به نتایج ارائه‌شده در شکل (8)، چون در روش شبیه‌سازی مونت کارلو در شروع محاسبات، از نقاط بیشتری استفاده می‌شود، این روش بیشترین زمان را در رسیدن به جواب بهینه دارد. روش نقطۀ داخلی به دلیل اینکه از ضرایب جریمه‌ای در حل استفاده می‌کند، دارای کمترین زمان اجرا است، پس سایر روش‌ها در بین روش نقطۀ داخلی و شبیه‌سازی مونت کارلو قرار دارند. اگر روش نقطۀ داخلی بهترین جواب ازنظر زمان اجرا در نظر گرفته شود، روش سه نقطه‌ای دارای 22/28 درصد، روش تخمین دو نقطه‌ای دارای 78/77 درصد، روش تبدیل بی بو 23/222 و روش شبیه‌سازی مونت کارلو دارای 12/2001 درصد افزایش زمان اجرا نسبت به روش نقطۀ داخلی است.
شکل (8): نمودار مقایسۀ نتایج زمان اجرا (بر حسب ثانیه) در روش‌های مورد مطالعه
با توجه به تحلیل تابع هزینۀ ارائه‌شده در شکل (9)، هزینۀ پخش توان بهینه بر حسب دلار بر ساعت در روش تخمین سه نقطه‌ای دارای کمترین مقدار نسبت به سایر روش‌ها است. سایر روش‌ها همانند شبیه‌سازی مونت کارلو، روش تخمین دو نقطه‌ای، روش تبدیل بی بو و روش نقطۀ داخلی به‌ترتیب دارای 49/0، 83/0، 11/1 و 50/7 درصد افزایش نسبت به روش تخمین سه نقطه‌ای است.
شکل (10) مقادیر میانگین هزینۀ پخش توان بهینۀ حاصل از پخش توان احتمالاتی با روش‌های شبیه‌سازی مونت کارلو و نقطۀ داخلی براساس تعداد دفعات تکرار مختلف را نشان می‌دهد. نقطۀ شروع در شکل (10) مقدار هزینۀ حاصل از مقدار قطعی پخش توان را نشان می‌دهد. با توجه به شکل، مقدار آن بیشتر از مقدار متناظر به هنگام منظورکردن اثر نیروگاههای بادی و خورشیدی است. این موضوع نشان‌دهندۀ کاهش هزینۀ بهره‌برداری شبکه به هنگام نفوذ انرژی‌های تجدیدپذیر است.
شکل (9): نمودار مقایسۀ نتایج تابع هزینۀ پخش توان بهینه (بر حسب دلار بر ساعت) در روش‌های مورد مطالعه
با مقایسه دو روش شبیه‌سازی مونت کارلو و نقطۀ داخلی، روش نقطۀ داخلی در نخستین مرحلۀ تکرار
(50 تکرار) به حالت همگرایی و جواب بهینه دست می‌یابد و از این تکرار به بعد جواب ثابت می‌ماند و برای این روش همان جواب بهینه هزینه خواهد بود. در روش شبیه‌سازی، مونت کارلو بعد از چندین مرحلۀ تکرار (80000 تکرار) به حالت همگرایی دست می‌یابد.
از دیدگاه دیگر، روش نقطۀ داخلی سریع‌ترین حالت همگرایی را دارد، اما جواب جواب بهینۀ آن بر حسب دلار بر ساعت بیشتر از روش دیگر است؛ در حالی که روش شبیه‌سازی مونت کارلو دیرتر به حالت همگرایی می‌رسد؛ درعوض جواب آن دارای مقدار بهتر و کمتری نسبت به روش نقطۀ داخلی است.
شکل (10): مقایسه روند همگرایی مقادیر میانگین تابع هزینه پخش بار بهینه در روش شبیه‌سازی مونت کارلو و روش نقطۀ داخلی
در صورتی که مقادیر همبستگی توان تولیدی منابع تجدیدپذیر واقع در شین‌های مختلف سیستم قدرت تغییر پیدا کند، اثرات این تغییر بر متغیرهای پخش توان بهینه کاملاً مشهود است و نتایج به‌صورت جدول (2) نشان داده می‌شوند:
جدول (2): مقایسۀ تأثیر تغییر ضریب همبستگی منابع تجدیدپذیر بر تغییرات توان اکیتو تولیدی (بر حسب مگاوات) در چهار روش محاسباتی با لحاظ ضریب همبستگی برابر 5/0 برای انرژی باد و خورشید و 7/0 برای بار
|
میانگین توان تولیدشده بر اثر تغییرات ضریب همبستگی |
شماره ژنراتور |
|||||||||
|
IPM |
UT |
3PEM |
2PEM |
MCS |
||||||
|
بعد |
قبل |
بعد |
قبل |
بعد |
قبل |
بعد |
قبل |
بعد |
قبل |
|
|
80/48 |
97/49 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
1 |
|
75/97 |
100 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
2 |
|
75/5 |
5 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
3 |
|
03/5 |
03/5 |
100 |
100 |
05/109 |
100 |
9/105 |
100 |
100 |
100 |
4 |
|
45/13 |
15 |
44/287 |
33/286 |
6/308 |
17/279 |
54/299 |
50/284 |
89/281 |
20/282 |
5 |
|
90 |
99/89 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
6 |
|
85/23 |
99/23 |
92/81 |
15/81 |
26/86 |
76/76 |
68/83 |
98/79 |
42/78 |
58/78 |
7 |
|
98/24 |
01/25 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
8 |
|
70/26 |
75/26 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
9 |
|
65/28 |
99/29 |
100 |
100 |
05/109 |
100 |
9/105 |
100 |
100 |
100 |
10 |
|
25/30 |
98/29 |
320 |
320 |
96/348 |
320 |
89/338 |
320 |
320 |
320 |
11 |
|
10/30 |
30 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
12 |
|
9 |
03/8 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
13 |
|
04/20 |
20 |
68/73 |
18/73 |
13/78 |
90/69 |
8/75 |
34/72 |
14/71 |
28/71 |
14 |
|
102 |
100 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
15 |
|
13/80 |
99/79 |
25/71 |
75/70 |
51/75 |
52/67 |
26/73 |
92/69 |
74/68 |
88/68 |
16 |
|
46/98 |
99/99 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
17 |
|
9 |
03/8 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
18 |
|
87/8 |
01/8 |
34/72 |
06/72 |
74/77 |
39/70 |
46/75 |
63/71 |
02/71 |
09/71 |
19 |
|
04/101 |
98/99 |
119 |
119 |
77/129 |
119 |
02/126 |
119 |
119 |
119 |
20 |
|
36/245 |
39/242 |
304 |
304 |
51/331 |
304 |
94/321 |
304 |
304 |
304 |
21 |
|
47/98 |
97/99 |
63/53 |
41/53 |
66/57 |
25/52 |
97/55 |
11/53 |
68/52 |
73/52 |
22 |
|
01/100 |
99/99 |
75/52 |
54/52 |
73/56 |
41/51 |
07/55 |
24/52 |
83/51 |
88/51 |
23 |
|
05/25 |
03/25 |
100 |
100 |
05/109 |
100 |
9/105 |
100 |
100 |
100 |
24 |
|
5 |
02/5 |
255 |
255 |
07/278 |
255 |
05/270 |
255 |
255 |
255 |
25 |
|
96/199 |
99/199 |
63/39 |
36/39 |
20/42 |
92/37 |
95/40 |
98/38 |
46/38 |
52/38 |
26 |
|
100 |
99/99 |
100 |
100 |
05/109 |
100 |
9/105 |
100 |
100 |
100 |
27 |
|
45/41 |
99/41 |
24/490 |
490 |
43/534 |
40/485 |
95/518 |
63/490 |
12/490 |
19/490 |
28 |
|
45/41 |
99/41 |
63/260 |
29/260 |
96/282 |
52/258 |
75/274 |
52/259 |
18/259 |
25/259 |
29 |
|
36/99 |
99/99 |
01/38 |
93/37 |
15/41 |
51/37 |
95/39 |
82/37 |
66/37 |
68/37 |
30 |
|
33/32 |
33 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
31 |
|
01/5 |
03/5 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
32 |
|
5 |
5 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
33 |
|
98/4 |
03/5 |
20/67 |
96/66 |
33/72 |
59/65 |
21/70 |
60/66 |
10/66 |
16/66 |
34 |
|
12/100 |
99/99 |
22/68 |
03/68 |
64/73 |
95/66 |
49/71 |
75/67 |
35/67 |
40/67 |
35 |
|
26/26 |
25 |
100 |
100 |
05/109 |
100 |
9/105 |
100 |
100 |
100 |
36 |
|
73/50 |
50 |
15/53 |
91/52 |
07/57 |
67/51 |
4/55 |
58/52 |
52 |
18/52 |
37 |
|
100 |
100 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
38 |
|
15/10 |
01/10 |
104 |
104 |
41/113 |
104 |
14/110 |
104 |
104 |
104 |
39 |
|
04/25 |
99/24 |
16/258 |
02/258 |
281 |
29/257 |
88/272 |
83/257 |
56/257 |
54/257 |
40 |
|
75/40 |
99/41 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
41 |
|
20/20 |
98/19 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
42 |
|
41 |
40 |
100 |
100 |
05/109 |
100 |
9/105 |
100 |
100 |
100 |
43 |
|
100 |
98/99 |
100 |
100 |
05/109 |
100 |
9/105 |
100 |
100 |
100 |
44 |
|
25/41 |
40 |
61/262 |
48/262 |
86/285 |
74/261 |
6/277 |
28/262 |
02/262 |
05/262 |
45 |
|
81 |
80 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
46 |
|
64/216 |
72/216 |
72/45 |
65/45 |
58/49 |
25/45 |
14/48 |
54/45 |
4/45 |
41/45 |
47 |
|
99/129 |
98/129 |
74/50 |
67/50 |
06/55 |
28/50 |
46/53 |
57/50 |
42/50 |
44/50 |
48 |
|
99/99 |
100 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
49 |
|
14/17 |
17 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
50 |
|
06/40 |
40 |
91/39 |
87/39 |
36/43 |
63/39 |
10/42 |
81/39 |
72/39 |
73/39 |
51 |
|
23/40 |
40 |
71/55 |
66/55 |
59/60 |
44/55 |
84/58 |
61/55 |
52/55 |
53/55 |
52 |
|
85/79 |
94/79 |
74/64 |
12/64 |
89/67 |
12/60 |
84/65 |
09/63 |
64/61 |
81/61 |
53 |
|
96/19 |
20 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
54 |
از جدول (2) می‌توان نتیجه گرفت با لحاظ ضرایب همبستگی، مقادیر توان‌های تولیدی، دستخوش تغییر می‌شوند و مقادیر تابع هدف تغییرات محسوسی را تجربه می‌کنند. با توجه به شرایط واقعی حاکم بر سیستم قدرت و تغییرات دما، سرعت باد و میزان تابش که تغییر اجباری میزان بار مصرف‌کنندگان را در پی خواهد داشت و سبب تغییر میزان انرژی حاصل‌پذیر از منابع انرژی باد و خورشید است، لحاظ ضرایب همبستگی واقعی براساس داده‌های همزمان اندازه‌گیری‌شده مربوط به این منابع، نقش حیاتی و مهمی در برنامه‌ریزی و بهره‌برداری سیستم قدرت خواهد داشت.
برای جمع‌بندی نتایج به‌دست‌آمده، در این بخش از مقاله تمامی اطلاعات به‌دست‌آمده از چهار روش مورد مطالعه به‌صورت جدول (3) خلاصه شدهاند.
جدول (3): مزایا و معایب روش‌های مورد مطالعه
عنوان روش |
مزایا |
معایب |
روش شبیه‌سازی مونت کارلو |
1- جواب بهینه تابع هدف (هزینۀ سوخت) در حد مطلوب است.2- پروفیل ولتاژ ژنراتورها دارای تغییرات زیادی است، اما در بازۀ حداقل و حداکثر مجاز قرار دارد. 3- چون این روش براساس تکرار به جواب نهایی می‌رسد، منحنی همگرایی برای حالت پخش توان بهینه رسم می‌شود. |
1- برای شروع محاسبات به نقاط اولیه زیاد نیاز است.2- برای رسیدن به جواب بهینه، به مدت زمان طولانی نیاز است. 3- برای حصول همگرایی به تکرارهای متعددی نیاز دارد. 4- ازنظر ریاضی، روابط پیچیده‌تری دارد. 5- با افزایش مقیاس شبکۀ مورد مطالعه، زمان لازم برای حصول به جواب بهینه افزایش می‌یابد. |
روش تخمین سه نقطه‌ای |
1- دارای نقاط اولیه محدود برای شروع محاسبات است.2- دارای بهترین جواب بهینه برای تابع هدف است. 3- زمان محاسبات به دلیل ساده‌بودن روابط ریاضی پایین است. 4- پروفیل ولتاژ ژنراتورها دارای تغییرات زیادی است، اما در بازۀ حداقل و حداکثر مجاز قرار دارد. |
1- به دلیل اینکه برای رسیدن به جواب نهایی از روش تکرار استفاده نمی‌کند، امکان رسم منحنی همگرایی برای این روش وجود ندارد. |
روش تبدیل بی بو |
1- جواب بهینه برای تابع هدف در حد مطلوب است.2- دقت جواب نهایی برای توان تولیدی ژنراتورها مشابه روش شبیه‌سازی مونت کارلو است. 3- پروفیل ولتاژ تولیدی ژنراتورها دارای تغییرات زیادی است، اما در کل این تغییرات در بازه حداقل و حداکثر قرار دارند. |
1- چون برای حصول جواب نهایی از روش تکرار استفاده نمی‌کند، امکان رسم منحنی همگرایی برای این روش وجود ندارد.2- زمان انجام محاسبات در این روش، در حد متوسط |
روش نقطۀ داخلی |
1- با توجه به اصلاح روابط ریاضی دارای سریع‌ترین حالت همگرایی در همان تکرارهای اولیه است.2- زمان انجام محاسبات در این روش بسیار پایین است. |
1- جواب بهینه در این روش بسیار بالا است. |
در این مقاله، برای حل عدم‌قطعیت‌های به‌وجودآمدۀ ناشی از وجود تولیدات نیروگاههای بادی و سلول‌های خورشیدی، از روش‌های احتمالاتی برای حل پخش توان بهینه استفاده می‌شود. این روش‌های احتمالاتی شامل شبیه‌سازی مونت کارلو، روش‌های تخمین نقطه‌ای، روش تبدیل بی بو و روش نقطۀ داخلی‌اند. در این مقاله، علاوه بر عدم‌قطعیت‌های نیروگاههای بادی و خورشیدی، عدم‌قطعیت تغییرات بار در نظر گرفته شده است.
این روش‌های حل، نقاط نمونۀ مناسب از عدم‌قطعیت متغیرهای ورودی را انتخاب می‌کند، سپس پخش توان بهینۀ احتمالاتی را با درجۀ بالایی از دقت و درستی انجام می‌دهند. روش‌های پیشنهادی روی شبکۀ استاندارد 118 شینه IEEE مطالعه شدهاند. با توجه به نتایج به‌دست‌آمده به‌صورت زیر دسته‌بندی می‌شوند:
در کل از میان روش‌های مورد مطالعه در شبکه 118 شینه IEEE براساس نتایج به‌دست‌آمده، روش تخمین سه نقطه‌ای دارای یک جواب بهینه ازنظر زمان، پروفیل ولتاژ و هزینۀ بهینه نسبت به سایر روش‌ها است.
ضمایم
فهرست علایم مربوط به معادله‌های (19) تا (21) در جدول (4) تعریف شده است.
حالت بهینگی روش کروش کان تاکر مرتبۀ اول برای عملکرد تابع لاگرانژ معادله (50) مربوط به روش نیوتون رافسون در سیستم غیرخطی به‌صورت معادلات (71) به دست می‌آید.
جدول (4): معرفی نمادهای استفاده‌شده در روابط ریاضی
|
توضیح |
نماد |
|
بردار توان اکتیو بار مصرف شده (بر حسب بر واحد) |
|
|
بردار توان اکتیو تولیدشده ژنراتورها (بر حسب بر واحد) |
|
|
بردار توان اکتیو تولیدشده با DER (بر حسب بر واحد) |
|
|
بردار توان راکتیو بار (بر حسب بر واحد) |
|
|
بردار توان راکتیو مصرف‌شده / تولیدشده با DER (بر حسب بر واحد) |
|
|
بردار توان راکتیو ژنراتورها (بر حسب بر واحد) |
|
|
بردار دامنۀ ولتاژ شین (بر حسب بر واحد) |
|
|
بردار عدم‌قطعیت متغیرهای ورودی |
|
|
بردار عدم‌قطعیت متغیرهای خروجی |
|
|
بردار زاویۀ ولتاژ شین (بر حسب رادیان) |
|
|
(71) |
|
با توجه به معادلات نیوتون، جهت برای هر متغیر در معادلات (72) تعریف شده است:
|
(72) |
|
مراحل و با جایگزینی معادلات (71) در (72) به‌ترتیب در معادلات (73) نشان داده می‌شوند.
|
(73) |
|
ساختار ماتریس (71) با جایگزینی معادلات (73) در (72) دوباره به‌صورت معادلات (72) بازنویسی می‌شود.
[1]تاریخ ارسال مقاله: 10/01/1398
تاریخ پذیرش مقاله: 18/04/1398
نام نویسندۀ مسئول: حمدی عبدی
نشانی نویسندۀ مسئول: ایران – کرمانشاه – بلوار طاق بستان – دانشگاه رازی – دانشکدۀ فنی و مهندسی
[1] Optimal Power Flow - OPF
[2] Monte Carlo Simulation - MCS
[3] Gram Charlier Series
[4] Point Estimate Method - PEM
[5] Probabilistic Optimal Power Flow - POPF
[6] Unscented Transformation (UT) Method
[7] Interior Point Methods - IPM
[8] Hong
[9] Harr
[10] Morales
[11] Weibull Distribution
[12] Rayleigh Distribution
[13] Beta Distribution
[14] Normal Distribution
)