%編寫懲罰函數 pen.m
function y=pen(x)
s=0;
for i=1:1:length(x)
   s=s+x(i)
end
y=min(0,min(x)).^2+(1-s).^2;

%編寫交易費用函數h.m
function y=h(x,u,p)
if (x>u)
   y=x*p
else
   if (x>0)
      y=u*p
   else
      y=0;
   end
end

%編寫目標函數f.m
function y=f(x)
global w;
global M;

r=[0.28,0.21,0.23,0.25,0.05];
q=0.01*[2.5,1.5,5.5,2.6,0];
p=0.01*[1,2,4.5,6.5,0];
u=[103,198,52,40,0];

penalty=1e-25;
t=0;
for i=1:1:4
   t=t+h(x(i)*M,u(i),p(i));
end

t=t-M*sum(r.*x);
y=w*t-(1-w)*max(x.*q)*M+penalty*pen(x);

for k=1:1:11
   w=(k-1)*0.1;
   s=5000;
   for j=1:1:50
      c=rand(1,5);
      d=sum(c);
      x0=c/d;
      x=fminsearch('f',x0);
      if(f(x)<s)
         s=f(x);
         xc=x;
      end
   end
   t=0;
   for i=1:1:4
      t=t+h(xc(i)*M,u(i),p(i));
   end
      FR(k)=f(xc);
   net(k)=t+M*sum(r.*xc);
   risk(k)=M*max(xc.*q);
   xn(k,:)=xc;
end
figure(1)
plot(w,xn(:,1),'*',w,xn(:,2),'c',w,xn(:,3),'s',w,xn(:,4),'v',w,xn(:,5),'0');