نوشته‌ها

کسب درآمد اینترنتی

چگونه درآمد مستمر از سایت متلبی داشته باشیم

با سلام

وقتتون بخیر

اگر به فکر یک درآمد مداوم با کمترین زحمت هستید , راه کار خوبی برای شما در نظر گرفته ایم

همانطور که میدانیم سایت متلبی دارای رنک و رتبه خیلی خوبی در موتورهای جستجو است

ما این بستر آماده را برای همکاری بیشتر با شما بصورت زیر در اختیارتان قرار می دهیم:

آموزش های خود را اعم از تخصصی و غیر تخصصی, یا متلب و غیر متلب, تهیه و بر روی سایت متلبی بفروش برسانید

فایلهای آموزشی خود را رمزگذاری میکنید و با پنل کاربری که به شما داده میشود خودتون محصولاتتان را بفروش میرسانید

پس از هر خرید, پیامک به شما داده میشود و شما با مراجعه به پنل کاربری خود, از خرید مشتری آگاه و پسورد را برای مشتری پیامک می کنید.

با اینکه در روند فروش شما سایت متلبی هیچگونه دسترسی به فایلهای آموزشی شما ندارد. بصورت مکمل به شما تعهد داده می شود که سایت و مدیریت آن هیچگونه بهره برداری و فروشی از فایلهای شما نخواهد داشت.

شاید با خودتان بگویید اگر فرد دیگری فایل آموزشی را خرید و در اینترنت منتشر یا برای فروش قرار داد چه میشه کرد؟ سایت متلبی دارای مجوز هست و اگر کسی این کار را انجام دهد میتوان از ان فرد یا سایت شکایت قانونی کرد

 نیازی به پرداخت هزینه به عنوان اشتراک و غیره نیست

——————

آموزش میتواند در قالب ویدیو یا ورد باشد

در آموزش های ویدیویی یا ورد باید در ابتدا فهرستی از مطالب کلیدی که قرار هست به اون پرداخته بشه قرار بدین

در مورد آموزش ویدیویی ابتدا یک فایل ویدیو از فهرست کار مطابق بالا تهیه و برای مدیریت سایت متلبی ارسال بفرمایید تا اگر نیاز به اصلاح در فرمت و ساختار ویدیو هست , اطلاع رسانی بشه.

در فایل ویدیو یا فایل آموزش هیچ اطلاعاتی از شما نباید وجود داشته باشد و از لوگوی سایت به صورت فیلم زیر استفاده شود:

نرم افزار ضبط دسکتاپ حرکات موس را واضح نمایش دهد و خروجی آن حجم پایینی اشغال کند. پیشنهاد ما:

zd soft screen recorder

bandicam

دانلود فیلم

محدودیتی در حجم فایل های تهیه شده از سوی شما نیست. فقط زمان فشرده سازی و رمزگذاری روی فایلها, بایستی به قسمتهای زیر 1گیگ تقسیم بندی شوند.

این صفحه آپدیت می شود

لطفا نظرات و پیشنهادهای خودتان را به اطلاع مدیریت سایت برسانید

باسپاس

شبیه سازی با متلب

تشخیص فونم ها(لب خوانی) با SVM در متلب

میخواهیم در این پست آموزشی تشخیص واج ها یا فونم ها را با متلب انجام دهیم. برای این منظور داده ها را برای هم آوا(واج یا فونم) 4بار تهیه و ذخیره می کنیم. این تعداد هر چه بیشتر باشد الگوریتم دقیق تر عمل می کند. چون در اینجا داده های ما کم هست از الگوریتم SVM که خروجی بهتری را می دهد استفاده می کنیم.

در این آزمایش برای هر آوا 4بار با اشخاص مختلف و با گرفتن تصویر ویدویی از لب خوانی آنها دیتای اولیه را تهیه می کنیم. بعد نیاز هست تا آنها را با عدد مشخص کنیم و به این صورت که پشت سر هم باشند. یعنی اگر 10فونم داریم. باید 40ویدیو باشد که اولی 1 و آخری 40 نام گذاری شده است

بهتر هست سریعتر به همراه برنامه متلب آن آموزش متلب را آغاز کنیم. در توضیحات زیر منظور از عکس همین ویدیو های تهیه شدست:

قسمت ابتدایی برا پاک کردن صفحه، حافظه و اشکال

clc

clear all

close all

[n p]=uigetfile;

گرفتن مسیر عکسها

numImg=input(‘Enter The Image Numbers’);

گرفتن تعداد عکسها

TraImg=input(‘Enter The Image Number for each class’);

گرفتن اینکه هر کلاس چند عکس را شامل می­باشد

imgbank=cell(numImg,1);

تعریف بانک تصویر

T1=cell(numImg,1);

T2=cell(numImg,1);

متغیرهای گرفتن طول و عرض دهان

for ii=1:numImg

خواندن عکسها و ذخیره در متغیر خاص

برای هر عکس

a=num2str(ii);

اعداد عکسهارا به متن تبدیل کردن

n=strcat(a,’.avi’);

مخلوط کردن عکس با پسوند avi

img=read(mmreader([p n]));

خواندن عکس با گرفتن مسیر و نام عکس

imgbank{ii,1}=(img);

ذخیره عکس در بانک تصاویر

[rr cc zz ff]=size(img);

imgg=cell(ff,1);

ایجاد متغیری برای نگهداری فریمها

خواندن و استخراج تصویر لب بر اساس پیکسلها

for jj=1:ff

برای هر فریم

imgh=imgbank{ii,1}(:,:,:,jj);

خواندن از بانک تصاویر

imgbw1=(imgh(:,:,1)>130 & imgh(:,:,1)<210 & imgh(:,:,2)>40 & imgh(:,:,2)<90 & imgh(:,:,3)>30 & imgh(:,:,3)<90);

پیکسهای خاصی را که مربوط به لب است جدا کردن

imgbw2=(imgh(:,:,1)>100 & imgh(:,:,1)<150 & imgh(:,:,2)>45 & imgh(:,:,2)<75 & imgh(:,:,3)>45 & imgh(:,:,3)<75);

پیکسهای خاصی را که مربوط به اطراف لب است جدا کردن

imgt=imgbw1+imgbw2;

ترکیب آنها

استخراج مسک برای حذف جاهای اضافی

[rr cc]=size(imgt);

b=zeros(rr,cc);

for iii=31:rr-30

for jjj=76:cc-75

b(iii,jjj)=sum(sum(imgt(iii-30:iii+30,jjj-75:jjj+75)));

با بلوکی به بزرگ 60 در 150 دنبال جایی هستیم که تصویر بیشترین 1 را دارد که موقعیت همان لب است

end

end

[aa bb]=find(b==max(max(b)));

بیشترین یک را پیدا کردن

c=zeros(size(b));

g1=min(aa);

g2=max(aa);

g3=min(bb);

g4=max(bb);

if g1<41

اگر minاز 41 کمتر باشد منظورکنار تصویر لب داشته باشیم

g1=41;

مینیمم را 41 بگیر

end

if g2>rr-40

مانند بالا

g2=rr-40;

end

if g3<121

g3=121;

end

if g4>cc-120

g4=cc-120;

end

c(g1-40:g2+40,g3-120:g4+120)=1;

imgt1=imgt.*c(1:rr,1:cc);

مسک را در تصویر ضرب و لب را استخراج کن

حذف نویز با مساحتهای کم

cc1 = bwconncomp(imgt1);

stats = regionprops(cc1, ‘Area’);

idx = find([stats.Area] > 500);

imgt2 = ismember(labelmatrix(cc1), idx);

اگر چیزهای دیگری باشند که مساحتشان کوچک باشد حذف کن

imgg{jj,1}=imgt2;

تصویر لب را در ماتریس ذخیره کن

end

t=[];

for jj=2:ff

u=imgg{jj,1}-imgg{jj-1,1};

استخراج حرکت دو تصویر پشت سر هم زمانی

[a b]=find(u==1);

اگر حرکت داشته باشیم

t=[t;[max(a)-min(a) max(b)-min(b)]];

عرض و طول این حرکات را ذخیره کن

end

T1{ii}=t(:,1);

عرضها در یک سلول

T2{ii}=t(:,2);

طولها در یک سلول

end

feature1=[];

for ii=1:numImg

برای هر تصویر

a=T1{ii};

b=T2{ii};

عرض و طولها را در دو متغیر قرار بده

feature1=[feature1 ; [std(a) std(b) median(a) median(b) a(fix(length(a)./2)-2) a(fix(length(a)./2)-1) a(fix(length(a)./2)) a(fix(length(a)./2)+1) a(fix(length(a)./2)+2) b(fix(length(a)./2)-2) b(fix(length(a)./2)-1) b(fix(length(a)./2)) b(fix(length(a)./2)+1) b(fix(length(a)./2)+2) ]];

فضای ویژگیها را بر اساس مدین، انحراف معیار و مقادیر وسط دو ماتریس بساز

%feature1=[feature1 ; [std(a) std(b) median(a) median(b) ]];

end

tt=max(feature1);

ماکزیمم ویژگیها محاسبه

feature=[];

for ii=1:numImg

feature=[feature ;feature1(ii,:)./tt];

همه ویژگیها به ماکزیمم تقسیم تا نرمالیزه شوند

end

out=[];

for ii=1:numImg./TraImg

out=[out ii.*ones(1,TraImg)];

تهیه ماتریس خروجیها بصورت 1*n

end

out1=zeros(size(19,76));

for ii=1:numImg

out1(out(ii),ii)=1;

تهیه ماتریس خروجیها بصورت m*n(مناسب شبکه عصبی)

end

%net=newff(feature’,out1,[13],{‘purelin’ ‘tansig’});

%net=train(net,feature’,out1);

%sim(net,feature(1,:)’)

[result] = multisvm(feature,out’,feature(1,:))

دسته بندی با SVM

به همین راحتی میتوان تشخیص فونم را با متلب و با الگوریتم SVM انجام داد.

فانکشنی را نیز برای این برنامه باید در نظر بگیرید به صورت زیر:

function [result] = multisvm(TrainingSet,GroupTrain,TestSet)

u=unique(GroupTrain);
numClasses=length(u);
result = zeros(length(TestSet(:,1)),1);

%build models
for k=1:numClasses
%Vectorized statement that binarizes Group
%where 1 is the current class and 0 is all other classe
G1vAll=(GroupTrain==u(k));
if G1vAll==ones(size( G1vAll))
G1vAll=[ G1vAll(1:end-1); 2];
end

models(k) = svmtrain(TrainingSet,G1vAll);
end

%classify test cases
for j=1:size(TestSet,1)
for k=1:numClasses
if(svmclassify(models(k),TestSet(j,:)))
break;
end
end
result(j) = k;
end

حالا کافیست برای اجرای برنامه هر بار یک فونم را بصورت زیر بخواهید تا تشخیص دهد:

[result] = multisvm(feature,out’,feature(7,:))

این دستور به این معنی هست که فونم 7ام مربوط به کدام کلاس هست که الگوریتم SVM این کار را انجام می دهد

شبیه سازی با متلب

آموزش انجام پروژه درس کنترل مدرن با متلب

به نام خدا

در این آموزش متلب قصد داریم تمرین درس کنترل مدرن را بررسی و در متلب شبیه سازی کنیم:

تمرین شماره یک بصورت زیر است:

modern-control1

برای حل این تمرین از مثال حل شده فصل 6 کتاب کنترل مدرن دکتر علی خاکی صدیق کمک می گیریم

در این تمرین تابع تبدیل حلقه باز سیستم (بدون سیستم جبران ساز) به صورت زیر است

modern-control2

از جایی که 0.5=ξ و3=nω  در این مسئله مشخص می باشد و با توجه به تابع تبدیل استاندارد که به صورت زیر است

modern-control3

با جایگذاری مقادیر داده شده  مسئله  در تابع تبدیل بالا تابع تبدیل برابر می شود  با

modern-control4

که قطب های آن به صورت زیر بدست می آید

modern-control5

در واقع ما با طراحی این کنترل کننده باید به این قطب ها برسیم

برای براورده ساختن مشخصه های مسئله وبرای تغییر مسیر نمودار مکان ریشه ها به گونه ای که از قطب های حلقه-بسته (ذکر شده در بالا)عبور کند، جبران سازی طراحی خواهیم کرد. برای تعیین جبران ساز پیش فاز زاویه لازم φ جهت افزودن به مجموع زوایای قطبها و صفرهای حلقه-باز با یکی از قطبهای غالب حلقه –بسته مطلوب را پیدا می کنیم.مجموع این زوایا و زاویه φ باید  باشد.در سیستم فعلی، زاویه  در قطب حلقه بسته مطلوب عبارتست از

modern-control6

از اینرو جبران ساز پیش فاز باید زاویه را فراهم آورد. موقعیت  قطب و صفری که این زاویه را فراهم می آورند،منحصر بفرد نیست. یک جبران ساز مناسب (که با سعی و خطا بدست آمده است)عبارتست از

modern-control7

بهره   از شرط دامنه تعیین می شود به عبارت دیگر از شرط دامنه داریم

modern-control8

با استفاده از شرط دامنه در قطب مورد نظر داریم

modern-control9

بهره  برابر است با

modern-control10

پس به طور کلی اطلاعات بدست آمده برابر است با

http://www.matlabi.ir

ابتدا برای مشاهده مکان هندسی تابع تبدیل حلقه – بسته بدون کنترل کننده به صورت زیر عمل می کنیم

modern-control12

s=tf(‘s’);

gs=10/(s*(s+1));

sys1 = feedback(gs,+1);

rlocus(sys1)

 مکان هندسی حلقه- باز بدون کنترل کننده

مکان هندسی حلقه- باز بدون کنترل کننده

مکان هندسی حلقه- بسته بدون کنترل کننده

مکان هندسی حلقه- بسته بدون کنترل کننده

برای بدست آوردن قطب های تابع تبدیل استاندارد به روش زیر عمل می کنیم

modern-control4

a=[1 3 9];

b=roots(a);

modern-control5

در ادامه برای تعیین زاویه مورد نظر از برنامه زیر استفاده شده

s=-1.5+2.59i;

z=10/(s^2+s);

r = abs(z);

theta = atan2(imag(z),real(z));

anglout = radtodeg(theta);

zaviye=180-anglout;

یعنی کنترل کننده ما باید زاویه بالا را جبران کند تا به جواب مطلوب برسیم  که با بدست آوردن مقادیر t1 وt2 مورد نظر از طریق سعی وخطا به صورت زیر به زاویه مطلوب می رسیم

t1=3.7222; %40.8934

t2=0.45;

z1=((t1*s)+1)/((t2*s)+1);

r1 = abs(z1);

theta1 = atan2(imag(z1),real(z1));

anglout1 = radtodeg(theta1);

  K بدست می آید با

