Printable View
شما تازه توی این برنامه تعداد numplayersرو مشخص کردین. توی برنامه های قبلی که استفاده می شد از کجا عددش رو میگرفت؟ جایی وارد کرده بودین؟
طبق شکلهایی که توی مقاله هست باید واسه هر supplier یه راندمان وجود داره. یعنی یه عدد تقسیم بر عدد دیگه. اینجا من اگه صورت رو تقسیم بر مخرج کنم فقط یه عدد بدست میاد. اینطوری که نمیتونم پلات کنم!
ضمن اینکه صورت کسر، نقاط تعادل نش هست که عددش 0.66872 هست و اون اعداد 0.035774 بردار payoff هست ظاهرن!
و اینکه چرا همه ی عددها یکسان در میاد؟
برنامه را که گفتم اجرا کن حداقل ببین جوابش مثل جواب اون نمودار شده یا نه. اگر هم اشکالی وجود داشت معلوم می شد. چرا اینقدر چک و چونه می زنی؟ :))
من نمیتونم این را اجرا کنم وگرنه جوابش را الان با نمودار اینجا میذاشتم
مثلا numplayers در برنامه قبلی numsupplier بوده
برای نمودار هم که معلمه . این کار را 9 بار اجرا کن یعنی
یکبار numplayers را برابر 2 قرار بده
یکی برای 3
یکی برای 4
یکی برای 5
...
یکی برای 10
در نتیجه 9 تا عدد داریم که میشه باهاش نمودار کشید البته با یک حلقه ساده هم می تونیم انجام بدیم
عدد 0.66872 چیه؟ میزان تولیده
عدد 0.035774 چیه؟ پول به دست امده هست
در فرمول efficiency میزان پول یا سود به دست اومده تقسیم میشه نه میزان تولید
اینکه یکسان به دست میاد هم علتش اینه که ما در ارائه ورودی ها برای هیچ تولید کننده ای ارجحیتی نسبت به بقیه قائل نشدیم که انتظار داشته باشیم میزان تولید یا سود او با دیگران متفاوت باشه و هدف ما تنها به دست آورد میزان سودی بوده که از طریق تولید کننده ها در تعادل نش به دست میاد
ببخشید که انقدر سوال می پرسم. چون می خوام بدونم که چیکار می کنی :لبخندساده:
جوابش در زیر پیوست شد.
مگه این خودش از numplayer دو تا 10 اجرا نمیکنه؟ پس چرا دوباره یه حلقه براش تعریف کنیم؟
http://uplood.ir/DD4W
دیگه دارم مجبور می شم یه پک کامل متلب رو برای شما بصورت اشانتیون بفرستم تا روی سیستمت نصب کنی :لبخندساده:
این برنامه فقط برای 10 انجام شده حلقه ای در کار نیست اون حلقه هم اگه دقت کنی کامنت شده
برنامه قبلی یه اشتباه داشت باید مخرج برابر 4 و خورده ای بشه
function ret = fffu(q)
alpha=.3; beta=1.2; k=.38; D=6.3;
ret=-sum(alpha.*q.^2+k.*(D-sum(q)).*q);
end
clc; close all; clear;
numplayers=10
lb=zeros(1,numplayers);
ub=ones(1,numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5;,[],[],[],[],lb,ub);
disp('makhraj=')
disp(-fx1)
برای پلات کردنش حلقه for هم گذاشتم :
clc; close all; clear;
interfac
global Q;global k;global D;global alpha;
alpha=.3; beta=1.2; k=.38; D=6.3;
efficiency = zeros(1,10)
for numplayers=2:10
lb=zeros(1,numplayers);
ub=ones(1,numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5,[],[],[],[],lb,ub);
disp('makhraj=')
disp(-fx1)
disp('Soorat=')
disp(sum(c))
efficiency(numplayers) = sum(c)/(-fx1)
end
plot(2:10,efficiency(2:10).*(2:10)), grid
ylabel('Symmetric NE efficiency');
xlabel('Number of Suppliers in Microgrid')
جواب efficiency
efficiency =
Columns 1 through 7
0 0.0664 0.0443 0.0332 0.0266 0.0221 0.0190
Columns 8 through 10
0.0166 0.0148 0.0133
http://www.uplooder.net/img/image/66...c5dcb8/fig.png
پست بالایی رو تازه دیدم . اینطوری شد :
clc; close all; clear;
interfac
global Q;global k;global D;global alpha;
alpha=.3; beta=1.2; k=.38; D=6.3;
efficiency = zeros(1,10)
for numplayers=2:10
lb=zeros(1,numplayers);
ub=ones(1,numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5,[],[],[],[],lb,ub);
disp('makhraj=')
disp(-fx1)
disp('Soorat=')
disp(sum(c))
efficiency(numplayers) = sum(c)/(-fx1)
end
plot(2:10,efficiency(2:10).*(2:10)), grid
ylabel('Symmetric NE efficiency');
xlabel('Number of Suppliers in Microgrid')
جواب :
No active inequalities.
makhraj=
4.0937
Soorat=
0.3577
efficiency =
Columns 1 through 7
0 0.0925 0.0767 0.0762 0.0799 0.0824 0.0842
Columns 8 through 10
0.0855 0.0866 0.0874
نمودار :
http://www.uplooder.net/img/image/67...315ef/fig2.png
اصلا یه کاری کنیم
در ابتدا نمی خواهیم صورت را به دست بیاریم
تنها کاری که شما می کنید اینه که در فایل interfac بدون حلقه for فقط مقدار numplayers را از 2 تا 10 تغییر بده و برای هر بار که اجرا کردی یکی از اعدادی که در قسمت
with payoffs [ [0.035774
را کپی می کنی در یک جایی مثل نوت پد میذاری
بنابراین نوت پد میشه شامل 9 تا خط که خط اولش یک عدد مثل 0.035774 هست
لطفا این عدد ها را اینجا بذار
این مرحله اول که انجام شد و صورت را به دست آوردیم می ریم سراغ مرحله دوم یعنی به دست آوردن مخرج
clc; close all; clear;
makhraj = zeros(1,10)
for numplayers=2:10
lb=zeros(1,numplayers);
ub=ones(1,numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5,[],[],[],[],lb,ub);
disp('makhraj=')
disp(-fx1)
makhraj(numplayers) =(-fx1)
end
اگه بخوام اون payoff ها رو بدست بیارم کافیه interfac رو توی command window بنویسم.
ولی نمیدونم numplayer رو تو کجا گذاشتی که تغییر بدم. توی کدوم اسکریپت هست؟
خط 30 از فایل interfac که به این صورته
numplayers=10
را تغییر بده
payoff های دو تا 10 :
[1.934] , [1.554] , [1.174] , [0.794] , [0.414] , [0.077012] , [0.057655] , [0.044775] , [0.035774]
حالا اینها را که به دست آوردی بذار داخل یه بردار به نام a و ضرب در 2:10 کن
a.*(2:10)
تا مخرج به دست بیاد
صورت را هم که برنامه اش را گذاشتم
اینطوری شد :
clc; close all; clear;
a=[1.934 ,1.554 ,1.174 ,0.794 , 0.414 ,0.077012 ,0.057655 ,0.044775 ,0.035774];
soorat=a.*(2:10)
makhraj = zeros(1,10)
for numplayers=2:10
lb=zeros(1,numplayers);
ub=ones(1,numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5,[],[],[],[],lb,ub);
disp('makhraj=')
disp(-fx1)
makhraj(numplayers) =(-fx1)
efficiency(numplayers) = soorat(numplayers-1)/(-fx1)
end
plot(2:10,efficiency(2:10).*(2:10)), grid
ylabel('Symmetric NE Efficiency');
xlabel('Number of Suppliers in Microgrid')
نمودار :
http://www.uplooder.net/img/image/83...2f7ae/fig3.png
این الآن درسته؟
تو این برنامه که سایز صورت و مخرج یکی نیست!
مهندس وقتی ما payoff تغذیه کننده های مختلف رو بدست آوردیم دیگه نیازی نیست که ما بیایم مثلن برای تغذیه کننده های با تعداد 5 تا رو ضرب در 5 بکنیم.(این خودش payoff پنج تا تغذیه کننده هست!) منظورم این خطه :
soorat=a.*(2:10)
من الآن فقط خود بردار payoff رو بعنوان صورت در نظر گرفتم نمودارش اینطوری شد :
http://www.uplooder.net/img/image/84...c6eeb/fig5.png
معقول تره به نظرم!
اگه شکل (2) مقاله رو نگاه کنید دقیقن مثل نمودارهای upperbound و simulated شد.
علتش اینه که توی برنامه ای که گذاشتی دو بار ضربدر 2:10 کردی
a.*2:10
و
plot(2:10,efficiency(2:10).*(2:10)), grid
که لازمه موقع plot دیگه در 2:10 ضرب نشه
پس برنامه کلی به این صورت شد با نمودار دومی که گذاشتم؟
clc; close all; clear;
payoff=[1.934 ,1.554 ,1.174 ,0.794 , 0.414 ,0.077012 ,0.057655 ,0.044775 ,0.035774];
soorat=payoff.*(2:10)
% soorat=payoff;
makhraj = zeros(1,10) ;
efficiency = zeros(1,10);
for numplayers=2:10
lb=zeros(1,numplayers);
ub=ones(1,numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5,[],[],[],[],lb,ub);
disp('makhraj=')
disp(-fx1)
makhraj(numplayers) =(-fx1);
efficiency(numplayers) = soorat(numplayers-1)/makhraj(numplayers);
end
خب ما تا حالا از جدول 3 و 4 استفاده نکردیم. اینطوری باشه که این واسه همه ی مسائل که یکی در میاد!
سایز صورت و مخرج یکینیست. چیکار کنم؟
تا اینجا که پیش رفتیم مربوط به رقابت بین تولید کننده ها بود از این به بعد که در جدول 3و 4 اومده مربوط به همکاری بین اونهاست که فکر کنم از shapley استفاده کرده ( فکر کنم مقاله را به صورت کامل مطالعه نکردی)
مطالعه کردم ولی فکر میکردم که صورت کسر باید با نقطه تعادل نش بدست بیاد و مخرج کسر با گیم تئوری مثل shapley !
الآن باید چیکار کنیم؟
این قسمت رو که دیگه توی همون مقاله اولی بدست آوردیم دیگه. داریم کدش رو. با همون mat-tug
شاید با اون قبلی متفاوت باشه باید دقیقا مقاله را بخونیم چون خوم هم این قسمتش را دقیقا مطالعه نکردمنقل قول:
نوشته شده توسط coronaa
ما توی این برنامه تونستیم تازه مطابق شکل(4)، نمودار simulated efficiency رو رسم کنیم. نمودار upperbound و lowerbound هم از همین برنامه باید بدست بیاد؟
قبلن توی اولین کد یه همچین نموداری بدست آوردی که شبیه lowerbound هست. همینه یا اون کد اشتباس؟
http://www.uplooder.net/img/image/93...768b7/fig5.png
مهندس این کدی هست که قبلن واسه تقسیم بندی نمرات انجام داده بودی.نقل قول:
نوشته شده توسط rahnema1
قسمت rules اگه بخوام اینو روی کاغذ بصورت if و then بنویسم چطوریه؟
به این صورته ؟
if (input is Low) then (out1 is Low)
if (input is High) then (out2 is High)
این اشتباه بودنقل قول:
نوشته شده توسط coronaa
اینو یادم رفته فقط میدونم متلب میتونه همین قاعده را با یک دستور ( که یادم نیست) به صورت if then بنویسهنقل قول:
نوشته شده توسط coronaa
حالا بیا یه چیز من بهت یاد بدم :))نقل قول:
نوشته شده توسط rahnema1
این دوتا دستور باید زیر برنامه نوشته بشه. بعد بصورت GUI اینو بالا میاریم و داخلش قوانین نوشته شده.
writefis(a,'myfile.fis');
fuzzy('myfile.fis');
مهندس نمودار upperbound و lowerbound از روی همین برنامه بدست نمیآد؟
این کد زیر رو به برنامه اضافه کردم. نمودار relative error برای upperbound بدست اومد. شکل(4)
error(numplayers) = abs(ub(numplayers-1)-efficiency)/efficiency;
http://www.uplooder.net/img/image/16...upperbound.png
ولی lowerbound بدست نمیآد. فقط یه خط راست میشه. شاید ورودی lowerbound که همه رو صفر در نظر گرفتیم اشتباه باشه!
http://www.uplooder.net/img/image/60...d_relative.png
شما نظرتون چیه؟
upperbound و lowerbound در ص 425 و در Corollary 2 توضیح داده شده
ببینم چی کار می کنی
من الآن روی مقدار lowerbound مشکل دارم. فکر میکنم صفر نباید باشه!
نمودار error of upperbound درسته؟ تائید کردی؟
یه مشکل دیگه ای هم که دارم. توی برنامه اصلی سایز صورت و مخرج یکی نیست!
اینها اشتباهه ظاهرا پست من را نخوندی . باید از fmincon یا fminbnd استفاده کنینقل قول:
نوشته شده توسط coronaa
ببینیم چی کار می کنی
من از این برنامه استفاده کردم. شما اون برنامه اولی رو گفتی اشتباهه!
clc; close all; clear;
payoff=[1.934 ,1.554 ,1.174 ,0.794 , 0.414 ,0.077012 ,0.057655 ,0.044775 ,0.035774];
soorat=payoff.*(2:10)
makhraj = zeros(1,10) ;
efficiency = zeros(1,10);
error =zeros(1,10);
for numplayers=2:10
lb=zeros(1,numplayers);
ub=ones(1,numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5,[],[],[],[],lb,ub);
disp('makhraj=')
disp(-fx1)
makhraj(numplayers) =(-fx1);
efficiency(numplayers) = soorat(numplayers-1)/makhraj(numplayers);
error(numplayers) = abs(ub(numplayers-1)-efficiency)/efficiency;
end
plot(2:10,efficiency(2:10)), grid
ylabel('Symmetric NE Efficiency');
xlabel('Number of Suppliers in Microgrid')
axis([2 10 0 1.1])
% ----------------------------------------------
figure
plot(2:10,error(2:10)), grid
ylabel('relative error of the upperbound');
xlabel('Number of Suppliers in Microgrid')
% legend('cos_x','sin_x')
axis([2 10 0 1.1])
این lb,ub برابر با upperbound و lowerbound نیست که
خواهشا و لطفا به همون صفحه مراجعه کن
ای بابا ینی این مقادیر payoff و چیزهای دیگه رو بدست آوردیم همه اشتباه بود؟
پس نمودارا چرا درست در اومد؟؟؟
چرا اون صفحه را مطالعه نمی کنی؟؟نقل قول:
نوشته شده توسط coronaa
اصلا به حرفهام گوش نمی دی
اون چیزی که قبلا به دست آوردیم درست بود و lb و ub که 0 و 1 بود و ثابت بودند و نیازی به محاسبه اونها نیست
اما اینجا حدبالا و پایین برای effciency را می خواهیم به دست بیاریم نه برای supply
چشم! فقط منو نزن :))
من یه لحظه داشتم ناامید می شدم. فکر کردم همه اشتباهه.
باشه من میخونمش. ببینم میتونم کاری انجام بدم یا نه.
ممنونم.
سلام دوست عزیز :لبخندساده:
من یه فانکشن ساختم به اسم Lmax و یکی هم به اسم Jmax ، طبق روابط 16 و 17 . بصورت زیر :
function ret = Lmax(L)
global k;global D;global alpha;
ret=-((alpha-k)./m).*L.^2+k.*D.*L;
end
function ret = Jmax(J)
global k;global D;global alpha;
ret=-(((alpha./m-k).*J.^2)+K.*D.*J);
end
بعد اومدم توی برنامه اصلی حد بالا و پایین راندمان رو بصورت شرطی قرار دادم.(بالای روابط 16)
حالا دیگه نمیدونم مقادیر fminbnd و fmincon رو چطوری مشخص کنم.
لاندا از کجا باید بدست بیاد؟
clc; close all; clear;
global q;global k;global D;global alpha ;global lambda;
alpha=.3; beta=1.2; k=.38; D=6.3; m=10;
L=sum(q);
J=sum(q);
for numsupplier=2:10
if lambda > k*D/((m+1)*k-2*alpha)
lobnd=((m*k^2*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/(((alpha-k).*Lmax^2/m)+k*D*Lmax);
upbnd=((m*k^2*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/((alpha/m-k)*Jmax^2+k*m*D*Jmax);
end
if lambda < k*D/((m+1)*k-2*alpha)
lobnd= (m^2.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k).*Lmax.^2)+k*m*d*Lmax);
upbnd= (m.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k*m).*Jmax.^2)+k*m*d*Jmax);
end
[x1, fx1]=fmincon(@fff,Q,[],[],[],[],lobnd,upbnd);
for tekrar=1:400
for supplier=1:numsupplier
[x2, fx2]=fminbnd(@fff,0,1);
Q(supplier)=x2;
end
end
efficiency(numsupplier)=fx2/fx1;
end
plot(2:10,efficiency(2:10).*(2:10)), grid
ylabel('Symmetric NE efficiency');
xlabel('Number of Suppliers in Microgrid')
در ضمن، ظاهرن با روابط 40 و 41 و 42 میشه بازده رو بدست آورد. هم رابطه ی V_NE مشخصه و هم رابطه ی V_opt
اینجا دیگه انگار payoff هم نمیخواد!
چه دلیلی داشته که از این حلقه ها استفاده کردی؟
for tekrar=1:400
for supplier=1:numsupplier
آخه گفته بودی از fminbnd و fmincon باید استفاده کنی فکر کردم مثل اون برنامه اولی باید نوشته بشه!نقل قول:
نوشته شده توسط rahnema1
الآن این قسمت برنامه رو حذف کنم؟
for tekrar=1:400
for supplier=1:numsupplier
[x2, fx2]=fminbnd(@fff,0,1);
Q(supplier)=x2;
end
end
efficiency(numsupplier)=fx2/fx1;
بخشهایی را خوب پیش رفتی
اولا توضیح بده اون صفحه چیکار می خواهد بکنه و مراحلش را بگو
ثانیا در هلپ متلب در مورد دو تابع fmincon و fminbnd مطالعه کن و بگو هر کدوم برای چه کاری مناسبه و اینکه میشه یکی را به جای دیگری استفاده کرد؟
ما میخوایم حد بالا و پایین راندمان رو بدست بیاریم. یعنی مشخص کنیم که راندمانی که بدست میاد از این مقادیر نباید تجاوز کنه.
توی این کرانهای بالا و پایین، دو متغیر Lmax و Jmax هم دارم که باید مشخص بشن. یعنی خود این دوتا هم باید اول ماکسیمم بشن تا توی کرانها قرار بگیرن.
L و J مجموع مقادیر انرژی تغذیه کننده ها هستن که ما اول باید این دوتا رو با روابط 16 و 17 بیشینه کنیم و بعد توی کرانهای بالا و پایین قرار بدیم.
ضمن اینکه خود این راندمای که میخواد بدست بیاد بصورت شرطی باید باشه. یعنی اگه لاندا از اون مقدار بزرگتر بود ، کران بالا و پایین راندمان به این صورته و اگر لاندا کوچکتر از اون مقدار بود ، راندمان باید توی این بازه قرار بگیره. یعنی حد بالا و پایین یا همون upperbound و lowerbound آن به این صورت میشه.
ـــــــــــــــــــــــــ ـــــــــــــــــــــــــ ـــــــــــــــــــــــــ ــــــــــــــــــــــــ
در مورد توابع fmincon و fminbnd ، خدمت شما عرض کنم که هر دو برای مسائل غیر خطی مقید شده استفاده میشن. ولی fminbnd برای یک متغیر استفاده شده و fmincon برای توابع چند متغیره استفاده میشه. پس نمیشه جای همدیگه هم استفاده کرد!
به نظرم توی این مسئله هم باید از fmincon استفاده کنیم.
شاید بشه برای مسائل یک متغیره هم از fmincon استفاده کرد فقط به جای چند متغیر از یک متغیر استفاده می کنیم ( من که متلب ندارم شما امتحان کن ببین میشه یا نه)نقل قول:
نوشته شده توسط coronaa
سوال: در اینجا یک تابع داریم که در اون تابع باید همزمان دو تا متغیر max بشن
یا اینکه دو تا تابع داریم که برای هرکدوم باید یک متغیر max بشه؟
ببینم چی کار می کنی نمودار را بکش
طبق روابط 16 و 17 دو تا تابع داریم که باید ماکزیمم بشن و هر تابع هم یک متغیر جدا داره. L و J دوتا متغیر هستن که L برای کران پایین و J برای کران بالا استفاده میشه.
مهندس من سوادم در همین حده! خب همون برنامه ای که بالا گذاشتم در همین حد بلدم. ران که میکنم یه سری متغیرها مثل لاندا رو نمیشناسه.
باید چیکار کنم؟
یه همچین چیزی
بعضی ها را کد خودت کپی کردم ( ضمانتش با خودت)
alpha=.3; beta=1.2; k=.38; D=6.3;lambda=1;
for m=2:10
parameterfunL=@(L,apha,k,m,D) (alpha-k)/m*(L^2)+k*D*L;
parameterfunJ=@(j,apha,k,m,D) (alpha/m-k)*(j^2)+k*D*j;
fl=@(L) parameterfunL(L,apha,k,m,D);
fj=@(j) parameterfunJ(j,apha,k,m,D);
Lmax=fminbnd(@fl,0,m*lambda);
Jmax=fminbnd(@fj,0,m*lambda);
if lambda > k*D/((m+1)*k-2*alpha)
lobnd=((m*k^2*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/(((alpha-k).*Lmax^2/m)+k*D*Lmax);
upbnd=((m*k^2*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/((alpha/m-k)*Jmax^2+k*m*D*Jmax);
end
if lambda < k*D/((m+1)*k-2*alpha)
lobnd= (m^2.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k).*Lmax.^2)+k*m*d*Lmax);
upbnd= (m.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k*m).*Jmax.^2)+k*m*d*Jmax);
end
end
برای خط 9 و 10 اشاره گر رو برداشتم. خطا میداد :
Lmax=fminbnd(fl,0,m*lambda);
Jmax=fminbnd(fj,0,m*lambda);
خب حالا این upbnd و lobnd رو بدست میاده. یک عدد هست. همین درسته یا باید یک بردار بده؟
بالایی اشتباه شدم. به ازای m های مختلف up,lo بدست میاد ولی نمیدونم چطوری باید توی یه بردار بذارم!
حل شد.
الآن به دوتا بردار تبدیل شد ولی عدداش یه طوریه!
upbnd
1.0e+004 *
Columns 1 through 7
0 1.6698 1.5303 0.7856 0.7421 0.2988 0.0941
Columns 8 through 10
0.0000 0.0000 0.0000
lobnd
1.0e+004 *
Columns 1 through 7
0 3.3397 4.5908 3.1424 3.7103 1.7930 0.6588
Columns 8 through 10
0.4396 0.3148 0.2507
مهندس دوتا برنامه رو با هم ترکیب کردم و بجای lp و ub از برنامه جدید استفاده کردم تا ببینم نمودارش چطوری میشه.
این کارم اشتباهه؟
فانکشن fffu
function ret = fffu(q)
alpha=.3; beta=1.2; k=.38; D=6.3;
ret=-sum(alpha.*q.^2+k.*(D-sum(q)).*q);
end
برنامه اصلی :
clc; close all; clear;
alpha=.3; beta=1.2; k=.38; D=6.3;lambda=1;
for m=2:10
parameterfunL=@(L,alpha,k,m,D) (alpha-k)/m*(L^2)+k*D*L;
parameterfunJ=@(j,alpha,k,m,D) (alpha/m-k)*(j^2)+k*D*j;
fl=@(L) parameterfunL(L,alpha,k,m,D);
fj=@(j) parameterfunJ(j,alpha,k,m,D);
Lmax=fminbnd(fl,0,m*lambda);
Jmax=fminbnd(fj,0,m*lambda);
if lambda > k*D/((m+1)*k-2*alpha)
lobnd(m)=((m*(k^2)*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/(((alpha-k).*Lmax^2/m)+k*D*Lmax);
upbnd(m)=((m*(k^2)*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/((alpha/m-k)*Jmax^2+k*m*D*Jmax);
end
if lambda < k*D/((m+1)*k-2*alpha)
lobnd(m)= (m^2.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k).*Lmax.^2)+k*m*D*Lmax);
upbnd(m)= (m.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k*m).*Jmax.^2)+k*m*D*Jmax);
end
end
fprintf('Lowerbound\n ');
disp(lobnd)
fprintf('\n---------------------------------------\n');
fprintf('Upperbound \n');
disp(upbnd)
% ---------------------------------------------------
% ---------------------------------------------------
payoff=[1.934 ,1.554 ,1.174 ,0.794 , 0.414 ,0.077012 ,0.057655 ,0.044775,0.035774];
soorat=payoff.*(2:10);
makhraj = zeros(1,10) ;
efficiency = zeros(1,10);
error =zeros(1,10);
for numplayers=2:10
lb=lobnd(1:numplayers);
ub=upbnd(1:numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5,[],[],[],[],lb,ub);
disp('makhraj=')
disp(-fx1)
makhraj(numplayers) =(-fx1);
efficiency(numplayers) = soorat(numplayers-1)/makhraj(numplayers);
error(numplayers) = abs(ub(numplayers-1)-efficiency)/efficiency;
end
plot(2:10,efficiency(2:10)), grid
ylabel('Symmetric NE Efficiency');
xlabel('Number of Suppliers in Microgrid')
axis([2 10 0 1.1])
% ----------------------------------------------
figure
plot(2:10,error(2:10)), grid
ylabel('relative error of the upperbound');
xlabel('Number of Suppliers in Microgrid')
% legend('cos_x','sin_x')
axis([2 10 0 1.1])
% ----------------------------------------------
یه پیغام خطا هم میده :
Exiting due to infeasibility: 1 lower bound exceeds the corresponding upper bound.
makhraj=
??? In an assignment A(I) = B, the number of elements in B and
I must be the same.
Error in ==> plot5 at 42
makhraj(numplayers) =(-fx1);
ظاهرن توی مرحله یک ، fx1 پوچ در میاد!
چیکار باید بکنم؟
مهندس من وقتی برای خط پایینی مقدار lb و ub رو بصورت زیر انتخاب میکنم، واسه اولین مرحله( یعنی m برابر 2)، fx1 رو پوج میده. این کارم اشتباهه؟
lb=lobnd(1:numplayers);
ub=upbnd(1:numplayers);
[x1, fx1]=fmincon(@fffu,ones(1,numplayers)*.5,[],[],[],[],lb,ub);
میخواستم طبق شکل (4)، نمودار حد بالا و پایین رو بدست بیارم.
خواهش میکنم یه کمکی بکن.
ین یکی را امتحان کن
alpha=.3; beta=1.2; k=.38; D=6.3;lambda=1;
ll=zeros(1,10);
uu=zeros(1,10);
for m=2:10
parameterfunL=@(L,alpha,k,m,D) -(alpha-k)/m*(L^2)+k*D*L;
parameterfunJ=@(j,alpha,k,m,D) -(alpha/m-k)*(j^2)+k*D*j;
fl=@(L) parameterfunL(L,alpha,k,m,D);
fj=@(j) parameterfunJ(j,alpha,k,m,D);
Lmax = fmincon(@fl,m*lambda/2,[],[],[],[],0,m*lambda);
Jmax = fmincon(@fj,m*lambda/2,[],[],[],[],0,m*lambda);
if lambda > k*D/((m+1)*k-2*alpha)
lobnd=((m*k^2*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/(((alpha-k).*Lmax^2/m)+k*D*Lmax);
upbnd=((m*k^2*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/((alpha/m-k)*Jmax^2+k*m*D*Jmax);
end
if lambda < k*D/((m+1)*k-2*alpha)
lobnd= (m^2.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k).*Lmax.^2)+k*m*D*Lmax);
upbnd= (m.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k*m).*Jmax.^2)+k*m*D*Jmax);
end
ll(m)=lobnd;
uu(m)=upbnd;
end
من توی خط 12 و 13 اشاره گرها رو بر میدارم. وگرنه خطا میده:
جواب :
ll =
0 Inf Inf Inf Inf Inf Inf Inf Inf Inf
>> uu
uu =
1.0e+015 *
Columns 1 through 7
0 Inf Inf 2.2085 Inf Inf Inf
Columns 8 through 10
0.0542 Inf Inf
دو تا فایل برای این تابعها درست کن
function ret=fj(j)
global m;
alpha=.3; beta=1.2; k=.38; D=6.3;lambda=1;
ret= -(alpha/m-k)*(j^2)+k*D*j;
end
function ret=fl(L)
global m;
alpha=.3; beta=1.2; k=.38; D=6.3;lambda=1;
ret= -(alpha-k)/m*(L^2)+k*D*L;
end
بعدش هم این یکی را اجرا کنچ
global m;
alpha=.3; beta=1.2; k=.38; D=6.3;lambda=1;
ll=zeros(1,10);
uu=zeros(1,10);
for m=2:10
Lmax = fmincon(@fl,m*lambda/2,[],[],[],[],0,m*lambda);
Jmax = fmincon(@fj,m*lambda/2,[],[],[],[],0,m*lambda);
if lambda > k*D/((m+1)*k-2*alpha)
lobnd=((m*k^2*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/(((alpha-k).*Lmax^2/m)+k*D*Lmax);
upbnd=((m*k^2*D^2*(k-alpha))/(((m+1)*k-2*alpha))^2)/((alpha/m-k)*Jmax^2+k*m*D*Jmax);
end
if lambda < k*D/((m+1)*k-2*alpha)
lobnd= (m^2.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k).*Lmax.^2)+k*m*D*Lmax);
upbnd= (m.*(alpha.*lambda.^2+k.*(D-m*lambda).*lambda))/(((alpha-k*m).*Jmax.^2)+k*m*D*Jmax);
end
ll(m)=lobnd;
uu(m)=upbnd;
end
جواب :
http://uplood.ir/53lR