A=[200 180 140]; B=[7 6.3 6.8]; G=[0.008 0.009 0.007]; BL=[0.000218 0.000093 0.000028 0.000093 0.000228 0.000017 0.000028 0.000017 0.000179]; Pd=170; Ng=3; lambda=7; dp=10; dlambda=.01; Pmin=[10 10 10]; Pmax=[85 80 70]; P=[0 0 0]; while abs(dp) >= .001 SS=2*BL*P'; for i=1:Ng; for j=1:Ng; S(i) = SS(i)-(2*BL(i,j)*P(i)); P(i)=(lambda*(1-S(i)))-B(i)/2*(G(i)+lambda*BL(i,i)); if P(i)>Pmax(i) P(i)=Pmax(i); elseif P(i)<Pmin(i) P(i)=Pmin(i); end end end pt=sum(P); PL=P(i)*BL(i,j)*P(j); for i=Ng; dp=P(i)-PL-Pd; end if dp>0 lambda=lambda-dlambda; elseif dp<0 lambda=lambda+dlambda; end end Pd P lambda Total_Cost=sum(A+B.*P+G.*P.^2)