k=1/abs(z1*z);

با دستورات زیر مکان هندسی ریشه های سیستم حلقه بسته را بدون کنترل کننده وهمینطور همراه با کنترل کننده در یک شکل رسم می شود

q=tf(‘s’);

gs=10/(q*(q+1));

gc=((k*(t1*q+1))/(t2*q+1));

sys1 = feedback(gs,1);

gg=series(gc,gs);

sys2=feedback(gg,1);

subplot(2,1,1); rlocus(sys1)

title(‘Without the controller’)

ylabel(‘halghe baste’)

subplot(2,1,2); rlocus(sys2)

title(‘With the controller’)

ylabel(‘halghe baste’)

که به صورت زیر است

modern-control15

مکان هندسی ریشه های سیستم حلقه- بسته را بدون کنترل کننده وهمینطور همراه با کنترل کننده

 

همینطور اگر بخواهیم مکان هندسی حلقه-باز را رسم کنیم به شکل زیر است

figure

subplot(3,1,1); rlocus(gs)

title(‘tabe tabdil’)

ylabel(‘halghe baz’)

subplot(3,1,2); rlocus(gc)

title(‘jobran saz’)

ylabel(‘halghe baz’)

subplot(3,1,3); rlocus(gg)

title(‘tabe tabdil va jobran saz’)

ylabel(‘halghe baz’)

         modern-control16

مکان هندسی حلقه-باز بدون کنترل کننده وهمینطور همراه با کنترل کننده

 

رسم پاسخ پله سیستم بدون کنترل کننده وبدون کنترل کننده  و همینطور اطلاعات پاسخ پله آن به صورت زیر است به این صورت انجام می شود

et=stepinfo(sys1)

et1=stepinfo(sys2)

figure

subplot(2,1,1)

plot(step(sys1))

subplot(2,1,2)

plot(step(sys2))

modern-control17

modern-control18

پاسخ پله سیستم بدون کنترل کننده وبدون کنترل کننده

 

در پیوست برنامه متلب این تمرین وجود دارد که مرحله به مرحله توضیحاتی که داده شد در آن انجام شده و همچنین مکان هندسی حالات مختلف نیز رسم شده است.

 

.

شبیه سازی با متلب

شناسایی سیستم غیر خطی ربات بازوی مسطح دو درجه آزادی توسط شبکه عصبی

در این آموزش متلب میخواهیم شناسایی سیستم غیر خطی ربات بازوی مسطح دو درجه آزادی توسط شبکه عصبی مصنوعی را با متلب انجام دهیم.

معادلات دینامیکی این ربات به صورت زیر است:

matlabi_robot2dof

ابتدا باید داده های لازم را برای آموزش, اعتبارسنجی, و آزمایش شبکه عصبی را تهیه کرد. برای این منظور می بایست ورودی های خاصی را به سیستم مورد نظر اعمال و خروجی آن که زوایای تتا1 و تتا2 است, را اندازه گیری کرد.

قبل از هر چیز نیاز هست تا پروژه متلب را دانلود کنید

matlab-file

براي اين کار ابتدا بايد سيستم را در نرم افزار matlab  شبيه سازي نمود تا بوسيله آن داده هاي لازم توليد شوند:

براي شبيه سازي برنامه t2.m  نوشته شده است:

clc

clear

%define intional value

L1=0.7;

L2=0.5;

m1=4;

m2=3;

g=9.81;

Tau1=1.5;

Tau2=2;

sim(‘fcn.slx’);

 

در اين برنامه تابع fcn.slx فراخواني مي شود

پروژه متلب

با اجراي اين برنامه مقادير Theta1 ,Theta2  به عنوان حروجي به صورت زير بدست مي آيد:

matlabi_robot2dof2 matlabi_robot2dof3

با توجه به شبيه سازي ديده مي شود که نتايج شبيه سازي شده داراي مقادير ثابتي نبوده و براي رسيدن به مقادير ثابت از کنترل کننده pid  استفاده مي کنيم

براي شبيه سازي برنامه t3.m  نوشته شده است:

clc

clear

%define intional value

L1=0.7;

L2=0.5;

m1=4;

m2=3;

g=9.81;

Tau1=1.5;

Tau2=2;

sim(‘fcn1.slx’);

 

در اين برنامه تابع fcn1.slx فراخواني مي شود

پروژه متلب

با اجراي اين برنامه مقادير Theta1 ,Theta2  به عنوان حروجي به صورت زير بدست مي آيد:

پروژه متلب matlabi_robot2dof6

اضافه نمودن شبکه عصبي:

پروژه متلب

با توجه به شبيه سازي بالا مقادير ورودي و خروجي براي شبکه عصبي آماده شده است و برنامه t4.m  را مي نوسيم:

 

clc

clear

%define intional value

L1=0.7;

L2=0.5;

m1=4;

m2=3;

g=9.81;

Tauc1=1:0.1:2;

Tauc2=3:0.1:4;

for i=1:size(Tauc1,2)

    Tau1=Tauc1(i);

    Tau2=Tauc2(i);

sim(‘fcn1.slx’);

Teta1(i)=Theta1.data(end);

Teta2(i)=Theta2.data(end);

end

Teta1;

Teta1;

% define inputs and outputs for traning

x1 =Tauc1 ; % inputs tau1

x2 =Tauc2 ; % inputs tau2

y1 =Teta1 ;% outputs  Thetadd1

y2 =Teta2 ;% outputs  Thetadd2

x=[x1;x2]

y=[y1;y2]

net = feedforwardnet(10);

net = train(net,x,y);

view(net)

y = net(x);

perf = perform(net,y,x);

 

 

 

توصيحات برنامه :

1- دادن مقادير اوليه

%define intional value

L1=0.7;

L2=0.5;

m1=4;

m2=3;

g=9.81;

 

 

 

2-فرض کرديم مقادير گشتاور 1 و گشتاور 2 به صورت زير تغيير مي کند با  يک حلقه مقادير نتايج را بعد از شبيه سازي در محيط سيمولينک با فراخواني تابع fcn1.xls  به دست آورده و ذحيره مي کنيم:

Tauc1=1:0.1:2;

Tauc2=3:0.1:4;

for i=1:size(Tauc1,2)

    Tau1=Tauc1(i);

    Tau2=Tauc2(i);

sim(‘fcn1.slx’);

Teta1(i)=Theta1.data(end);

Teta2(i)=Theta2.data(end);

end

Teta1;

Teta1;

3- تعريف متغيير ها براي  آموزش شبکه عصبي:

% define inputs and outputs for traning

x1 =Tauc1 ; % inputs tau1

x2 =Tauc2 ; % inputs tau2

y1 =Teta1 ;% outputs  Thetadd1

y2 =Teta2 ;% outputs  Thetadd2

x=[x1;x2]

y=[y1;y2]

 

4- ايجاد شبکه عصبي :

net = feedforwardnet(10);

net = train(net,x,y);

view(net)

y = net(x);

perf = perform(net,y,x);

 

اجراي برنامه t4.m  مشاهده خروجي ها:

پروژه متلب شبکه عصبی

با توجه به شکل بالا شبکه عصبي 2 ورودي 2 خروجي 2 لايه اصلي که لايه پنهان آن داراي 10 لايه مي باشد استفاده شده است.

صفحه اصلي مربوط به شبکه عصبي فوق:

پروژه متلب

با کليک بر روي گزينه performance   شکل زير بدست مي آيد:

شبکه عصبی

با کليک بر روي گزينه traning state   شکل زير بدست مي آيد:

شبکه عصبی متلبی

با کليک بر روي گزينه erore histogram   شکل زير بدست مي آيد:

آموزش متلب

با کليک بر روي گزينه    Regression شکل زير بدست مي آيد:

آموزش متلب

 

 

 

 

 

 

 

شبیه سازی با متلب

مدل فرایند گاورنر 4درجه آزادی جهت کنترل سرعت توربین

مدل فرایند گاورنر جهت کنترل سرعت توربین به فرم زیر درنظر گرفته شده است. همچنین برای جبرانسازی افت فرکانس سیستم از یک دروپ کنترل کننده انتگرالی برای حذف خطای حالت ماندگار استفاده می گردد.

معادلات سیستم در متلب

کد نوشته شده در نرم افزار متلب به صورت زیر می باشد. (Parameters.m)

 

clc

% Dynamic Model Parameters:

 

Tt = 0.3;

Th = 0.1;

Fs = 60;

R  = 0.05*Fs;

H  = 3;

D  = 1.2125;

 

% Dynamic Model;

s = tf(‘s’);

Gdroop = (-1/(R*s))

Tfgov  = (1/(Th*s+1));

Tftur  = (1/(Tt*s+1));

Tfgen  = (Fs/(2*H*s+D));

Plant  = Tfgen*Tftur*Tfgov*Gdroop;

t=0:0.1:20;

 

با تعریف ورودی پله و سیگنال سینوسی به صورت زیر نتایج به فرم زیر قابل نمایش است.

 

t=0:0.1:20;

u=sin(t);

step(Plant);

y = lsim(Plant,u,t);

 پاسخ پله سیستم

شکل 1- پاسخ پله سیستم

 پاسخ سیستم به ورودی سینوسی

شکل 2 – پاسخ سیستم به ورودی سینوسی

بر اساس روش های کنترل کلاسیک در اینجا PID به طراحی کنترل کننده مناسب جهت پایدارازی سیستم و ردیابی ورودی مرجع پرداختیم که کنترل کننده حاصل به شکل زیر می باشد. اساس تنظیم این روش مبتنی بر طراحی زیگلر- نیکولز می باشد.

 شبیه سازی سیستم با کنترل کنندهPID

شکل 3- شبیه سازی سیستم با کنترل کنندهPID

 پاسخ سیستم به ووردی پله در حضور کنترل کننده PID

شکل 4 – پاسخ سیستم به ووردی پله در حضور کنترل کننده PID

 پاسخ سیستم به ووردی سینوسی در حضور کنترل کننده PID

شکل 5 – پاسخ سیستم به ووردی سینوسی در حضور کنترل کننده PID

جهت طراحی مناسبترین کنترل کننده ، جبران ساز پیش فاز – پس فاز را بکار برده ایم. همچنین جهت تنظیم پارامترهای این کنترل کننده از IMC استفاده کرده ایم. با توجه به نتایج نشان داده شده در شکل زیر کنترلر به عنوان تنظیم کننده عملرد بهتری نسبت به حالت ردیابی دارد.

شبیه سازی سیستم با جبران ساز پیش فاز - پس فاز

شکل 6- شبیه سازی سیستم با جبران ساز پیش فاز – پس فاز

 پاسخ سیستم به ووردی پله در حضور جبران ساز پیش فاز - پس فاز

شکل 7 – پاسخ سیستم به ووردی پله در حضور جبران ساز پیش فاز – پس فاز

 پاسخ سیستم به ووردی سینوسی در حضور جبران ساز پیش فاز - پس فاز

شکل 8 – پاسخ سیستم به ورودی سینوسی در حضور جبران ساز پیش فاز – پس فاز

در این قسمت فضای حالت در متلب به صورت زیر تعریف می گردد:

% Convert Plant to State Space ;

Pss = ss(Plant);

[A,B,C,D] = ssdata(Pss)

پروژه متلب

شکل 9- شبیه سازی سیستم برای طرحی فیدبک حالت

در این حالت با توجه به طراحی انجام شده در فضای حالت در مقایسه با توابع تبدیل توانایی اعمال شرایط اولیه و بررسی پاسخ سیستم نسبت به شرایط اولیه در حالت ورودی صفر را داریم. این کنترل کننده به لحاظ عملکردی دارای سرعت پاسخ بیشتری می باشد.

همچنین با تنظیم بهره فیدبک حالت کنترلی هرچه قدر گین های بزرگتری انتخاب کنیم قطب های حلقه بسته سیستم در مکان دورتری قرار خواهد گرفت و در نتیجه سیستم پایدار تر خواهد شد.

 پاسخ سیستم به شرایط اولیه در حضور کنترل فیدبک حالت کامل

شکل 10 – پاسخ سیستم به شرایط اولیه در حضور کنترل فیدبک حالت کامل

در مرحله بعد بهینه سازی پارامترهای کنترل کننده با تعریف توابع هدف مختلفی در متلب قابل انجام است. اکثر روشهای بهینه سازی کلاسیک مبتنی بر سیگنال خطای سیستم می باشد . در روش ITAE با در نظرگیری انتگرال  قدرمطلق خطا پارامترهای k را طوری محاسبه می نماییم که مقدار تابع هدف مورد نظر حداقل گردد. علت استفاده از قدر مطلق در این معیار جلوگیری از صفر شدن اثر سطوح مثبت و منفی می باشد.

با توجه به فایل شبیه سازی که در شکل نمایش داده شده است و  تابع تعریف شده در فایل ITAE.m با اجرای دستور زیر فرایند بهینه سازی ضرایب PID آغاز می گردد.

>> global P I

P = var(1);

I = var(2);

[t,x,y] = sim(‘PIDcontrolTIAE’,50);

[zn,fval,exitflag]=fminsearch(‘ITAE’,[8,8],optimset(‘MaxIter’,100))

در مرحله بعد بهینه سازی پارامترهای کنترل کننده با تعریف توابع هدف مختلفی قابل انجام است. اکثر روشهای بهینه سازی کلاسیک مبتنی بر سیگنال خطای سیستم می باشد . در وش ITAE با در نظریری انتگرال قدرمطلق خطا پارامترهای k را طوری محاسبه می نماییم که مقدار تابع هدف مورد نظر حداقل گردد. علت استفاده از قدر مطلق در این معیار جلوگیری از صفر شدن اثر سطوح مثبت و منفی می باشد. با توجه به فایل شبیه سازی که در شکل نمایش داده شده است و تابع تعریف شده در فایل ITAE.m با اجرای دستور زیر فرایند بهینه سازی ضرایب PID آغاز می گردد. >> global P I P = var(1); I = var(2); [t,x,y] = sim('PIDcontrolTIAE',50); [zn,fval,exitflag]=fminsearch('ITAE',[8,8],optimset('MaxIter',100))

شکل 11- شبیه سازی سیستم برای تنظیم ضریب PID به روش ITAE

 پاسخ سیستم به ورودی پله در حضور کنترلPID و کاربرد روش ITAE

شکل 12 – پاسخ سیستم به ورودی پله در حضور کنترلPID و کاربرد روش ITAE

میخواهیم بر اساس الگوریتم رگولاتور بهینه خطی درجه دو LQR در متلب بهینه ترین گین ها را جهت فیدبک نمودن حالت های سیستم طراحی نماییم. از آنجایی که هرچه K را بزرگتر انتخاب کنیم قطبهای حلقه بسته در مکان دورتری از محور موهومی قرار می گیرند اما این مشکل نیز به وجود خواهد آمد که سیگنال های کنترل متناسب با این بهره های بزرگ ، بزرگ خواهند شد و با توجه به محدودیت عملگرهای سیستم امکان ناپایداری وجود خواهد داشت.

