main.m

电网机制设计

clc,clear
global beta r
global Tlimit
global alpha n miu miuo lambda
lambda = 0.75;
alpha = 1.2;
beta  = 0.5;
r     = 0.01;
n     = 10 ;
capa = 12;
miu   = n.^(alpha-1)/(n.^(alpha-1)+2);
miuo   = n.^(alpha-1)/(n.^(alpha-1)+1);

% find the Tlimit

findoptim
% define the global paramters

miu, miuo
% the right of the certain inequation 
rightvalue = sqrt(pi).*miu./(2.*n.^(alpha-1).*r.^(3/2))

% define the temp paramters 
a = 0;% 积分下限
b = Tlimit;% 积分上限
m = 800;
tol = 0.0001;

% computer 
[R,quad,err,h]=romber('quadequ',a,b,m,tol);
quadvalue = quad


% cost
capa = 12;% A
capdt = 0.9;% D（T）
Ilow = n.*(2.*r.^(3/2)./(sqrt(pi).*miu.*capa)).^(1./(alpha-1));%开a-1次方是^(alpha-1)吗？
%Ilow = n.pow((2.*r.^(3/2)./(sqrt(pi).*miu.*capa)),1/(alpha-1))
bb = r.*Tlimit;
[R1,quad1,err1,h1]=romber('quadequo',a,bb,m,tol);
thequadone = quad1;
cost = capa.*Ilow.^alpha.*(miu.*sqrt(pi)./(2.*n.^(alpha-1))-thequadone)./(r.^(3/2).*capdt);
Iinterval = linspace(Ilow,5*Ilow,500);% 
costinterval = capa.*Iinterval.^alpha.*(miu.*sqrt(pi)./(2.*n.^(alpha-1))-thequadone)./(r.^(3/2).*capdt);
plot(Iinterval,costinterval)
xlabel('投资（美金）'); ylabel('最优预估价（美金）');
title('国家电网公司的最优预估价随成本的变化')
% legend(' I=100')

aa = 0;% 下限
bb = 100;% 上限
delta = 0.00001;
[c,err,yc]=bisect('highequ',aa,bb,delta);
Icost=c;
findequationshigh