Economical Optimization w/o G limits.

eos.m

clc; clearvars; 

lambda = 6;
a = [0.004 0.006 0.009];
b = [5.3 5.5 5.8];
c = [500 400 200];
Pd = 800;
error = 100;
tolerance = 0.001;
k = 0;

while(error > tolerance)
    k = k + 1;
    P(k,:) = (lambda - b)./(2*a);
    error = abs(sum(P(k,:))-Pd);
    deltaLambda = error/sum(1./(2*a));
    lambda = lambda + deltaLambda;
end
Ci = sum(a.*P(k,:).^2+b.*P(k,:)+c);

fprintf('Total load: \t\t%.3f MW\n', Pd);
disp('Load shared by each unit (in MW): ');
disp(P(k,:));
fprintf('Incremental cost: \t$ %.3f \n', lambda);
fprintf('Total cost: \t\t$ %.3f /h\n', Ci);