بنابراین با روش طراحی کنترلر بهینه بهترین بهره ممکن جهت جایابی قطبهای حلقه بسته سیستم را پیدا می کنیم.

 پاسخ سیستم به شرایط اولیه با در نظ گرفتن کنترل کننده LQR

شکل 13-  پاسخ سیستم به شرایط اولیه با در نظ گرفتن کنترل کننده LQR

با توجه به طراحی کنترل کننده فیدبک حالت سیستم نیاز به اندازه گیری همه حالت های سیستم جهت فیدبک نمودن دارد. اما در عمل تنها تعداد محدودی از خروجی ها قابل اندازه گیری هستند بنابراین با طراح رویت گر حالت به تخمین حالت های سیستم می پردازیم.

شبیه سازی رویتگر کامل

شکل 14 – شبیه سازی رویتگر کامل

مقایسه حالت های تخمینی (پایین) با حالت واقعی سیستم (بالا)

شکل 15- مقایسه حالت های  تخمینی (پایین)  با حالت واقعی سیستم (بالا)

با توجه به طراحی بهره رویتگر قطبهای حلقه بسته سیستم دینامیکی رویتگر می توانند در هر نقطه موردنظر ما قرار گیرند. در ایجا نیز مانند طراحی فیدبک حالت انتخاب بهره های بزرگ باعث دور شدن قطبهای رویتگر از محور موهومی و سریع شدن پاسخ خمین رویتگر خواهند شد.اما با توجه به وجود نویز اندازه گیری سنسورها بزرگ بودن مقدار بهره رویتگر باعث تقویت این نویز و کاهش عملکرد حالت های تخمین زده شده می گردد. بنابراین در اینجا نیز انتخاب بهینه بهره رویتگر مساله مهمی می باشد.

مقایسه حالت های تخمینی (پایین) با حالت واقعی سیستم (بالا) - پاسخ واگرا

شکل 16- مقایسه حالت های  تخمینی (پایین)  با حالت واقعی سیستم (بالا) – پاسخ واگرا

هدف از طراحی رویگر کاهش یافته ، کاهش مرتبه دینامیکی رویتگر با اجتناب کردن از تخمین حالت یکی از متغیرها که اندازه گیری آن با سنسور و هزینه پایین و بدون نویز مقدور می باشد است. اما در سیستم در نظرگرفته شده برای دینامیک گاورنر و با توجه به اینکه تنها یک خروجی قابل اندازه گیری برای اعمال به رویتگر و بهبود تخمین حالت های سیستم با بروز کردن مقادیر داریم بنابراین استفاده از رویتگر کاهش یافته به معنای اجتناب کردن از اعمال ورودی سنسور به رویتگر می باشد و این نکته مانع از بهبود عملکرد رویتگر با بروز نمودن پارامترها است.

وجود نویز در سنسور با توجه به فیدبک نمودن خروجی اندازه گیری شده به رویت گر باعث تخمین حالت های سیستم همراه با مقداری نویز و متعاقبا کاهش عملکرد سیستم رویتگر و در نتیجه سیستم کنترل کننده فیدبک حالت خواهد شد.

جهت دانلود برنامه متلب آموزش داده شده اینجا را کلیک کنید

governer_matlabi

 

 

 

 

 

 

 

 

شبیه سازی با متلب

کنترل ولتاژ شبکه توزیع با تغییر تپ چنجر و خازن گذاری در متلب

سیستم توزیع انرژی الکتریکی به همراه مصرف کننده های عمده و جزئی از طریق سیستم انتقال به ولتاژ فشار قوی متصل است . سطح ولتاژ در سیستم توزیع پایین است و در نتیجه اندازه جریان ها زیاد می باشد به همین دلیل تلفات اهمی در سیستم های توزیع در مقایسه با سیستم های انتقال از اهمیت بیشتری برخوردار است . مسئله کاهش تلفات و بهبود کارایی تحویل انرژی الکتریکی سیستم قدرت عمدتاً به بخش های توزیع الکتریکی برمی گردد .

هرساله هزینه هفتگی برای جبران سازی تلفات انرژی الکتریی در سیستم قدرت از طریق تولید مازاد برمصرف صرف می شود . آمار نشان می دهد تلفات انرژی برق در شبکه های انتقالی نیرو طی سالهای گذشته 3 الی 4 درصد و در بخش توزیع 11 الی 12 درصد است و این نیز دلیلی دیگر برای اهمیت وتوجه به مسئله تلفات در سیستم توزیع می باشد . لذا بایستی با شناخت عوامل ایجاد کننده تلفات در جهت حذف این عوامل و کاهش تلفات قدمهای موثری برداشته شود . یکی از منابع مهم ایجاد تلفات ، توان راکتیوهای جاری در خطوط می باشد

تنظيم كننده ولتاژ نـوع پلـه اي شـامل يـك اتـوترانس و يـكدستگاه تغيير تپ زير بار مي باشـد كـه در يـك واحـد كامـل ساخته شده است . تغييـر ولتــاژ از طــريق تغييـر تـپ هـاي اتوترانس به دست مي آيد. تنظيم كننده هاي گامي اسـتاندارد داراي محدوده تنظيم 10% ( محدوده 10% ± آنها به ترتيب 32 يا 16 گام) مـي باشـند.بعضــي از واحــدها شامـل يـك كليد معكوس كننده مـي باشـند كـه قـادر اسـت محدوده تنظيم ولتاژ 10% ± را به وجود آورد. تنظـيم كننـده هاي ولتاژ گامي جديد داراي محـدوده تنظيـم ولت ـاژ 10% مي باشنـد كه به گام هاي 5 و 1 درصد تقسيم شده انـد. تنظيم كننده هاي ولتاژ گامي (اتوبوستر) با اتوترانس ها تفاوت دارند يعني آن ها به جاي اين كه بر حسـب كيلـو ولـت آمپـركلي مدار درجه بندي شوند برحسـب ميـزان كيلوولـت آمپـرتنظيمي درجه بندي شده اند. براي مثال ، يـك تنظـيم كننـده ولتاژ گامي سه فاز با قدرت 750 كيلو ولت آمپر و با محـدودهتنظيمي 10% ± در يك مدار سه فـاز بـا قـدرت نـامي 7500 كيلوولت آمپر مورد استفاده قرار مي گيرد. ولتاژ يك فاز را به ميزان 10% مقدار نامي مثبت يا منفي تغيير خواهد داد. تنظـيم كننده هاي ولتاژ برحسب آمپر نيز درجه بنـدي مـي شـوند دركاربرد تنظيم كننده جريان نامي آن بايستي برابـر يـا بيشـتر ازجريان خط باشد.

دو نوع تنظيم كننده پله اي وجود دارد نـوع پسـت ونوع شبكه توزيع كه اغلـب مواقـع ، تنظـيم كننده هاي ولتاژ توزيع كننده هاي ولتاژ تغذيـه كننـده توزيـع ناميده مي شوند به صورت يك فازه با قـدرت نـامي 5/12 تـا167 كيلوولت آمپر و ولتاژ نامي 4/14 كيلو ولت وكمتر وجوددارد. واحدهاي كوچكتر كه اغلب به نوع خط موسوم هستند، بــراي نصــب روي تيرهاســاخته مــي شــوند. واحــدهاي بزرگتربراي نصب روي تير مناسب نيستند و در داخـل پسـتها و يا روي سكو قرار داده مي شوند. تمام واحـدها كامـلا” خـودكار و شامـل تمـام لـوازم مربوطه مي باشند .

درصورتي كه از تنظيم كننده هاي نوع توزيع در يك محـدوده تنظيـم كاهـش ي استفـاده شـود ، ظرفيـت آنهــا را مـي تـوان افزايش داد . باكاهش 5% ± محدوده تنظيم ، ظرفيت تـا160% افزايش مي باشد. جـدول زیر افـزايش ظرفيـت مجـاز را بـراي كاهش هاي مختلف در محدوده تنظيم بـه طـور فهرسـت واربيان مي نمايد. برتري افزايش ظرفيت كـه در يـك محــدوده تنظيـم كاهش ي به دست مـي آيـد قبـل از ايـن كـه رشـد بـارمستلزم يك تغيير وضعيت باشد عبارت از افزايش زمـان كـارتنظيم كننده است. يعني اين كه تـامين يـك محـدوده تنظـيمكوچكتر مجاز امكان پذير است .

 

توانائي عبور جريان

اتصال كوتاه برحسب واحد

جريان نامي

مداوم برحسب واحد

بزرگترين محدوده تنظيمي كه بايد به كار برد
0/1*25 براي 2 ثانيه 1/00 10% تنظيم درجهت مثبت

در نگاهی گذرا به نظر می رسد که خازن یک وسیله بسیار ساده یعنی دو صفحه فلزی است که با یک ماده عایقی دی الکتریم از هم جدا شده اند . خازن بخش متحرک ندارد و با اعمال فشار الکتریکی کار می کند . اما در واقع خازن قدرت یک وسیله پیچیده و کاملاً فنی است که در آن مواردی التریک بسیار نازک و فشارهای التریکی قوی با شیوه های فرآیندی کاملاً پیچیده به یکدیگر مرتبط اند . در گذشته بیشتر خازن های قدرت را از دو برگ آلومینیوم خالص می ساختند که دست کم با سه لایه کاغذ گرافت با اشباع شیمیایی از هم جدا می شوند .

خازن های قدرت در 30 سال اخیر بهبود زیادی یافته اند که هم به علت بهبود موادی الکتریک و با استفاده بهتر از آنها و هم به علت بهبود در تکنیکهای فرآیند است. اندازه خازنها از بهینه kVar 25-15 به پهنه ی kVar300-200 افزایش یافته است . ( گروه خازنها در اندازه های 300 تا kVar 1800 عرضه می شوند).]1[.

امروزه خازنهای قدرت در دسترس شرکتهای برق رسانی نسبت به 30 سال پیش بازده بیشتر و هزینه برکیلو وار کمتر دارند به طور کلی ، امروزه بیش از گذشته به خازنها توجه می شود که تا حدی به علت اضافه شدن بعد جدیدی در تحلیها یعنی اقتصاد متغیر است . در پاره ای شرایط ، حتی تعویض خازنهای قدیمی تر برمبنای ارزیابی اتلاف کمتر خازنهای جدید توجیه می شود . تکنولوژی خازنها ، به سمت طرحهای بسیار کم اتلافی متحول شده که ناشی ازساخت فیلمی (پوسه ای ) است . در نتیجه ، شرکتهای برق رسانی می توانند انتخاب خود را براساس ارزیابی اقتصادی اتلاف تکنولوژیهای خازنی موجود انجام دهند]2[.

خازنهای موازی

خازنهای موازی یعنی خازنهایی که به موازات خطها بسته می شوند ، فراوان در سیستم های توزیع به کار می روند . خازن های موازی ، توان یا جریان واکنشی را تامین می کنند تا مولفه ی ناهمفاز جریان مورد نیاز یک بار القایی را جبران کند از جهتی ، خازنهای موازی با کشیدن جریان پیش فازی که بخشی یا همه مولفه پس فازی جریان بار القایی را در نقطه نصب ، خنثی می کند مشخصه آن را اصلاح می کنند . بنابراین ، خازن موازی همان اثر خازن سنکرون یعنی ژنراتور یا موتور سنکرون (همزمان ) پر تحریک را دارد .با بکارگیری خازن موازی برای فیدر ، می توان جریان بار را کم کرد و ضریب توان مدار را بهبود بخشید . در نتیجه افت ولتاژ بین سرفیدر و بار کاهش می یابد ولی خازنهای موازی اثری برجریان یا ضریب توان مدار بعد از نقطه نصب خود ندارند .

 

افت ولتاژ را در فیدر، یا در خطهای انتقال کوتاهی که ضریب توان پس فاز دارند می توان چنین تقریب زد :

(2-1 )

matlab

R = کل مقاومت مدار فیدر ،

matlab= کل رئکتانس القایی مدار فیدر ،

matlab= مولفه اکتیو جریان ، A

matlab= مولفه راکتیو جریان که 90 درجه نسبت به ولتاژ عقبتر است ، A

اگر مانند شکل 2-1 ب خازن در محل مصرف (بار ) نصب شود ، افت ولتاژ حاصل را می توان تقریباً چنین حساب کرد :

(2-2)matlab

که در آن ،

matlabمولفه راکتیو جریان که 90 درجه نسبت به ولتاژ جلوتر است ، A

تفاوت بین افت ولتاژهای حاصل از معادله های (2-1) و (2-2) برابر افزایش ولتاژ ناشی از نصب خازن است که می تواند چنین بیان شود :

(2-3 )matlab

توجیه اقتصادی تنظیم ولتاژ و خازنگذاری در شبکه توزیع

جبران توان راکتیو در سیستم توزیع باعث کاهش انتقال جریان راکتیو از نیروگاهها به سیستم توزیع می گردد که در پی آن تلفات سیستم قدرت کاهش می یابد . هنگامی که توان راکتیو را فقط در نیروگاهها تامین کنیم ، اندازه همه اجزای سیتم ایمنی ژنراتور ، ترانسفورماتور ، خط انتقال و توزیع ، تجهیزات کلید خانه و حفاظت باید متناسبا افزایش یابد . خازنها می تواند با کاهش تقاضای توان راکتیو در پشت خود تا ژنراتور ، این شرایط را تعدیل کنند . جریان خط از پشت سر خازن تا تجهیزات تولید کاهش می یابد . در نتیجه اتلافها و بار راکتیو خطهای توزیع ، ترانسفورماتورهای پست ها و خطهای انتقال کاهش می یابد . نصب خازنها می تواند قابلیت ژنراتور و پست را برای بار اضافی ، بسته به میزان عدم تصحیح ضریب توان سیستم ، تا حداقل 30 درصد ، افزایش دهد و می تواند قابلیت تک تک مدارها را از دیدگاه افت ولتاژ ، تقریبا 30 تا 100 درصد افزایش دهد ]1[. از طرفی ، کاهش جریان ترانسفرماتور ، تجهیزات توزیع و خطها ، بار این دستگاههای با کیلو ولت آمپر محدود را کم می کند و متعاقباً نصب وسایل جدید را به تأخیر می اندازد .لذا ، منافع اقتصادی ایجاب می کند که خازنها در سیستم توزیع نصب شوند جهت بررسی توجیه اقتصادی خازنگذاری مسائلی همچون کاهش تلفات انرژی ، آزادسازی ظرفیت ، کاهش تلفات پیک و سود حاصل از خازنگذاری را بیان می کنند . در مقابل هزینه خرید ، نصب و نگهداری تجهیزات کنترل خازنها به عنوان هزینه های سرمایه گذاری مطرح می شوند .

آزاد سازی ظرفیت و تصحیح ضریب توان پست های توزیع

سیستم برق رسانی نوعی از ماههای تابستان ، بارراکتیو با ضریب توان 80 درصد دارد . بنابراین ،دربارهای توزیع نوعی ، جریان ،مطابق شکل 2-2 نسبت به ولتاژ پس فاز دارد . کسینوس زاویه بین جریان و ولتاژ تغذیه را ضریب توان مدار می گویند . اگر مولفه های اکتیو راکتیو جریان را در ولتاژ سربار ضرب کنیم . رابطه حاصل را می توان روی مثلثی به نام مثلث توانها مطابق شکل 2-2 ب نشان داد . شکل 2-2 ب رابطه موجود بین کیلو وات ، کیلو ولت آمپر و کیلو وار را نشان می دهد . توجه کنید که با افزودن خازن می توان مولفه توان راکتیو توان ظاهری بار را کم یا کلاً حذف کرد .

 

شکل 2-2 ( الف ) نمودار فاز برداری (ب) مثلث توانهای یک بار توزیع نوعی

شکل های 2-3 و 2-4 ، افزایش مولفه توان راکتیو را به ازای هر 10 درصد تغییر ضریب توان نشان می دهد . بنابراین شکل 3-2 ، حتی ضریب توان 80 درصد ، آنقدر بزرگ است که باعث افزایش 25 درصد توان ظاهری ( کیلو ولت آمپر) خط می شود . در این ضریب توان kVar 75 خازن ، برای حذف 75 کیلو وار مؤلفه راکتیو لازم است .

شکل 2-3 نمایش چگونگی افزایش توانهای ظاهری و راکتیو بر حسب ضریب توان بار ، ضمن ثابت نگه داشتن توان اکتیو .

شکل 2-4 نمایش تغییر توانهای اکتیو و راکتیو بر حسب ضریب توان یا بار ، ضمن ثابت نگه داشتن توان ظاهری

چنانکه پیشتر گفتیم ، تولید توان راکتیو در یک نیروگاه و رساندن آن به بارهایی که در فاصله دوری از آن قرار دارند، از نظر اقتصادی عملی نیست، اما می تواند به آسانی با نصب خازن در مراکز بار انجام شود.

شکل 2-5 ، تصحیح ضریب توان یک سیستم معین را نمایش می دهد . مطابق این شکل ، خازنها توان راکتیو پس فازی را منبع می کشند ، یعنی توان راکتیو پس فازی را برای بار فراهم می کنند . فرض می کنیم یک بار با توان اکتیو p ، توان راکتیو پس فاز Q و توان ظاهری S در ضریب توان پس فاز زیر تغذیه می شود .

شکل 2-5 نمایش تصحیح ضریب توان

هنگامی که خازن را به توزی مکان بار نصب کنیم ، ضریب توان از بهبود می یابد

بنابراین مطابق شکل 2-5 توان ظاهری و توان راکتیو (با تامین توان راکتیو، ) به ترتیب از به و از به کاهش یافته اند . بدیهی است که کاهش جریان واکنشی باعث کاهش کل جریان می شود ودر نتیجه اتلاف توان را پایین می آورد لذا تصحیح توان ، با پایین آوردن ظرفیت کیلو ولت آمپری وکاهش اتلافهای توان در تجهیزات بین نقطه ی نصب خازن ها و نیروگاه که شامل خطهای توزیع ، ترانسفورماتورهای پست و خطهای انتقال است ، باعث صرفه جویی در هزینه های سرمایه گذاری و سوخت می شود . ضریب توان اقتصادی ، ضریب توانی است که در آن ، منافع اقتصادی ناشی از افزودن خازنهای موازی ، درست با هزینه خازنها برابر شود . در گذشته ، ضریب توان اقتصادی ، حدود 95 درصد بود . امروزه هزینه های زیاد نیروگاه و سوخت ، ضریب توان اقتصادی را به سمت یک سوق داده است. ولی هر چه ضریب قدرت به یک نزدیکتر شود ، تاثیر خازنها در بهبود ضریب توان ، کاهش کیلو ولت آمپر انتقالی خطر ، افزایش ظرفیت بار ، یا کاهش اتلافهای مس خط با کاهش جریان خط ، به شدت کاهش می یابد . بنابراین هزینه نصب خازنهای لازم برای رساندن ضریب توان به یک ، هرچه به یک نزدیکتر شویم زیادتر است .

کاهش تلفات و تلفات پیک

با تأمین توان راکتیو در محل مصرف مشترکین ، جریان عبوری از خط کاهش یافته لذا تلفات خط کاهش می یابد . هم چنین با داشتن طول پریود مورد مطالعه میزان کاهش تلفات انرژی نیز محاسبه می گردد .حال با در نظر گرفتن هزینه هر کیلو وات ساعت تولید انرژی ، سود ناشی از کاهش تلفات انرژی در اثر خازنگذاری مشخص می شود . کاهش تلفات در زمان پیک شبکه منفعت خوبی به همراه دارد . با کاهش تلفات در پیک نیروگاهها از حدود نا می خود فاصله می گیرند در نتیجه نیاز به تولید کمتر می شود . ضمن اینکه با افزایش مشترکین جدید احداث نیروگاههای جدید به تعریق افتاده و در هزینه احداث نیروگاه جدید نیز صرفه جویی می شود .

بهبود پروفیل ولتاژ

درآمد شرکتهای برق رسانی به علت افزایش ولتاژ ناشی از خازنگذاری در سیستم که موجب افزایش مصرف انرژی می شود بالا می رود . این امر به ویژه برای مصرف کننده های خانگی درست است . افزایش مصرف انرژی به ماهیت وسایل بکار رفته بستگی دارد . مثلا مصرف انرژی لامپها با مربع ولتاژ افزایش می یابد . به عنوان مثال جدول2-1 درصد افزایش انرژی کیلو وات ساعتی اضافی را بر حسب نسبت میانگین ولتاژ پس از افزودن خازنها به میانگین ولتاژ پیش از افزودن خازنها نشان داده است .

میانگین ، پس از خازنگذاری V

میانگین پیش از خازنگذاری V

افزایش kWh

%

00/1 0
05/1 8
10/1 16
15/1 25
20/1 34
25/1 43
30/1 52

جدول 2-1

بنابراین منافع در آمد ناشی از کیلو وات ساعت افزوده مصرف انرژی را می توان چنین محاسبه کرد :

(2-8)

= درآمد سالانه افزوده ناشی از افزایش kWh مصرف انرژی ، سال / واحد پول

= kWh افزایش مصرف انرژی، %

BEC= مصرف انرژی ابتدایی kWh ( یامبنا ) ، سال / kWh

= هزنیه انرژی ، kWh / واحد پول

از جمله مزایای مهم دیگر اتوبوستر و خازنگذاری درشبکه توزیع و فوق توزیع می توان به آزاد سازی ظرفیت فیدرها و تجهیزات وابسته به آن ، تعویق یا حذف هزینه های سرمایه گذاری ناشی از بهبود ویا توسعه سیستم و آزاد سازی ظرفیت پست توزیع نام برد. علاوه بر این خازنگذاری سبب آزاد سازی ظرفیت تولید و انتقال در شبکه قدرت نیز می گردد . به طور خلاصه ، رگولاتور ولتاژ و خازنها وسیله بسیار موثری برای کاهش هزینه های صنعت برق رسانی باتوجه به افزایش پیوسته هزینه های سوخت و نیروگاه است. هرگاه شرکتهای برق رسانی بتوانند سرمایه گذاری جدید نیروگاهی را به تعویق اندازند یا حذف کنند و نیازمندیهای انرژی را کاهش دهند سود می برند . بنابراین خازنها به کمینه کردن هزینه های بهره برداری کمک می کنند و برق رسانی و مصرف کنندگان جدید را با کمترین سرمایه گذاری در سیستم ممکن می سازند . امروزه شرکت های برق رسانی آمریکا تقریباً به ازای هرkW2 ظرفیت تولید نصب شده kVar 1 خازنهای قدرت نصب کرده اند تا از منافع اقتصادی آن بهره مند شوند.

تعاریف و تقسیم بندی‌های مختلفی ازناپایداری ولتاژ ارائه شده است. در تقسیم‌بندی اول دقیقاًّ مثل ناپایداری زاویه بار این ناپایداری به دو دسته اغتشاش کوچک و بزرگ تقسیم می‌شود. در تحلیل اغتشاش کوچک از مدل خطی و در بزرگ از مدل غیر خطی استفاده می‌کنند.

https://upload.wikimedia.org/wikipedia/fa/thumb/b/bc/Un12.png/220px-Un12.png

انواع ناپایداری

یک سیستم قدرت در یک نقطه‌کار مشخص اگر پس از یک اغتشاش کوچک ولتاژ به همان مقدار قبل از اغتشاش باز گردد و یا در همسایگی آن آرام گیرد پایدار ولتاژ اغتشاش کوچک نامیده می‌شود. با توجه به تئوری پایداری سیستم‌های خطی این در صورتی امکانپذیر است که مدل خطی شده سیستم حول نقطه‌کار مورد نظر مقادیر ویژه‌ای با قسمت منفی داشته باشد. یک سیستم قدرت در یکنقطه‌کار مشخص چنانچه ولتاژ سیستم پس از اغتشاش بزرگ در نقطه تعادل بعد از خطا قرار گیرد پایدار اغتشاش بزرگ نامیده می‌شود با توجه به تئوری پایداری سیستم‌های خطی این در صورتی امکانپذیر است که مقدار خطای ولتاژ در ناحیه همگرایی نقطه تعادل باقی بماند. یک سیستم قدرت در یک نقطه‌کار مشخص و تحت یک اغتشاش معین چنانچه نقطه تعادل بعد از خطای آن به زیر حد قابل قبول برسد به فروپاشی ولتاژ می‌انجامد و این زمانی است که ولتاژ به مرز 0.7 پریونیت کاهش یابد.

توضیحات برنامه main1

clc;clear ;close all

اطلاعات باس و خطوط شبکه 33 باسه ieee

Image result for ieee 33 bus system

% basemva = 1; accuracy = 0.001; accel = 1.8; maxiter = 100;

% IEEE 30-BUS TEST SYSTEM (American Electric Power)

% Bus Bus Voltage Angle —Load—- ——-Generator—– Static Mvar

% No code Mag. Degree MW Mvar MW Mvar Qmin Qmax +Qc/-Ql

busdata=[

1 1 1.0 0.0 0.000 0.000 0.0 0.0 0 0 0

2 0 1.0 0.0 0.100 0.060 0.0 0.0 0 0 0

3 0 1.0 0.0 0.090 0.040 0.0 0.0 0 0 0

4 0 1.0 0.0 1.120 0.080 0.0 0.0 0 0 0

5 0 1.0 0.0 0.060 0.030 0.0 0.0 0 0 0

6 0 1.0 0.0 1.060 0.020 0.0 0.0 0 0 0

7 0 1.0 0.0 0.200 0.100 0.0 0.0 0 0 0

8 0 1.0 0.0 0.200 0.100 0.0 0.0 0 0 0

9 0 1.0 0.0 0.060 0.020 0.0 0.0 0 0 0

10 0 1.0 0.0 0.060 0.020 0.0 0.0 0 0 0

11 0 1.0 0.0 0.045 0.030 0.0 0.0 0 0 0

12 0 1.0 0.0 0.060 0.035 0.0 0.0 0 0 0

13 0 1.0 0.0 0.060 0.035 0.0 0.0 0 0 0

14 0 1.0 0.0 0.120 0.080 0.0 0.0 0 0 0

15 0 1.0 0.0 0.060 0.010 0.0 0.0 0 0 0

16 0 1.0 0.0 0.060 0.020 0.0 0.0 0 0 0

17 0 1.0 0.0 0.060 0.020 0.0 0.0 0 0 0

18 0 1.0 0.0 0.090 0.040 0.0 0.0 0 0 0

19 0 1.0 0.0 0.090 0.040 0.0 0.0 0 0 0

20 0 1.0 0.0 0.090 0.040 0.0 0.0 0 0 0

21 0 1.0 0.0 0.090 0.040 0.0 0.0 0 0 0

22 0 1.0 0.0 0.090 0.040 0.0 0.0 0 0 0

23 0 1.0 0.0 0.090 0.050 0.0 0.0 0 0 0

24 0 1.0 0.0 0.420 0.200 0.0 0.0 0 0 0

25 0 1.0 0.0 0.420 0.200 0.0 0.0 0 0 0

26 0 1.0 0.0 0.060 0.025 0.0 0.0 0 0 0

27 0 1.0 0.0 0.060 0.025 0.0 0.0 0 0 0

28 0 1.0 0.0 0.060 0.020 0.0 0.0 0 0 0

29 0 1.0 0.0 0.120 0.070 0.0 0.0 0 0 0

30 0 1.0 0.0 0.200 0.100 0.0 0.0 0 0 0

31 0 1.0 0.0 0.150 0.070 0.0 0.0 0 0 0

32 0 1.0 0.0 0.210 0.100 0.0 0.0 0 0 0

33 0 1.0 0.0 0.060 0.040 0.0 0.0 0 0 0];

% Line code

% Bus bus R X 1/2 B = 1 for lines

% nl nr p.u. p.u. p.u. > 1 or < 1 tr. tap at bus nl

linedata=[01 2 0.0922 0.0470 0 1

02 03 0.4930 0.2512 0 1

03 04 0.3661 0.1864 0 1

04 05 0.3811 0.1941 0 1

05 06 0.8190 0.7070 0 1

06 07 0.1872 0.6188 0 1

07 08 0.7115 0.2351 0 1

08 09 1.0299 0.7400 0 1

09 10 1.0440 0.7400 0 1

10 11 0.1967 0.0651 0 1

11 12 0.3744 0.1298 0 1

12 13 1.4680 1.1549 0 1

13 14 0.5416 0.7129 0 1

14 15 0.5909 0.5260 0 1

15 16 0.7462 0.5449 0 1

16 17 1.2889 1.7210 0 1

17 18 0.7320 0.5739 0 1

02 19 0.1640 0.1565 0 1

19 20 1.5042 1.3555 0 1

20 21 0.4095 0.4784 0 1

21 22 0.7089 0.9373 0 1

03 23 0.4512 0.3084 0 1

23 24 0.8980 0.7091 0 1

24 25 0.8959 0.7071 0 1

06 26 0.2031 0.1034 0 1

26 27 0.2842 0.1447 0 1

27 28 1.0589 0.9338 0 1

28 29 0.8043 0.7006 0 1

29 30 0.5074 0.2585 0 1

30 31 0.9745 0.9629 0 1

31 32 0.3105 0.3619 0 1

32 33 0.3411 0.5302 0 1];

پخش بار شبکه فوق بدون تاثیر خازن و تپگزاری

[Vm1,Loss1]=runpf(busdata,linedata);

تعیین مکان باس خازن و ظرفیت خازن

capaBus=[4 7 13];

capaVar=[0.8 0.4 0.2];

و پخش بار مجدد شبکه تغییریافته اعمال خازن فوق به شبکه قدرت اصلی

busdata(capaBus,end)=capaVar;[Vm2,Loss2]=runpf(busdata,linedata);

تعیین شماره خط و درصد تپ ترانس در ورودی خط انتقال

TapLine=[10 26];

TapSet=[ 0.98 0.98];

و انجام پخش بار اعمال تپ فوق الذکر برشبکه تغییریافته

linedata(TapLine,end)=TapSet;[Vm3,Loss3]=runpf(busdata,linedata);

رسم نمودار های ولتاژ برای سه حالت فوق

subplot(8,1,1:5)

plot(Vm1,’-o’,’linewid’,2)

hold on

plot(Vm2,’-o’,’linewid’,2)

plot(Vm3,’-o’,’linewid’,2)

% set(gca,’ylim’,[0.8 1.2])

xlabel(‘Bus Index’)

ylabel(‘Voltage (PU)’)

hold off

legend(‘Normal’,’With Capa’,’With Capa+Tap’)

grid on

رسم نمودار های تلفات برای سه حالت فوق

subplot(8,1,7:8)

bar([Loss1 Loss2 Loss3])

set(gca,’xticklabel’,{‘Normal’,’With Capa’,’With Capa+Tap’})

set(gca,’ylim’,[0.4 0.7])

ylabel(‘Losses (PU)’)

نتایج برنامه

در حالت عادی ولتاژ تا 0.91 پریونیت کاهش یافته و با افزودن خازن در باس 4و 7و 13 پروفیل ولتاژ بهبود یافته و با تنظیم تپ ترانس در خطوط 10 و 26 پروفیل ولتاژ در باسهای انتهایی شبکه بهبود یافته است و تلفات هم به نسبت کاهش یافته است.

البته میتوان بسهولت در برنامه مکان و ظرفیت خازن و تپ ترانس را تغییر داد و نتایج را با اجرای برنامه مشاهده کرد

این پروژه یکسری فانکشن های کمکی دارد که توسط برنامه mains1 فراخوانی میشود و بخودی خود قابلیت اجرا ندارند

 

 

 

 

 

 

 

 

 

شبیه سازی با متلب

ماشین های القایی متقارن در متلب

نظریه ماشین های القایی متقارن در متلب

  analysis of electric machinery exercises in Matlab

برای شبیه سازی یک موتور القایی در متلب روشهای گوناگونی می تواند مطرح باشد, از آنجایی که معادلات حاکم بر این موتورها در اختیار است, می توانیم با انتگرال گیری عددی به جوابهای قابل قبولی برسیم.

در ابتدا با استفاده از الگوریتم انتگرال گیری رانگ کوتاه معادلات موتور القایی دلخواهی را در فضای مرجع دلخواه با متلب شبیه سازی کرده و نتایج را مشاهده می کنیم.

clear all
clc
J=.089;
rr=.816;
t=0;
i=0;
dt=.0001;
P=4;
Tl=11.9/2;
Te=0;
tet0=0;tetr0=0;tet=0;tetr=0;
rs=.435;
Wr0=0;
Xlr=.754;
Xls=.754;
Wb=120*pi;
Xm=26.13;
Xss=Xls+Xm;
Xrr=Xlr+Xm;
D=(Xss*Xrr)-(Xm^2);
vqr=0;vdr=0;vor=0;
Q0=[0;0;0;0;0;0];
Qqs=0;Qds=0;Qos=0;Qqr=0;Qdr=0;Qor=0;
while t<=3
Vas=220*(.8165)*cos(Wb*t);
Vbs=220*(.8165)*cos((Wb*t)-(2*pi/3));
Vcs=220*(.8165)*cos((Wb*t)+(2*pi/3));
Ks=(2/3)*[cos(tet0) cos((tet0)-(2*pi/3)) cos((tet0)+(2*pi/3));sin(tet0) sin((tet0)-(2*pi/3)) sin((tet0)+(2*pi/3));.5 .5 .5];
V1=Ks*[Vas;Vbs;Vcs];
vqs=V1(1,1);vds=V1(2,1);vos=V1(3,1);
Q=[Qqs;Qds;Qos;Qqr;Qdr;Qor];
V=[vqs;vds;vos;vqr;vdr;vor];
Wr=Wr0;
W=Wr;
A=[rs*Xrr/D (W/Wb) 0 (-rs*Xm/D) 0 0;(-W/Wb) rs*Xrr/D 0 0 -rs*Xm/D 0;0 0 rs/Xls 0 0 0;-rr*Xm/D 0 0 rr*Xss/D (W-Wr)/Wb 0;0 (-rr*Xm/D) 0 (-(W-Wr)/Wb) rr*Xss/D 0;0 0 0 0 0 rr/Xlr];
B=[1/Wb 0 0 0 0 0;0 1/Wb 0 0 0 0;0 0 1/Wb 0 0 0;0 0 0 1/Wb 0 0;0 0 0 0 1/Wb 0;0 0 0 0 0 1/Wb];
Brev=inv(B);
Arev=-inv(B)*A;
Q01=Arev*Q+Brev*V;
Te=(P/2)*(3/2)*(Xm/(D*Wb))*((Q(1,1)*Q(5,1))-(Q(2,1)*Q(4,1)));
Wr01=(P/2)*((Te-Tl)/J);
tet01=W;
tetr01=Wr;
Wr=Wr0+(Wr01*dt/2);
tet=tet0+(tet01*dt/2);
tetr=tetr0+(tetr01*dt/2);
Q=Q0+(Q01*dt/2);

Q02=Arev*Q+Brev*V;
Te=(P/2)*(3/2)*(Xm/(D*Wb))*((Q(1,1)*Q(5,1))-(Q(2,1)*Q(4,1)));
Wr02=(P/2)*((Te-Tl)/J);
Wr=Wr0+(Wr02*dt/2);
W=Wr;
tet02=W;
tetr02=Wr;
tet=tet0+(tet02*dt/2);
tetr=tetr0+(tetr02*dt/2);
Q=Q0+(Q02*dt/2);

Q03=Arev*Q+Brev*V;
Te=(P/2)*(3/2)*(Xm/(D*Wb))*((Q(1,1)*Q(5,1))-(Q(2,1)*Q(4,1)));
Wr03=(P/2)*((Te-Tl)/J);
Wr=Wr0+(Wr03*dt);
W=Wr;
tet03=W;
tetr03=Wr;
tet=tet0+(tet03*dt);
tetr=tetr0+(tetr03*dt);
Q=Q0+(Q03*dt);

Q04=Arev*Q+Brev*V;
Te=(P/2)*(3/2)*(Xm/(D*Wb))*((Q(1,1)*Q(5,1))-(Q(2,1)*Q(4,1)));
Wr04=(P/2)*((Te-Tl)/J);
Wr0=Wr0+(Wr01+2*Wr02+2*Wr03+Wr04)*(dt/6);
W=Wr0;
tet04=W;
tetr04=Wr;
tet0=tet0+(tet01+2*tet02+2*tet03+tet04)*(dt/6);
tetr0=tetr0+(tetr01+2*tetr02+2*tetr03+tetr04)*(dt/6);
Q0=Q0+(Q01+Q02*2+2*Q03+Q04)*(dt/6);
Qqs=Q0(1,1);Qds=Q0(2,1);Q0s=Q0(3,1);Qqr=Q0(4,1);Qdr=Q0(5,1);Q0r=Q0(6,1);
Te=(P/2)*(3/2)*(Xm/(D*Wb))*((Q(1,1)*Q(5,1))-(Q(2,1)*Q(4,1)));
I=(1/D)*[Xrr 0 0 -Xm 0 0;0 Xrr 0 0 -Xm 0;0 0 D/Xls 0 0 0;-Xm 0 0 Xss 0 0;0 -Xm 0 0 Xss 0;0 0 0 0 0 D/Xlr]*Q;
iqs=I(1,1);ids=I(2,1);ios=I(3,1);iqr=I(4,1);idr=I(5,1);ior=I(6,1);
Kr=(2/3)*[cos((tet0-tetr0)) cos(((tet0-tetr0))-(2*pi/3)) cos(((tet0-tetr0))+(2*pi/3));sin((tet0-tetr0)) sin(((tet0-tetr0))-(2*pi/3)) sin(((tet0-tetr0))+(2*pi/3));.5 .5 .5];
Irotor=inv(Kr)*[iqr;idr;ior];
iar=Irotor(1,1);ibr=Irotor(2,1);icr=Irotor(3,1);
Ks=(2/3)*[cos(tet0) cos((tet0)-(2*pi/3)) cos((tet0)+(2*pi/3));sin(tet0) sin((tet0)-(2*pi/3)) sin((tet0)+(2*pi/3));.5 .5 .5];
Istator=inv(Ks)*[iqs;ids;ios];
ias=Istator(1,1);ibs=Istator(2,1);ics=Istator(3,1);
i=i+1;
time(i)=t;
Iar(i)=iar;
Tee(i)=Te;
Ias(i)=ias;
WW(i)=(30*Wr0)/(2*pi);
t=t+.0001;
if t>=1.5
Tl=11.9;
end
end
plot(time,WW);xlabel(‘t’) ;ylabel(‘Wr’)
grid on
figure
plot(time,Ias);xlabel(‘t’) ;ylabel(‘ias’)
grid on
figure
plot(time,Iar);xlabel(‘t’) ;ylabel(‘iar’)
grid on
figure
plot(time,Tee);xlabel(‘t’) ;ylabel(‘Te’)
grid on
figure
plot(WW,Tee);xlabel(‘wr’) ;ylabel(‘Te’)
grid on

با توجه به متن برنامه متلب در ابتدا موتور را بی بار راه اندازی می کنیم و نتیجه ای که از خروجی متلب مشاهده می کنیم:

در این شکل می بینیم که سرعت گردش روتور از صفر تا 1800rpm افزایش می یابد و بعد ثابت می ماند, چون به سرعت نامی خود رسیده است. Wb

جریان استاتور با متلب

در شکل فوق جریان استاتور را در طول زمان اجرای برنامه متلب می بینیم.

جریان روتور با متلب

شکل بالا تغییرات جریان روتور را به نمایش گذاشته است

گشتاور موتور القایی با متلب

گشتاور موتور القایی با متلب

در شکل بالا گشتاور الکترومغناطیسی بعد از نوسانات ابتدایی به صفر می رسد, حالا اگر بخواهیم گشتاور را بر حسب سرعت در متلب داشته باشیم به صورت زیر به دست می آید:

matlabi_elghayi55

بعد از بررسی ویژگی ماشین القایی مادامی که بی بار کار می کند میتواند بار نامی مکانیکی را روی محور ماشین القایی قرار داد. برای اینکار در برنامه متلب چند خط زیر را اضافه می کنیم:

if t>=1.5
Tl=11.9;
end

با اضافه کردن این قسمت از برنامه که در بالا با رنگ سبز نیز مشخص شده, این اتفاق در متلب می افتد و نتایج به صورت زیر است:

راه اندازی موتور القایی با متلب

راه اندازی موتور القایی با متلب

همانطور که از شکل فوق مشخص هست در حالتی که بار روی محور ماشین القایی قرار می گیرد سرعت 1800rpm نمی شود و کمی افت دارد.

جریان استاتور با متلب

جریان استاتور با متلب

جریان روتور با متلب

جریان روتور با متلب

دلیل نوسانات در جریان روتور و استاتور قرار گرفتن بار نامی روی  محور ماشین است

matlabi_elghayi4

همانطور که مشاهده می شود گشتاور مغناطیسی در یک مقدار غیر از صفر ثابت می ماند.

matlabi_elghayi5

 

بزودی بخش دیگری از تحلیل رفتار موتور القایی را در متلب در سایت متلبی خواهید دید.

شبیه سازی با متلب

مدلسازی موتور شنت با متلب

با استفاده از مشخصات و پارامترهای یک موتور شنت شبیه سازی را انجام می دهیم:

clear all
clc
ra=.6;
rf=240;
laa=.012;
lff=120;
laf=1.8;
Ta=laa/ra;
Tf=lff/rf;j=1;bm=0;bl=.2287;
disp(‘type your simulink program name and then type return’)
keyboard
subplot(311),plot(y(:,1),y(:,2)),title(‘ia,A’),grid
subplot(312),plot(y(:,1),y(:,3)),title(‘Te,n.m’),grid
subplot(313),plot(y(:,1),y(:,4)),title(‘if,A’),grid

در این برنامه متلب از دستور keyboard استفاده شده چرا که وقتی برنامه به این خط از برنامه متلب می رسد متوقف می شود و در حالت standby باقی می ماند تا ما اسم برنامه سیمولینک مورد نظر خود را می نویسیم, بعد از اجرای برنامه سیمولینک به کامنت ویندوز متلب برمی گردیم و return را می نویسیم و اینتر میزنیم. حالا دستورات plot شکلهای خروجی را رسم می کنند. بلوک دیاگرام موتور شنت به صورت زیر است:

matlabi_shant motor1

در شکل فوق گشتاور بار را صفر در نظر گرفته ایم که در ادامه غیر صفر بودن آنرا نیز بررسی خواهیم کرد.

matlabi_motor shant2

یکی از مشهودترین اتفاق در شکل بالا جریان راه اندازی بالا در لحظه اولیه هست, بخاطر وجود چنین جریان مخربی است که موتور را بصورت پله ای راه اندازی می کنند. وقوع گشتاور بالا در لحظه راه اندازی نیز بخاطر وجود جریان بالاست. حالا در ادامه فرض می کنیم که بعد از چند لحظه که موتور راه افتاد بار مکانیکی را روی محور ماشین قرار می دهیم. قرار دادن بار را روی ماشین در متلب با step انجام می دهیم. step را طوری تنظیم می کنیم که در لحظه مثلا 4ثانیه بار نامی که از تقسیم توان نامی بر سرعت نامی بدست می آید روی محور قرار گیرد. مقدار گشتاور بار معادل 29.195 بدست می آید. من عمدا مقدار گشتاور را دوبرابر قرار دادم تا exaggerate کرده باشم و شکلها واضح تر باشند. در حقیقت تاثیر بار را بهتر خواهیم دید. نتایج عبارت اند از:

matlabi_motor shant3

تغییرات در لحظه 4ثانیه بوضوح مشهود است.مدلسازی موتور شنت با متلب را نیز انجام دادیم 🙂

 

شبیه سازی با متلب

شبیه سازی رفتار یک رله مغناطیسی با متلب

در این بخش از آموزش های متلب قصد داریم تا یک سیستم الکترومکانیکی با میدان واسط مغناطیسی را با متلب تحلیل کنیم که مدل ابتدایی یک رله الکترومغناطیسی می تواند باشد که البته فقط خطی حرکت می کند.

برای اینکار از دومعادله دیفرانسیل جهت شبیه سازی با متلب استفاده می کنیم که یکی از معادلات حاکم بر قسمت الکتریکی سیستم و دیگری حاکم بر قسمت مکانیکی آن است. بایستی برای شبیه سازی از روش انتگرال گیری استفاده کنیم و در اینجا از روش رنگ کوتاه در متلب استفاده کردیم. برای اجرای انتگرال گیری از یک فانکشن در متلب استفاده کردیم که به صورت زیر است:

function z = out(s)
global t
i1=s(1);
x=s(2);
v=s(3);
we=s(4);
f=s(5);
%……………………………………..
xo=.003;
r=10;
k=2667;
u=5;
if t>.15
u=0;
end
m=.055;
d=4;
c=6.293e-5;
%……………………………………..
di1=(47.67)*(u-r*i1);
dx=v;
dv=((-c/2)*((i1/x)^2)-(d*v)-k*(x-xo))/m;
dwe=(u*i1-r*i1^2);
df=(-.5)*c*((i1/x)^2)*v;
%……………………………………..
z(1)=di1;
z(2)=dx;
z(3)=dv;
z(4)=dwe;
z(5)=df;
%……………………………………..
end

با توجه به ام فایل بالا می بینیم که s بعنوان یک بردار تعریف شده است, ضمن اینکه چون معادله مکانیکی حاکم بر سیستم مشتق دوم نسبت به زمان دارد دو حالت برای آن تعریف کرده ایم و اما برای معادله الکتریکی i را یک حالت در نظر گرفتیم.

برنامه اصلی در متلب نیز به صورت زیر است:

clear all
clc
global t
%%…………………………
dt=.0001;
t=0;
i1=0;
x=3e-3;
v=0;
we=0;
f=0;
s=[i1,x,v,we,f];
n=1;
c=6.293e-5;
%………………………….
while (t<.25)
A=out(s);
B=out(s+A*dt/2);
C=out(s+B*dt/2);
D=out(s+C*dt);
s=s+(A+B*2+C*2+D)*dt/6;
%………………………
if (t<.025)
s=[0,3e-3,0,0,0];
end
%………………………
n=n+1;
t=t+dt;
ii1(n)=s(1);
xx(n)=s(2);
vv(n)=s(3);
wee(n)=s(4);
ff(n)=s(5);
tt(n)=t;

end

figure
plot(tt,xx);xlabel(‘t’) ;ylabel(‘X’);
figure
plot(tt,ii1);xlabel(‘t’) ;ylabel(‘i1’);
figure
plot(tt,vv);xlabel(‘t’) ;ylabel(‘v’);
figure
plot(tt,wee);xlabel(‘t’) ;ylabel(‘we’);
figure

plot(tt,ff);xlabel(‘t’) ;ylabel(‘f’);

نتایج را در زیر مشاهده می کنیم:

matlabi_electromagnetic1

در شکل بالا تغییرات x را مشاهده می کنیم و در ادامه تغییرات جریان را می بینیم:

matlabi_electromagnetic2 matlabi_electromagnetic3 matlabi_electromagnetic4 matlabi_electromagnetic5

شبیه سازی با متلب

پیش بینی سری زمانی به کمک شبکه عصبی در متلب

خلاصه مبحث آموزشی :

پیشگویی سری زمانی یک مساله مهم کاربردی است که در بسیاری از زمینه ها نظیر پیش بینی وضع هوا، اقتصاد ، کنترل و … استفاده می شود . در این پست ما از هوش مصنوعی (شبکه عصبی ، الگوریتم ژنتیک و الگوریتم فازی) برای مدل کردن سیستم ( مساله سری زمانی ) استفاده می کنیم و این روشها را با یکدیگر مقایسه می نماییم .

مساله سری زمانی به عنوان یک مساله شناسایی سیستم در نظر گرفته شده است ، درحالیکه ورودی های سیستم مقادیر گذشته سری زمانی اند و خروجی مطلوب مقدار آینده آن می باشد . شبکه های عصبی backpropagation ، RBF (Radial Basis Functions) و PNN (polynomial neural network ) مورد بررسی قرار گرفته و سپس مدل پیشرفته PNN یعنی GBSON (genetic-based self organizing network) که تلفیقی از الگوریتم ژنتیک و شبکه عصبی PNN می باشد مورد بررسی قرار گرفته است . در نهایت نیز یک مدل خاص از شبکه فازی برای تحلیل مورد استفاده قرار گرفته است .

مقدمه :

یک سری زمانی مجموعه ای از مشاهدات X(t) می باشد که هر کدام در یک زمان معین استخراج شده اند . یک سری زمانی گسسته ، سری ای است که در زمانهای گسسته نمونه برداری شده است . سری زمانی پیوسته سری ای است که با ثبت پیوسته مشاهدات در یک دوره زمانی بدست آمده است.[1]

مثالی از یک سری زمانی در شکل (1) آمده است .

برای تحلیل داده های سری زمانی باید آن را به اجزایش تجزیه کرد . هر جزء به عنوان یک فاکتور که در سری زمانی تاثیرگذار است در نظر گرفته میشود . سه جزء اصلی سری زمانی مشخص شده اند .

1- Trend (تمایل) : گرایش به افزایش یا کاهش در دراز مدت را در یک سری زمانی نمایش می دهد .

2- Seasonality (فصلی) : رفتار تناوبی یک سری زمانی را در یک دوره معین از زمان نمایش می دهد .

3- Fluctuation (نوسان، تغییر) جزء بی قاعده سری زمانی می باشد . [1]

شکل (1) سری زمانی مکی گلاس

کاربردهای پیشگویی سری های زمانی را می توان در حوزه های تجارت ، اقتصاد ، کنترل موجودی ، پیش بینی وضع هوا ، پردازش سیگنال و کنترل و سایر زمینه ها پیدا کرد .

در این مقاله هوش مصنوعی برای حل مساله سری زمانی بکار گرفته شده است . مساله سری زمانی به عنوان یک مساله شناسایی سیستم در نظر گرفته شده است ، درحالیکه ورودی های سیستم مقادیر گذشته سری زمانی اند و خروجی مطلوب مقدار آینده آن می باشد . یک روش برای حل اینگونه مسائل استفاده از شبکه عصبی می باشد . از میان شبکه های عصبی شبکه های عصبی TDNN و PNN بیشتر برای تشخیص سری زمانی بکار می روند در این مقاله شبکه های عصبی backpropagation ، RBF و PNN مورد بررسی قرار گرفته است و سپس برای رفع اشکالات آن از الگوریتم ژنتیک استفاده شده است و در نهایت نیز یک سیستم فازی برای حل مساله بکار گرفته شده است .

  1. سری زمانی مکی گلاس:

سری زمانی مکی گلاس با استفاده از فرمول زیر بدست می آید :

Time series prediction using neural network

این سری برای T>17 حالت آشوبگونه و اغتشاشی از خود نشان می دهد و بنابراین در مسایل پیش بینی سری های زمانی چون پیش بینی نقاط بعدی این سری مشکل می شود از این سری جهت آزمایش استفاده می شود .[2] نمونه این سری برای نقطه شروع 8/0 X=و 25T= به صورت زیر می شود :

Time series prediction using neural network

2-شبکه های عصبی :

ایده شبکه عصبی در سال 1943 توسط مک کلوچ و پیتز ارایه شد. این دو نفر مدل ساده نرون را ایجاد کردند که این نرون قابلیت ساخت دروازه های AND و OR را داشت . سپس این مدل توسط دانشمندان دیگر توسعه پیدا کرد و از حالت تک نرون به شبکه ای از نرون ها توسعه یافت . برای استفاده از یک شبکه عصبی باید ابتدا شبکه را تعلیم داد. در سالهای گذشته ساختار های متنوعی برای شبکه های عصبی ارایه شده است. که از مهمترین آنها می توان به شبکه های بازگشتی و feedforward اشاره کرد . برای این شبکه ها الگوریتم های مختلف تعلیم ارایه شده اند که می توان آنها را به دو دسته نظارتی و بدون ناظر تقسیم کرد . از جمله کاربردهای شبکه های عصبی می توان به شناخت الگوهای مرئی (صوت،چهره،تصویر)، تقریب توابع، طبقه بندی، کنترل کننده ها و پیش بینی (وضع هوا،سریهای زمانی و…) اشاره کرد. در بخش های بعد از شبکه های عصبی گوناگون برای پیش بینی سری زمانی استفاده شده است که ابتدا نوع شبکه ، ساختار آن و سپس الگوریتم تعلیم شبکه توضیح داده شده است .

2-1 شبکه عصبی RBF :

شبکه های عصبی بطور وسیعی در مسایل سری های زمانی و پیش بینی بکار می روند ساختارها و مدلهای گوناگونی از این شبکه ها برای این گونه مسایل طراحی شده اند یکی از این شبکه های عصبی شبکه عصبی توابع پایه ای شعاعی می باشد . این شبکه یک لایه پنهان دارد که متشکل از نرونهایی است که تابع فعالیت آنها یک تابع گوسین ، چندجمله ای ، مکعبات و … می تواند باشد . پارامترهای شبکه عصبی RBF وزنهای بین لایه پنهان و لایه خروجی ، مراکز توابع فعالیت وشعاع توابع فعالیت می باشند . این شبکه برای مدلسازی سیستم سری زمتنی بکار برده شده است و نتایج آن در قسمت آزمایشات با بقیه شبکه ها مقایسه شده است .

الگوریتم تعلیم شبکه عصبی RBF :

ورودی الگوریتم:

حداکثر خطای مورد نظر ، الگوهای تعلیم و شعاع نرون های لایه مخفی

مراحل:

  1. l=1 تعداد یک نرون در لایه مخفی
  2. یک الگوی تعلیم را بطور تصادفی انتخاب کرده و آن را به عنوان مرکز نرون l قرار می دهیم .
  3. محاسبه وزنهای لایه دومTime series prediction using neural network
  4. خطای شبکه به تمامی الگوهای تعلیم را بدست می آوریمTime series prediction using neural network
  5. اگر خطا از خطای مورد نظر کمتر است توقف کن
  6. در غیر اینصورت l=l+1 به گام 2 برگرد .

مزایا و معایب شبکه RBF :

  • اشکالات RBF :
  1. این شبکه برای الگوهای تعلیم خطای نسبتا کمی دارد ولی برای الگوهای تست خطا زیاد می شود
  2. شعاع r باید ابتدا تعیین شود که با تغییر r وزنها خیلی تغییر می کنند .
  3. 2هر چه تعداد نرون ها بیشتر شود تعمیم پذیری شبکه (جواب آن به داده های تست) بدتر می شود .

2-2- شبکه عصبی Backpropagation :

  • الگوریتم شبکه عصبی Backpropagation :

الگوریتم دسته ای :

  1. وزنهای اولیه شبکه را در مقادیر تصادفی کوچک قرار می دهیم .
  2. برای تمام الگوهای تعلیم P=1,2,…,p مراحل انتشار مستقیم و معکوس زیر را تکرار کن

2-1) انتشار مستقیم : ورودی Xp را به شبکه اعمال و تمام خروجیها و خالص ورودیها (net)را حساب کن .

2-2) انتشار معکوس : ابتدا خطا را حساب می کنیم سپس را حساب مي كنيم . يعنی از لايه خروجي شروع كرده و اين محاسبات را انجام مي دهيم و براي لايه های دیگر نيز محاسبات را به صورت زير انجام مي دهیم:

2-3) مقادیر را ذخیره می کنیم .

2-4) حساسیت نسبت به تک تک الگوها را حساب کرده و همه را با هم جمع می کنیم

  1. وزنهای جدید را محاسبه می کنیم :
  2. به قدم 2 برگرد تا شرط پایان برقرار شود. (ساده ترین شرط توقف تعداد تکرار الگوریتم یا رسیدن الگوریتم به یکخطای مطلوب است )

الگوریتم تکراری نیز دقیقا مانند این روش است با این تفاوت که بعد از بدست آوردن خطای یک الگو وزنها را عوض می کنیم .البته الگوریتم تکراری همگرایی اش سریعتر است .

1476-4598-4-29-1-lسیستم فازی :

نظريه فازی برای اولين بار در سال 1960 توسط استاد ايرانی دانشگاه کاليفرنيا، پروفسور لطفی زاده مطرح شد. درکنترل فازی برای هر پديده يک تابع در نظر مي گيرند که ميزان تعلق آن را به يک مجموعه بيان مي کند .

در این شيوه سه مرحله وجود دارد. ابتدا مرحله فازی سازی است که مرحله تعريف مجموعه های فازی برای متغيرهای ورودی وخروجی است.برای تعريف اين مجموعه های فازی، بايد دانش اوليه ای از دامنه تعريف هر کدام از اين متغيرها داشته باشيم .

در مرحله بعد يا مرحله استنتاج،تعدادی قاعده فازی بوجود مي آوريم و در مرحله آخر که مرحله غيرفازی سازی نام دارد،با توجه به مقادیری که در مرحله استنتاج بدست آمده است مقداری حقيقی برای خروجی بدست مي آوريم

الگوریتم تعليم (روش جدول جستجو):

  1. مجموعه های فازیی را تعیین می کنیم که زوج های ورودی – خروجی را پوشش دهد
  2. از روی هر زوج ورودی – خروجی یک قاعده تولید می کنیم.
  3. یک درجه به هر قاعده تولید شده در گام دوم نسبت می دهیم و از بین قواعد متضاد قاعده ای که بالاترین درجه را دارا می باشد نگه می داریم. که به این ترتیب تعداد قواعد کاهش می یابد.
  4. پایگاه قواعد فازی را ایجاد می کنیم.
  5. ساخت سیستم فازی بر اساس پایگاه قواعد.

این طراحی بر اساس سیستم فازی با موتور استنتاج ضرب، فازی ساز منفرد و غیر فازی ساز میانگین مراکز به صورت زیر خواهد بود:

Time series prediction using neural network

اشکالات سیستم (روش جدول جستجو):

  1. در این روش توابع تعلق ثابت اند و به زوج های ورودی خروجی بستگی ندارند . بدین معنی که توابع تعلق با توجه به زوج های ورودی خروجی به طور بهینه ای به دست نمی آیند.
  2. پایگاه قواعد فازی در این روش ممکن است کامل نباشد که در نتیجه سیستم پاسخ خوبی به داده های تست نمی دهد.

GMDH (Group Method Data Handling) :

در اواخر دهه 1960 توسط Ivakhnenko به عنوان وسیله ای برای شناسایی ارتباطات غیرخطی بین متغیرهای ورودی و خروجی ارائه شد . الگوریتم GMDH یک ساختار بهینه از مدل را توسط نسلهای موفق توصیفات جزئی (PD) داده که به عنوان چندجمله ای های درجه دوم با دو ورودی در نظر گرفته شده اند را ارائه می کند . البته این الگوریتم مشکلاتی نیز دارد که بعدها این مشکلات توسط ارائه شبکه عصبی چند جمله ای تا حدودی مرتفع شده است .

تابع انتقال هرPD از فرمول زیر بدست می آید :

Time series prediction using neural network

اشکالات GMDH :

1- گرایش به ایجاد چندجمله ای کاملا پیچیده برای سیستمهای نسبتا ساده دارد .

2- به علت ساختار عمومی GMDH ( چند جمله ای دو متغیره درجه دوم ) ، گرایش به ایجاد یک شبکه (مدل) کاملا پیچیده وقتی که برای سیستمهای غیرخطی مراتب بالاتر بکار می رود دارد .

3- برای سیستمهای کمتر از سه ورودی ساختار مناسبی ندارد .

راه حل :

استفاده از شبکه عصبی چندجمله ای در این روش

GMDH-network

2-3- شبکه عصبی PNN (Polynomial neural network):

 

شبکه عصبی چندجمله ای الگوریتم جدیدی است که قابلیت انعطاف بسیار زیادی دارد به طوریکه هر گره در این شبکه هم از لحاظ تعداد ورودی (حداکثر سه ورودی) و هم از لحاظ تابع انتقال آن (حداکثر سه ورودی) می تواند متفاوت باشد .

مرتبه چندجمله ای ها در هر گره شبکه می تواند متفاوت (خطی ، درجه دوم ومکعب)باشد .

معماری این شبکه از ابتدا ثابت نیست (هم ساختار و هم پارامترهای شبکه ) اما کاملا بهینه می شود ، بطور مثال تعداد لایه ها می تواند با اضافه شدن لایه جدید اگر لازم باشد افزایش یابد .

الگوریتمهای آموزش شبکه‏های چندجمله‏ای

    • Group Method of Data Handling (GMDH)
    • Polynomial Network Training Routine (PNETTR)
    • Algorithm for Synthesis of Polynomial Network (ASPN)
    • همگی توسط A.G.Ivakhnenko معرفی شده‏اند.

 

ساختار کلی شبکه عصبی PNN به صورت شماتیک :

Time series prediction using neural network

الگوریتم تعلیم شبکه عصبیPNN :

ابتدا داده ها را به صورت زیر در می آوریم :

Time series prediction using neural network

برای محاسبه خروجی تخمینی ، ما برای هر جفت ورودی یک PD تشکیل می دهیم. تعداد PD ها بر اساس تعداد ورودیها می باشد .

مثلا اگر 4 ورودی داشته باشیم تعداد آنها 6 عدد می شود . پارامترهای مدل از روی حداقل خطای داده های آموزش بدست می آیند .

بعلاوه ما بهترین مدل را برای تشکیل لایه اول انتخاب می کنیم.

در نهایت ما PD های جدید را از روی متغیرهای میانه (Zm ها) که در تکرار جدید واقع شده اند ایجاد می کنیم .

سپس ما یک جفت از متغیرهای ورودی جدید را گرفته و عملیات را روی آن انجام می دهیم تا زمانیکه به معیار توقف برسیم .

هنگامی که لایه پایانی ساخته شد گره با بهترین عملکرد به عنوان خروجی در نظر گرفته می شود . بقیه گره های آن لایه حذف می شوند . این عملیات در لایه های قبلی تا لایه اول انجام می شود .

مراحل طراحی شبکه PNN :

  1. ابتدا متغیرهای ورودی را نرمالیزه کرده و مطابق جدول زیر آنها را تقسیم می کنیم .

X متغیرهای ورودی

خروجیY

x1 x2 ….. xN

Y

2-داده ها را به دو دسته تست و آموزش تقسیم می کنیم :

X

Y

x11 x21 ….. xN1

x12 x22 ….. xN2

.

.

.

x1ntr x2ntr ….. xNntr

Y1

Y2

.

.

.

Yntr

x11 x21 ….. xN1

x12 x22 ….. xN2

.

.

.

x1nte x2nte ….. xNnte

Y1

Y2

.

.

.

Ynte

n=nte+ntr داده آموزش برای تشکیل مدل PNN استفاده می شود . (محاسبه ثوابت هر PD گرههای موجود در هر لایه)سپس داده تست برای ارزیابی مدل PNN استفاده می شود .

  1. انتخاب یک ساختار برای PNN :

ساختار PNN براساس تعداد ورودی ها و مرتبه PD هر لایه انتخاب می شود . دو نوع ساختار PNN بنیادی و PNN اصلاح شده ، متمایز شده اند . هر کدام از آنها به دو دسته تقسیم می شوند : جدول زیر این دو ساختار را نشان می دهد .

  1. تعیین تعداد متغیرهای ورودی و مرتبه چندجمله ای تشکیل دهنده PD :

ما متغیرهای ورودی یک گره را ازبین N متغیر ورودی انتخاب می کنیم . تعداد کل PD های لایه حاضر بر طبق تعداد متغیرهای ورودی انتخاب شده از گره های لایه قبل تغییر می کند . در نتیجه تعداد گرهها برابر N!/(N-r)!r! می شود که r تعداد متغیرهای ورودی انتخاب شده است.

  1. تخمین ثوابت PD :

از رابطه زیر محاسبه می شوند :

Time series prediction using neural network

i شماره گره ، k شماره داده ntr تعداد داده آموزش ، n تعداد ورودیهای انتخابی ، m ماکسیمم درجه، n’ تعداد ثوابت تخمینی .

5- انتخاب PD با بهترین عملکرد :

هر PD توسط داده های آموزش و تست تخمین زده شده و ارزیابی می شود . سپس ما اینها را باهم مقایسه کرده و PD هایی را که عملکرد بهتری داشته اند جدا می کنیم . معمولا ما یک تعداد W از قبل تعیین شده از این PD ها را انتخاب می کنیم .

6-چک کردن شرط توقف :

الف- اگر خطای مرحله بعدی از خطای مرحله قبلی بزرگتر باشد .

ب-اگر تعداد لایه ها به تعداد از قبل تعیین شده توسط کاربر برسد.

7- انتخاب متغیرهای جدید برای لایه بعد.

خروجی های لایه قبل را به عنوان ورودی لایه جدید در نظر می گیریم.

مراحل 4 تا 7 را تا رسیدن به شرط توقف تکرار می کنیم .

Time series prediction using neural network

همان طور که دیده می شود ضرایب این شبکه به وسیله فرمول زیر به دست می آیند :

Time series prediction using neural network

این ضرایب به صورت زیر بدست می آیند ، همانطورکه می دانیم خروجی به صورت زیر بدست می آید:

حال برای محاسبه ضرایب مراحل زیر طی می شوند :

Time series prediction using neural network

هر چند این فرمول ضرایبی را می دهد که خطا را کاهش می دهد اما این خطا می نیمم خطا نیست و با تغییر ضرایب می توان به جوابهای بهتری دست یافت .

اما این ضرایب را چگونه می توان از میان انبوه اعداد یافت؟

یکی از راه حلها استفاده از الگوریتم ژنتیک برای بدست آوردن بهترین ضرایب جهت کاهش خطا می باشد .

حال به صورت زیر از الگوریتم ژنتیک برای بهبود ضرایب در شبکه عصبی PNN استفاده شده است :

الف – ورودیهای لایه اول را همان متغیرهای مستقل داده ورودی در نظر می گیریم .

ب – خروجی های هر لایه بهترین گره های بدست آمده توسط الگوریتم ژنتیک هستند که اینها خود ورودی های لایه بعد می شوند .

ج – ضرایب هر گره را به صورت جمعیت اولیه برای الگوریتم ژنتیک در نظر می گیریم و الگوریتم ژنتیک را برای یافتن بهترین ضرایب بکار می بریم . پس از بدست آوردن بهترین ضرایب این گره های با بهترین ضرایب را از بقیه جدا می کنیم و ورودی را به اینها می دهیم تا خروجی گره های لایه اول بدست آید . این خروجی ها را به عنوان ورودی لایه بعدی در نظر می گیریم . این کار را تا زمانی که به خطای مورد نظر برسیم ادامه می دهیم .

روشی که در این تحقیق استفاده شده است به صورت زیر می باشد :

ابتدا سه گره را در نظر می گیریم . سپس برای هر گره به ترتیب 20 جمعیت در نظر می گیریم . هر فرد به صورت زیر در نظر گرفته می شود .

هر فرد را به عنوان یک رشته 128 بیتی در نظر می گیریم که هر 32 بیت آن مربوط به یک ضریب می باشد . در اینجا 4 ضریب در نظر گرفته شده است .

بیت اول از هر 32 بیت به عنوان بیت علامت در نظر گرفته می شود (+ یا – ) . 2 بیت بعدی برای قسمت صحیح ضریب و مابقی بیت ها برای قسمت اعشاری در نظر گرفته می شوند .

Time series prediction using neural network

 

حال ابتدا با استفاده از روش تورنمنت (q=6) 6 فرد را انتخاب کرده دو نفر را به عنوان والد انتخاب می کنیم .

سپس عمل تقاطع را با ضریب Pc=0.95 و بعد آن عمل جهش را با ضریب Pm=.05 انجام می دهیم ، تا دو فرزند جدید ایجاد شوند . این کار را تا زمان تولید P=20 فرزند تکرار می کنیم.

حال با استفاده از روش Elitism (نخبه گرایی) بهترین فرد نسل قبل را نیز به نسل جدید منتقل می کنیم . این کار را به تعداد نسل های لازم (3500 نسل ) تکرار می کنیم .

سپس بهترین فرد نسل آخر را به عنوان ضرایب برای این گره انتخاب میکنیم و به گره بعد می رویم .

بعد از این که لایه اول تمام شد اگر به شرط توقف رسیدیم می ایستیم در غیر اینصورت به لایه بعد می رویم .

خروجی های لایه اول را به عنوان ورودی لایه بعد در نظر گرفته و همین مراحل را برای لایه بعد تا زمان رسیدن به پاسخ مطلوب ادامه می دهیم .

در کلیه محاسبات برای شبکه ها تعداد ورودی ها را 3 در نظر گرفته ایم . یعنی جمله چهارم را بر اساس سه جمله قبل مشخص می کنیم .

– نتایج بدست آمده از شبکه BP :

در اینجا از شبکه عصبی BP با ساختار 3L 15N 1L استفاده شده است که نتایج زیر بدست آمده است.

نمودار تعداد تکرار اموزش شبکه وخطای بدست امده برای داده ها:

نمودار جواب مطلوب و جواب شبکه به داده آموزش که به صورت زیر می باشد ولی به علت خطای بسیار ناچیز با چشم قابل مشاهده نیست .

خطای بدست آمده از داده آموزش:

نمودار جواب مطلوب و جواب شبکه به داده تست که به صورت زیر می باشد ولی به علت خطای بسیار ناچیز با چشم قابل مشاهده نیست .

Time series prediction using neural network

 

خطای داده تست:

Time series prediction using neural network

2- نتایج بدست آمده از شبکه RBF:

شکل زیر داده آموزش مطلوب و پاسخ شبکه به داده های ورودی را نشان می دهد.

Time series prediction using neural network

 

شکل زیر خطای شبکه نسبت به داده آموزش را نشان می دهد :

Time series prediction using neural network

نمودار جواب مطلوب و جواب شبکه RBF به داده تست که به صورت زیر می باشد :

Time series prediction using neural network

خطای داده تست:

3- نتایج بدست آمده از سیستم فازی (جدول جستجو):

شکل زیر داده آموزش مطلوب و پاسخ سیستم به داده های ورودی را نشان می دهد:

به علت خطای زیاد سیستم فازی با چشم نیز این خطا قابل مشاهده است .

شکل زیر خطای سیستم نسبت به داده آموزش را نشان می دهد : 1246/0 = Etr

نمودار جواب مطلوب و جواب سیستم فازی به داده تست که به صورت زیر می باشد :

Time series prediction using neural network

هر چند که در داده های تست سیستم فازی بهتر عمل کرده است اما باز هم خطا نسبت به شبکه عصبی زیاد است.

خطای داده تست: 032/0Ete =

Time series prediction using neural network

4- نتایج بدست آمده از شبکه عصبی PNN :

شکل زیر داده آموزش مطلوب و پاسخ سیستم به داده های ورودی را نشان می دهد:

Time series prediction using neural network

شکل زیر خطای شبکه نسبت به داده آموزش را نشان می دهد :Time series prediction using neural network

نمودار جواب مطلوب و جواب شبکه PNN به داده تست که به صورت زیر می باشد :

 

شکل زیر خطای شبکه نسبت به داده آموزش را نشان می دهد :

5- نتایج بدست آمده از شبکه عصبی GBSON:

شکل اصلی و شکل پیش بینی شده داده آموزش :

شکل زیر خطای شبکه نسبت به داده آموزش را نشان می دهد :

نمودار جواب مطلوب و جواب شبکه GBSON به داده تست :

خطای داده تست :

مقایسه عملکرد شبکه ها در پیش بینی سری زمانی

شبکه عصبی ژنتیک GBSON شبکه عصبی PNN شبکه عصبی RBF شبکه عصبی BP سیستم فازی
0023/0 00245/0 0011/0 0023/0 1246/0 خطای داده آموزش
0016/0 00186/0 0031/0 0017/0 032/0 خطای داده تست
00195/0 0022/0 0022/0 002/0 0639/0 میانگین خطا

نتیجه گیری :

از لحاظ سرعت شبکه عصبی PNN بهترین سرعت را چه در مرحله آموزش و چه در مرحله تست دارا می باشد .

بهترین پاسخ برای داده آموزش را شبکه عصبی RBF دارا می باشد و بدترین پاسخ مربوط به سیستم فازی می باشد .

بهترین پاسخ برای داده تست را شبکه عصبی- ژنتیک GBSON دارا می باشد و بدترین پاسخ مربوط به سیستم فازی می باشد . البته سرعت آموزش شبکه GBSON بسیار پایین است .

از لحاظ میانگین خطا نسبت به داده تست و آموزش نیز بهترین پاسخ مربوط به شبکه GBSON می باشد .

MFILE ها :

1- برنامه نوشته شده برای شبکه BP در مطلب :

clc

clear

load(‘E:\shenasaey projeh\x’);load(‘E:\shenasaey projeh\xtr’);load(‘E:\shenasaey projeh\xte’)

load(‘E:\shenasaey projeh\ytr’);load(‘E:\shenasaey projeh\yte’);

xtr=xtr’;

ytr=ytr’;

xte=xte’;

net.trainParam.epochs = 500;

net.trainParam.goal = 0;

net.performFcn=’mse’;

net = newff(xtr,ytr,[15]);

net = train(net,xtr,ytr);

y = sim(net,xtr);

y=[x(1:3) y];

figure(1)

plot(y)

hold on

plot(x(1:300),’k’)

etr=(x(1:300)-y).^2;

figure(2)

plot(x(1:300)-y,’r’)

Etr=sqrt(sum(etr)/300)

figure(3)

yte=sim(net,xte);

plot(yte,’r’)

hold on

plot(x(301:600),’b’)

e=(x(301:600)-yte).^2;

figure(4)

plot(x(301:600)-yte,’k’)

E=sqrt(sum(e)/300)

2- برنامه نوشته شده برای شبکه RBF :

clc

clear

load(‘E:\shenasaey projeh\x’)

load(‘E:\shenasaey projeh\xte’)

load(‘E:\shenasaey projeh\xtr’)

load(‘E:\shenasaey projeh\yte’)

load(‘E:\shenasaey projeh\ytr’)

xtr=xtr’;

ytr=ytr’;

xte=xte’;

yte=yte’;

net=newrbe(xtr,ytr);

Y=[];

for i=1:297

Y =[Y sim(net,xtr(:,i))];

end

Y=[x(1,1:3) Y];

plot(x(1,1:300),’b’)

hold on

plot(Y,’k’)

figure(2)

plot(x(1,1:300)-Y,’b’)

Etr=sqrt(sum((x(1,1:300)-Y).^2)/300)

Y=[];

for i=1:300

Y =[Y sim(net,xte(:,i))];

end

figure(3)

plot(yte,’k’)

hold on

plot(Y,’b’)

Ete=sqrt(sum((yte-Y).^2)/300)

figure(4)

plot(yte-Y,’b’)

  1. برنامه نوشته شده برای سیستم فازی :

clc

clear

T=30;

x=[1.2];

delta=1;

for i=2:T

x(i)=x(i-1)/1.1;

end

for i=T+1:round(600/delta)

x(i)=(1/1.1)*(x(i-1)+0.2*x(i-T)/(1+x(i-T)^10));

end

max_x=max(x)+0.2;

min_x=min(x)-0.2;

n=4;

sin=15;

stp=(max_x-min_x)/(sin-1);

y=min_x:0.01:max_x;

cent=min_x:stp:max_x;

% for i=1:sin

% k=trimf(y,[cent(i)-stp cent(i) cent(i)+stp]);

% end

so=15;

stp_o=(max_x-min_x)/(so-1);

cent_o=min_x:stp_o:max_x;

% for i=1:so

% k=trimf(y,[cent_o(i)-stp_o cent_o(i) cent_o(i)+stp_o]);

% end

for k=1:296

for t=1:4

for i=1:sin

mu_a(i)=trimf(x(k+t-1),[cent(i)-stp cent(i) cent(i)+stp]);

end

rule(k,t)=find(mu_a==max(mu_a));

rule_d(k,t)=max(mu_a);

end

for i=1:so

mu_b(i)=trimf(x(k+4),[cent_o(i)-stp_o cent_o(i) cent_o(i)+stp_o]);

end

B(k,1)=find(mu_b==max(mu_b));

B_d(k,1)=max(mu_b);

end

D=prod(rule_d’)’.*B_d;

for i=296:-1:1

for j=i-1:-1:1

if rule(i,:)==rule(j,:)

if D(j)<D(i)

D(j)=-inf;

else

D(i)=D(j);

B(i)=B(j);

D(j)=-inf;

end

end

end

end

for i=296:-1:1

if D(i)==-inf

D(i)=[];

rule(i,:)=[];

B(i)=[];

end

end

l_rule=length(rule);

hold on;

for u=1:596

suming_s=0;

suming_m=0;

pro=ones(1,l_rule);

y=[x(u) x(u+1) x(u+2) x(u+3)];

for k=1:l_rule

for t=1:4

pro(k)=pro(k)*trimf(y(t),[cent(rule(k,t))-stp cent(rule(k,t)) cent(rule(k,t))+stp]);

end

suming_s=suming_s+pro(k)*cent_o(B(k));

suming_m=suming_m+pro(k);

end

f(u+4)=suming_s/suming_m;

end

plot(x(1,1:300),’b’)

hold on

plot(f(1,1:300),’r’)

figure(2)

plot(x(1,1:300)-f(1,1:300),’b’)

Etr=sqrt(sum((x(1,1:300)-f(1:300)).^2)/300)

figure(3)

plot(x(1,301:600),’r’)

hold on

plot(f(1,301:600),’b’)

Ete=sqrt(sum((x(1,301:600)-f(1,301:600)).^2)/300)

figure(4)

plot(x(1,301:600)-f(1,301:600),’r’)

  1. برنامه نوشته شده برای شبکه PNN :

clc

clear

load(‘E:\shenasaey projeh\x’);load(‘E:\shenasaey projeh\xte’);load(‘E:\shenasaey projeh\yte’)

load(‘E:\shenasaey projeh\ytr’);load(‘E:\shenasaey projeh\xtr’)

c=[];

[sx,lx]=size(xtr);

y_col=[ytr;yte];

zte=[]; % khoroji laye 1

ztr=[];

xpd=cell(1,3);

ff=0;preE=10^10;

while 1

ff=ff+1

t=0;

for i=1:lx-1

for j=i+1:lx

t=t+1;

xpd{t}=Xpd([xtr(:,i) xtr(:,j)]);

end

end

for k=1:3

c(k,:)=inv(xpd{k}’*xpd{k})*xpd{k}’*ytr;

end

for kk=1:297

ztr=[ztr quad(xtr(kk,:),c)];

end

for kk=1:300

zte=[zte quad(xte(kk,:),c)];

end

z=[ztr’;zte’];

E_col=[];

for k=1:3

e=0;

for jj=1:597

e=e+(y_col(jj,1)-z(jj,k))^2;

end

u=e/597;

E_col=[E_col sqrt(u)];

end

xtr=ztr’;xte=zte’;

zz=z;

z=[];

if min(preE)<=min(E_col)

break

end

preE=E_col;

ztr1=ztr;zte1=zte;

ztr=[];zte=[];

[sx,lx]=size(xtr);

xpd=[];

end

Etr=[];

for k=1:3

e=0;

for jj=1:297

e=e+(y_col(jj,1)-zz(jj,k))^2;

end

u=e/297;

Etr=[Etr sqrt(u)]

end

Ete=[];

for k=1:3

e=0;

for jj=298:597

e=e+(y_col(jj,1)-zz(jj,k))^2;

end

u=e/300;

Ete=[Ete sqrt(u)]

end

y_col=[x(1:3)’; y_col];

figure(1)

plot(y_col(1:300,1))

hold on

ztr1=ztr1′;

zte1=zte1′;

ztr=[ztr1(:,1);zte1(:,1)];

ztr=[x(1:3)’;ztr];

plot(ztr(1:300,1),’r’)

s=(y_col-ztr);

figure(2)

plot(s)

figure(3)

plot(y_col(301:600,1))

hold on

plot(ztr(301:600,1),’r’)

s=(y_col(301:600,1)-ztr(301:600,1));

figure(4)

plot(s)

5- M-FILE مربوط به آموزش شبکه عصبی – ژنتیک GBSON :

clc

clear

load(‘C:\Documents and Settings\PC\Desktop\New Folder\Dtr’)

load(‘C:\Documents and Settings\PC\Desktop\New Folder\Dte’)

x=Dtr(:,1:3); % dadeye vorodi

xte=Dte(:,1:3);yte=Dte(:,4);

[sx,lx]=size(x);ytr=Dtr(:,4); % khoroji matlob

y_col=[ytr;yte];cs=[];c_f=[];P=20;l=128;

yyte=[]; % khoroji laye

Pc=.95;Pm=.1;q=6;

a1=cell(1,3);

c=cell(1,3);

c_laye=[];

xu=cell(1,3);

for laye=1:15

E=[];

ztr_col=[];ztr=[];

aa=[];a_new=[];

n=300;

for i=1:3

a1{i}=randint(P,l); %matrice valedin

end

for node=1:3

k=node;

a=a1{node};

for g=1:3500

E=[];

i=1;

while i<5

for ii=1:P

c_f(ii,i)=2.^(-1:-1:-29)*a(ii,32*(5-i)-28:(32*(5-i)))’;

cs(ii,i)=bi2de(a(ii,32*(5-i)-30:32*(5-i)-29),’left-msb’)+c_f(ii,i); %halgheye ijade sabetha

if a(ii,32*(5-i)-31)==1

cs(ii,i)=-cs(ii,i);

end

end

i=i+1;

end

for i=1:P %halgheye be dast avardane khoroji

xpd=Xpd(x,k);

for kk=1:300

ztr=[ztr;sum(xpd(kk,:).*cs(i,:))];

end

ztr_col(:,i)=ztr;

ztr=[];

end

syy=0;

for i=1:P %halgheye be dast avardane khata

for j=1:300

syy=syy+(ztr_col(j,i)-ytr(j,1))^2;

end

E=[E;sqrt((1/n)*(syy))];

syy=0;

end

if rem(g,50)==0

g

laye

E(1:3,:)

cs(1:3,:)

end

e=E;

fd=[];

ii=0;

while ii<3

ii=ii+1;

fd=find(e==min(e));

if length(fd)>1

for j=2:length(fd)

e(fd(j,1),1)=10^5*e(fd(j,1),1);

end

end

a_e(ii,:)=a(fd(1,1),:);

[sa,la]=size(a_e);

e(fd(1,1),:)=10^5*e(fd(1,1),:);

cr=[];

if sa>2

cr=[];

for k1=1:sa

for k2=k1+1:sa

if k1~=k2

k3=1;

cu=0;

while k3<128

if a_e(k1,k3)==a_e(k2,k3)

cu=cu+1;

end

k3=k3+32;

end

if cu==4

cr=[cr k2];

end

end

end

end

end

if cr~=0

if length(cr)==3

a_e(2:3,:)=[];

ii=ii-2;

end

if length(cr)==1

a_e(cr(1),:)=[];

ii=ii-1;

end

end

fd=[];

end

for i=1:P %halgheye fitnes

f(i,1)=1/(E(i,1));

end

sum_f=sum(f(1,:));

ss=0;

m=0;

for tt=1:P/2 %halgheye tolide farzandan

o=[];

while m<2 %halgheye entekhabe 2 valed

ff=[];

for i=1:q

rd=fix(P*rand)+1;

o=[o;a(rd,:)];

ff=[ff;f(rd,1)];

end

sum_ff=sum(ff(1,:));

r=sum_ff*rand;

for i=1:q

if ss<r

ss=ss+ff(i);

end

if ss>r

aa=[aa;o(i,:)];

ss=0;

r=sum_ff*rand;

[m,vv]=size(aa);

end

if m==2

ss=0;

break

end

end

end

a_old=aa;

ri=fix((l-2)*rand)+1;

rr=rand;

if rr<Pc %halgheye entekhabe cross over

for j=ri:l

a_old(1,j)=aa(2,j);

a_old(2,j)=aa(1,j);

end

ii=i+1;

else

ii=l;

end

for j=1:2 %halgheye mutation

for i=ii:l

zz=rand;

if zz<Pm

if a_old(j,i)==1

a_old(j,i)=0;

else

a_old(j,i)=1;

end

end

end

end

a_new=[a_new;a_old];

a_old=[];

aa=[];

m=0;

end %payane GA

a=a_new;

[tss,tse]=size(a_e);

for i=1:tss % halgheye jaygozinie behtarinhaye nasle ghabl ba nasle bad

if i<4

a(i,:)=a_e(i,:);

else

break

end

end

a_e=[];

a_new=[];

ztr=[];

end

E=[];

for i=1:P % halgheye mohasebeye khorojiha

xpd=Xpd(x,k);

for kk=1:300

ztr=[ztr;sum(xpd(kk,:).*cs(i,:))];

end

ztr_col(:,i)=ztr;

ztr=[];

end

syy=0;

for i=1:P % halgheye mohasebeye khataye khorojiha

for j=1:300

syy=syy+(ztr_col(j,i)-ytr(j,1))^2;

end

E=[E;sqrt((1/n)*(syy))];

syy=0;

end

fd=find(E==min(E));

ztr(:,1)=ztr_col(:,fd(1,1));

E=[];

xu{node}=ztr;

c{node}=cs(fd(1,1),:);

ztr=[];

ztr_col=[];

end

for im=1:3

syy=0;

xdh=xu{im};

for j=1:300

syy=syy+(xdh(j,1)-ytr(j,1))^2;

end

E=[E;sqrt((1/n)*(syy))];

if E<.002

break

end

syy=0;

end

for i=1:3

x(:,i)=xu{i};

c_laye(i,:)=c{i};

end

c_col{laye}=c_laye

end

 

M-FILE مربوط به تست شبکه عصبی – ژنتیک GBSON در متلب :

clc

clear

load(‘C:\Documents and Settings\PC\Desktop\New Folder\Dte’);

load(‘i:\arshad\Genetic\genetic project\Gbson\tornoment\c_colbitavan62tr8001’);

xte=Dte(:,1:3);yte=Dte(:,4);E=[];

zte=[];Ecol=[];

for ij=1:4

css=c_col{1,ij};

zte_col=[];

for i=1:3 %halgheye be dast avardane khoroji

xpd=Xpd(xte,i);

for kk=1:297

zte=[zte;sum(xpd(kk,:).*css(i,:))];

end

zte_col(:,i)=zte;

zte=[];

end

syy=0;

for i=1:3 %halgheye be dast avardane khata

for j=1:297

syy=syy+(zte_col(j,i)-yte(j,1))^2;

end

E=[E;sqrt((1/297)*(syy))];

syy=0;

end

xte=zte_col;

Ecol=[Ecol E];

E=[];

end

zte1=xte(:,3);

figure(1)

plot(yte)

hold on

plot(zte1,’r’)

s=(yte-zte1);

figure(2)

plot(s)

مراجع:

1.A Design of EA-based Self-Organizing Polynomial Neural Networks using

Evolutionary Algorithm for Nonlinear System Modeling

Dong-Won Kim and Gwi-Tae Park

2. Implementation of Artificial Intelligence in the Time Series Prediction Problem, L. EKONOMOU, S.SP. PAPPAS, September 22-24, 2006

3. Time Series PredictionUsingEvolving Polynomial Neural Networks, Amalia Foka, 1999

تهیه این آموزش متلب با زحمت زیادی صورت گرفته لطفا در فضای مجازی نشر ندهید